Изменения

Перейти к: навигация, поиск

EM-алгоритм

11 738 байт добавлено, 19:36, 4 сентября 2022
м
rollbackEdits.php mass rollback
== Определение ==
'''Алгоритм EM''' (англ. ''expectation- maximization'') {{---}} итеративный алгоритм поиска оценок максимума правдоподобия параметров для решения задачмодели, в ситуации, где некоторые переменные не являются наблюдаемымикогда она зависит от скрытых (ненаблюдаемых) переменных.
Алгоритм ищет параметры модели итеративно, каждая итерация состоит из двух шагов:
'''E(Expectation)''') шаг, в котором находится распределение {{---}} поиск наиболее вероятных значений скрытых переменных используя значение наблюдаемых переменных и текущего значения параметров.
'''M(MaximisationMaximization)''' шаг {{- пересчет --}} поиск наиболее вероятных значений параметров, находя максимум правдоподобия исходя из распределения скрытых переменных, для полученных на шаге E - шагезначений скрытых переменных.
== Задача разделения смеси распределений == EM алгоритм подходит для решения задач двух типов:
=== Общий алгоритм ===# Задачи с неполными данными.# Задачи, в которых удобно вводить скрытые переменные для упрощения подсчета функции правдоподобия. Примером такой задачи может служить кластеризация.
Необходимо описать плотность распределения функции на X как сумму k функций, которые можно рассматривать как элементы параметрического семейства функций <tex> p_j(x) = \phi(x;\theta_j)</tex>. Плотность распределения будет выглядеть как <br><tex>p(x) = \sum\limits_{iОсновной алгоритм =1}^k \omega_j p_j(x); \sum\limits_{i=1}^k w_j = 1; w_j >= 0 </tex> <br /> где <tex>\omega_j</tex>- априорная вероятность j компоненты распределения.Задача разделения смеси заключается в том, чтобы, имея выборку <tex>X^m</tex> случайных и независимых наблюдений из смеси <tex>p(x)</tex>, зная число <tex>k</tex> и функцию <tex>\phi</tex>, оценить вектор параметров <tex>\theta = (\omega_1,..,\omega_k,\theta_1,..,\theta_k)</tex>
=== Постановка задачи === Плотность распределения смеси имеет вид:<br/><tex>p(x) = \sum\limits_{j=1}^k w_j p_j(x)</tex>.<br/>Где <tex> \sum\limits_{j=1}^k w_j = 1; w_j \geq 0; p_j(x) = \phi(x;\theta_j)</tex> {{---}} функция правдоподобия <tex>j</tex>-ой компонеты смеси, <tex>w_j</tex> {{---}} априорная вероятность <tex>j</tex>-ой компоненты смеси.<br/> Перед нами стоит две задачи:<br/> # По заданной выборке <tex>X^m</tex> случайных и независимых наблюдений полученных из смеси <tex>p(x)</tex>, числу <tex>k</tex> и функции <tex>\phi</tex>, оценить вектор параметров <tex>\Theta = (w_1,..,w_k,\theta_1,..,\theta_k)</tex>.# Найти <tex>k</tex>. === Проблема ===  Задачи подобного рода мы умеем решать, максимизируя логармиф правдоподобия:<br><tex> Q(\Theta) = ln \prod\limits_{i=1}^mp(x_i) = \sum\limits_{i=1}^m ln\sum\limits_{j=1}^k w_j p_j(x_i; \theta_j) \longrightarrow \max\limits_{\Theta}</tex>.<br/> Но проблeма в том, что мы не знаем как аналитически посчитать логарифм суммы. Тут нам и поможет алгоритм EM. === Решение ===  Основная идея алгоритма EM заключается в том, что мы добавляем скрытые переменные такие, что:<br/> # Они могут быть выражены через <tex>\Theta</tex>.# Они помогают разбить сумму так: <tex>p (X, H|\Theta) = \prod\limits_{i=1}^k p (X|H, \Theta) p(H|\Theta)</tex>, где <tex>H</tex> {{---}} матрица скрытых переменных. Тогда алгоритм EM сводится к повторению шагов, указанных в [[#Определение|Определении]]. === E-шаг:===  <tex>p(x,\theta_j) = p(x)P(\theta_j | x) = w_jp_j(x)</tex>.<br />
Скрытые переменные представляют из себя матрицу <tex>p(x,\theta_j) H = p(xh_{ij})P(_{m \theta_j | x) = w_jp_j(x)times k}</tex> ,<br />Введем обозначение: где <tex> g_h_{ij} = P(\theta_j | x_i) </tex> это и будут скрытые параметры данной задачи {{-- апостериорная -}} вероятность того, что обучающий объект <tex> x_i </tex> получен из пренадлежит <tex>j</tex>-й компоненты ой компоненте.<br/>
По формуле Байеса справедливо равенство:<br />
<tex> g_h_{ij} = \frac{w_jp_j(x_i)}{p (x_i)} = \frac{w_jp_j(x_i)}{\sum\limits_{ts=1}^k w_t p_tw_s p_s(x_i)}</tex>.<br /> Также <tex>\sum\limits_{j=1}^k h_{ij} = 1</tex>.<br/> Таким образом при , зная значение значения вектора параметров <tex>\Theta</tex>, мы легко найти можем пересчитать значения скрытых переменных.<br/> === M-шаг ===  {{Теорема|statement=Если известны скрытые переменные, то задача минимизации <tex>Q(\Theta)</tex> сводится к <tex>k</tex> независимым подзадачам:<br/><center><tex>\theta_j = \arg\max\limits_{\theta}\sum\limits_{i=1}^m h_{ij}*\ln\phi(x_i;\theta)</tex>.</center>Оптимальные же веса считаются как:<br/><center><tex> w_j = \frac {1} {m} \sum\limits_{i=1}^m h_{ij}</tex>.</center>|proof=Посчитаем логарифм правдоподобия:<br><tex> Q(\Theta) = \sum\limits_{i=1}^m ln\sum\limits_{j=1}^k w_j p_j(x_i; \theta_j) \longrightarrow \max\limits_{\Theta}</tex>.<br/>При условии, что<tex> \sum\limits_{j=1}^k w_j = 1; w_j \geq 0</tex> имеет смысл рассматривать Лагранжиан задачи:<br/><tex> L(\Theta, X^m) = \sum\limits_{i=1}^m ln \biggl( \sum\limits_{j=1}^k w_j p_j(x_i) \biggr) - \lambda \biggl(\sum\limits_{j=1}^k w_j - 1 \biggr) </tex>.<br/>Приравняв нулю производную Лагранжиана по <tex>w_j</tex>, получим:<br/><tex>\frac{\partial L} {\partial w_j} = \sum\limits_{i=1}^m \frac{p_j(x_i)}{\sum\limits_{s=1}^kw_s p_s(x_i)} - \lambda = 0, j = 1..k</tex>.<br/>Умножим на <tex>w_j</tex> и просуммируем уравнения для всех <tex>j</tex>:<br /><tex>\sum\limits_{j=1}^k \sum\limits_{i=1}^m \frac{w_jp_j(x_i)}{\sum\limits_{s=1}^kw_s p_s(x_i)} = \lambda \sum\limits_{j=1}^kw_j</tex>.<br />А так как <tex>\sum\limits_{j=1}^k \frac{w_jp_j(x_i)}{\sum\limits_{s=1}^kw_sp_s(x_i)} = 1</tex> и <tex>\sum\limits_{j=1}^kw_j = 1</tex>, из чего следует <tex>\lambda = m</tex>.<br /> <tex>w_j = \frac{1}{m}\sum\limits_{i=1}^m \frac{w_jp_j(x_i)}{\sum\limits_{s=1}^kw_sp_s(x_i)} = \frac{1}{m}\sum\limits_{i=1}^m h_{ij}</tex>.<br /> Приравняв к нулю производную Лагранжиана по <tex>\theta_j</tex>, схожим способом найдем:<br /> <tex> \theta_j = \arg\max\limits_{\theta}\sum\limits_{i=1}^m h_{ij}*\ln\phi(x_i;\theta).</tex><br />  }} === Критерий остановки ===  Алгоритм EM выполняется до сходимости, но как нам определить, что сходимость наступила? Мы можем останавливаться, когда либо <tex>Q(\Theta)</tex>, либо <tex>H</tex> перестают сильно меняться. Но, обычно, удобней контролировать изменения значений скрытых переменных, так как они имеют смысл вероятностей и принимают значения из отрезка <tex>[0,1]</tex>. Поэтому один из возможных критериев остановки будет выглядеть так: <tex>\max\limits_{i,j} |h_{ij} - h_{ij}^{(0)}| > \delta</tex>. === Псевдокод ===   Input:<tex>X^m, k, \Theta^{(0)}</tex> Repeat '''E-step''': for all i = 1..m; j = 1..k: <tex>h_{ij} = \frac{w_j \phi(x_i; \theta_j)}{\sum\limits_{s=1}^k w_s \phi(x_i; \theta_j)}</tex> '''M-step''': for all j = 1..k: <tex>\theta_j = \arg\max\limits_{\theta} \sum\limits_{i=1}^m h_{ij}*ln \phi (x_i, \theta)</tex> <tex>w_j = \frac {1} {m} \sum\limits_{i=1}^m h_{ij}</tex> Until a stopping criterion is satisfied Return <tex>\Theta = (\theta_j, w_j)_{j=1}^k</tex> === Плюсы и минусы ===  Плюсы:<br/> * Сходится в большинтсве случаев.* Наиболее гибкое решение.* Существуют простые модификации, позволяющие уменьшить чуствительность алгоритма к шуму в данных. Минусы:<br/> * Чуствителен к начальному приближению. Могут быть ситуации, когда сойдемся к локальному экстремуму.* Число компонент <tex>k</tex> является [[Настройка_гиперпараметров|гиперпараметром]]. == Модификации ==  Базовый алгоритм EM является очень гибким для модификаций, позволяющих улучшить его работу. В этом разделе мы приведем краткое описание некоторых из них. === Generalized EM-algorithm ===  Осоновная идея этой модификации заключается в том, что на шаге M мы не будем пытаться найти наилучшее решение. Это применимо в случаях, когда максимизация <tex>Q(\Theta)</tex> является сликшом дорогой, поэтому нам достаточно сделать лишь несколько итераций, для того, чтобы сместиться в сторону максимума значения <tex>Q(\Theta)</tex>. Эта модификация имеет неплохую сходимость. === Stochastic EM-algorithm === Как уже было отмечено в [[#Плюсы_и_минусы|Плюсы и минусы]], базовый алгоритм чувствителен к начальному приближению и могут быть ситуации, когда алгоритм "застрянет" в локальном экстремуме. Для того, чтобы предотвратить это, будем на каждой итерации алгоритма случайно "встряхивать" выборку. В этой модификации у нас добавляется шаг S, на котором мы и будем "встряхивать" выборку. И на шаге M мы будем решать уже задачу максимуму невзвешенного правдоподобия. Эта модификация хороша тем, что нечуствиетльная к начальном приблежению.
Перейдем к M-шагу:== Пример. Разделение смеси Гауссиана ==
Посчитаем для аддитивности логарифм правдоподобия[[Файл:<br /><tex> Q(\Theta) = ln \prod\limits_{i=1}^mp(x_i) = \sum\limits_{i=1}^m ln\sum\limits_{j=1}^k w_j p_j(x_i) \longrightarrow max</tex> <br />при условии <tex>\sum\limits_{i=1}^k w_j = 1; w_j >= 0</tex> имеет смысл рассматривать лагранжиан задачи:<br /><tex>\frac{\partial L} {\partial w_j} = \sum\limits_{i=1}^m \frac{p_j(x_i)}{\sum\limits_{t=1}^kw_tp_t(x_i)} - \lambda = 0</tex><br />Gaussians2.png|right|thumb|400px|Несколько итераций алгоритма]]
Умножим на <tex>\omega_j</tex> Каноническим примером использования EM алгоритма является задача разделения смеси гауссиана. Данные у нас получены из нормального распределения. В этом случае параметрами функций ялвяются матожидание и просумируем уравнения для всех <tex>j</tex> дисперсия.<br />
<tex>\sumtheta = (w_1,..,w_k;\;\mu_1,..,\mu_k;\;\sigma_1,..,\limits_sigma_k)</tex> {j{---}} вектор параметров, <br/><tex>p_j(x) =1}^k N(x;\summu_j, \limits_sigma_j) = \frac1{i=1\sqrt{2\pi}\sigma_j}^m \exp \biggl(-\frac{w_jp_j(x_ix - \mu_j)^2}{2\sum\limits_{t=1}sigma_j^kw_tp_t(x_i)2} = \lambda \sum\limits_{j=1}^kw_jbiggr) </tex> {{---}} плотность распределения.<br />
Так как можно заменить порядок суммы и <tex>\sum\limits_{i=1}^m \frac{w_jp_j(x_i)}{\sum\limits_{t=1}^kw_tp_t(x_i)} = 1</tex> и <tex>\sum\limits_{j=1}^kw_j = 1</tex> из чего следует <tex>\lambda = m</tex> Посчитаем значения для каждого шага. <br />
<tex>\omega_j = \frac{1}{m}\sum\limits_{i=1}^m \frac{w_jp_j(x_i)}{\sum\limits_{t=1}^kw_tp_t(x_i)} = \frac{1}{m}\sum\limits_{i=1}^mg_{ij}</tex><br />E-шаг:
Приравняв к нулю лагранжиан по : <tex>h_{ij} = \frac{w_j N(x_i, \mu_j, \sigma_j)}{\theta_jsum\limits_{s=1}^k w_s N(x_i, \mu_s, \sigma_s)}.</tex> схожим способом найдем:<br />
<tex> \theta_j = arg \max\limits{\theta}\sum\limits_{i=1}^mg_{ij}ln(\phi(x_i;\theta)) </tex><br />M-шаг:
Таким образом на M-шаге необходимо взять среднее значение : <tex>w_j = \frac{1}{m} \sum\limits_{i=1}^m h_{ij}.</tex>: <tex>g_\mu_j = \frac {1} {mw_j} \sum\limits_{i=1}^m h_{ij}x_i.</tex> и решить : <tex> \sigma_j^2 = \frac {1} {mw_j} \sum\limits_{i=1}^m h_{ij}(x_i - \mu_j)^2, j = 1..k независимых оптимизационных задач.</tex>
== Использование в задаче кластеризации ==
=== Разделение смеси гауссиан ===[[Файл:Gaussianskmeans.pngjpg|right|thumb|200px| несколько итераций алгоритмаПример работы k-means]]Важным на практике примером является случай, когда параметрическое семейство - нормальные распределения. Параметрами функций будут являться матожидание и дисперсия.<br /><tex>\theta = (w_1,..,w_k;\;\mu_1,..,\mu_k;\;\sigma_1,..,\sigma_k)</tex> — вектор параметров, <br /><tex>p_j(x) = N(x;\mu_j, \sigma_j) = \frac1{\sqrt{2\pi}\sigma_j} \exp \biggl(-\frac{(x - \mu_j)^2}{2\sigma_j^2}\biggr) </tex>
== Как уже упоминалось в [[#Определение|Определении]], алгоритм EM подходит для решения задачи кластеризации. И одной из его имплементаций для этой задачи является алгоритм [[Алгоритм_k-Means|<tex>k</tex>-means как EM алгоритм ==Means]]. В этом алгоритме в качестве скрытых переменных выступают метки классов объектов. Параметрами же являются центроиды искомых классов. Тогда на шаге E мы относим объекты к какому-то одному классу на основе расстояний до центроид. А на шаге M мы пересчитываем центроиды кластеров, исходя из полученной на шаге E разметке.<br/>
Скрытыми переменными в данной задаче являются классыТакже стоит упомянуть алгоритм <tex>c</tex>-means<ref>[https://en.wikipedia.org/wiki/Fuzzy_clustering#Fuzzy_C-means_clustering C-means clustering, Wikipedia]</ref>. В нем качестве скрытых переменных выступают вероятности принадлежности объекта к которым относятся объекты для кластеризации. Сами же параметры это центры масс классовклассам. На шаге E - распределяются все объекты по классам исходя мы пересчитывем вероятности принадлежности объектов, иходя из расстояния от центрадо центроид. Шаг M, идейно, на шаге M находится оптимальное месторасположение центраостается без изменений.
Аналогично рассматривается и алгоритм c-means. Скрытые переменные здесь будут вероятности принадлежности к классам, которые находятся на E-шаге по расстоянию от центра. Центр так же рассчитывается на M-шаге исходя из скрытых переменных.
 
 
 
== Реализация на python ==
 
В пакете sklearn алгоритм EM представлен объектом GaussianMixture. Проиллюстрируем его работу на примере задачи кластеризации и сравним его с алгоритмом <tex>k</tex>-means:
[[Файл:em_clustering.png|thumb|600px|Результат выполнения программы]]
'''import''' numpy as np
'''import''' matplotlib.pyplot as plt
'''from''' sklearn '''import''' cluster, datasets, mixture
'''from''' sklearn.preprocessing '''import''' StandardScaler
'''from''' itertools '''import''' cycle, islice
np.random.seed(12)
<font color="green"># Создаем datasets с использованием стандартных sklearn.datasets</font>
n_samples = 2000
random_state = 170
noisy_circles = datasets.make_circles(n_samples=n_samples, factor=.5, noise=.05)
noisy_moons = datasets.make_moons(n_samples=n_samples, noise=.05)
blobs = datasets.make_blobs(n_samples=n_samples, random_state=8)
varied = datasets.make_blobs(n_samples=n_samples, cluster_std=[1.0, 2.5, 0.5], random_state=random_state)
<font color="green"># Создаем анизатропно разделенные данные</font>
X, y = datasets.make_blobs(n_samples=n_samples, random_state=random_state)
transformation = [[0.6, -0.6], [-0.4, 0.8]]
X_aniso = np.dot(X, transformation)
aniso = (X_aniso, y)
<font color="green"># Выставляем параметры для matplotlib.pyplot</font>
plt.figure(figsize=(9 * 2 + 3, 12.5))
plt.subplots_adjust(left=.02, right=.98, bottom=.001, top=.96, wspace=.05, hspace=.01)
plot_num = 1
defaul_n = 3
<font color="green"># Варьируем значение количества классов в зависимости от данных, ведь для нас это гиперпараметр</font>
datasets = [
(varied, defaul_n),
(aniso, defaul_n),
(blobs, defaul_n),
(noisy_circles, 2)]
for i_dataset, (dataset, n_cluster) in enumerate(datasets):
X, y = dataset
<font color="green"># Нормализация данных</font>
X = StandardScaler().fit_transform(X)
<font color="green"># Непосредственно наш алгоритм - Gaussian Mixture</font>
gmm = mixture.GaussianMixture(n_components=n_cluster, covariance_type='full')
<font color="green"># Для сравнения берем алгоритм - k-means</font>
two_means = cluster.KMeans(n_clusters=n_cluster)
clustering_algorithms = {
'Real distribution': None,
'Gaussian Mixture': gmm,
'k-Means': two_means
}
for name, algorithm in clustering_algorithms:
# Этап обучения
if algorithm is not None:
algorithm.fit(X)
# Применяем алгоритм
y_pred = y if algorithm is None else algorithm.predict(X)
# Рисуем результаты
plt.subplot(len(datasets), len(clustering_algorithms), plot_num)
if i_dataset == 0:
plt.title(name, size=18)
colors = np.array(list(islice(cycle(['#377eb8', '#ff7f00', '#4daf4a']), int(max(y_pred) + 1))))
plt.scatter(X[:, 0], X[:, 1], s=10, color=colors[y_pred])
plt.xlim(-2.5, 2.5)
plt.ylim(-2.5, 2.5)
plt.xticks(())
plt.yticks(())
plot_num += 1
plt.show()
 
Как и следовало ожидать, алгоритм EM работает на некоторых данных лучше чем k-means, однако есть данные, с которыми он не справляется без дополнительных преобразований.
== См. также ==
*[[Кластеризация]]
*[[Алгоритм_k-Means|Алгоритм k-Means]]
 
==Примечания==
<references />
== Источники информации ==
# Материалы лекции про кластеризацию курса "Машинное обучение" университета ИТМО, 2019 год
# [http://www.machinelearning.ru/wiki/images/6/6d/Voron-ML-1.pdf Математические методы обучения по прецедентам К. В. Воронцов]
# [http://www.machinelearning.ru/wiki/index.php?title=EM-%D0%B0%D0%BB%D0%B3%D0%BE%D1%80%D0%B8%D1%82%D0%BC Статья про EM-алгоритм на machinelearning.ru]
# [https://machinelearningmastery.com/expectation-maximization-em-algorithm/ A Gentle Introduction to Expectation-Maximization]
# [http://dendroid.sk/2011/05/09/k-means-clustering/ k-means]
[[Категория:Машинное обучение]][[Категория: Кластеризация]]
1632
правки

Навигация