на главную ] 

Восстановление смеси плотностей распределений с помощью EM-алгоритма.

Евгений Борисов

четверг, 26 декабря 2013 г.


В этой статье мы поговорим о задаче восстановления смеси плотностей распределений точек и о методе, который называется EM-алгоритм (expectation–maximization algorithm). Этот метод можно применять для решения разных практических задач, например для решения задачи кластеризации.

 
 

1 Введение

Рассмотрим множество точек X, поэтому набору точек требуется восстановить плотность их распределения. Пусть плотность имеет следующий вид

       ∑ k               ∑k
p(x) =     wjφ (x;𝜃j) ;      wj = 1  ; wj ≥  0
        j=1               j=1
где 𝜃j - параметры, φ - заданная функция, k - количество компонент смеси, wj - весовой коэффициент (априорная вероятность компоненты)

Задача: найти оптимальные 𝜃j,wj по простой выборке X и заданным k и φ.

 
 

2 Базовый вариант EM-алгоритма

Для решения описанной выше задачи можно использовать EM-алгоритм, который состоит из двух основных процедур: E-шаг - оценка параметров и М-Шаг - оптимизация параметров. Базовый вариант EM-алгоритма выглядит следующим образом.

EM(X,k,w,𝜃,δ)

  1. выбрать начальные значения w и 𝜃
  2. stepE: вычислить оценку параметров
    g′ij := gij

          ---wj-φ(xi,𝜃j)---
gij := ∑k   w  φ(x ,𝜃 )
         s=1  s   i  s

    i = 1,,m - номер объекта в выборке,
    j = 1,,k - номер компоненты смеси распределения

  3. stepM: оптимизировать параметры

              m
      -1 ∑
wj := m      gij
         i=1
                 ∑m
𝜃j := arg max    gij ln φ (xi;𝜃)
          𝜃   i=1
    (1)

  4. если max i,j|gij gij| > δ
    то переход на п.2

    где 0 < δ < 1 – максимальное значение изменения параметров (условие остановки EM)

  5. конец

Это базовый вариант метода. Результат его работы зависит от начальных значений параметров, здесь не вполне ясно как выбирать количество компонент k и начальные параметры 𝜃j. Для решения этих задач существует модификация базового алгоритма - EM-алгоритм с последовательным добавлением компонент.  

 
 

3 EM-алгоритм с последовательным добавлением компонент

Идея этой модификации ЕМ-алгоритма вполне описывается в названии. Мы последовательно добавляем компоненты к смеси, на каждом шаге выполняем базовый ЕМ и оцениваем качество результата по определённой процедуре. Этот вариант EM-алгоритма выглядит следующим образом.

  1. выбрать начальные значения параметров

    k = 1, w = 1

    параметры 𝜃 выбираются в зависимости от вида функции φ,
    например, если φ это многомерная нормальная (гауссовская) плотность
    то 𝜃1 = (w11, Σ1), где μ1 - мат.ожидание, Σ1 - матрица ковариаций на X.

  2. найти максимальное значение плотности p на X

    pmax  = mxa∈xX p (x )
           ∑ k
p(x) =     wjφ (x;𝜃j)
        j=1
  3. выбрать из X точки с плотностью ниже пороговой R pmax,

    где 0 < R < 1 – коэффициент правдоподобия

    U  = {x ∈ X |p(x) < R ⋅ pmax }
  4. если |U| < m0
    то переход на п.8

    где 0 < m0 < |X| – минимальная длинна выборки для обработки

  5. добавляем компоненту смеси

    k = k + 1

    mk  = |U|

    w  =  mk- ;  w  = w  (1 − w ) ;  j = 1,...k − 1
 k    m       j     j      k

    параметры для нормальной (гауссовской) плотности φ

          w  ∑
μk = --k     xi
     mk  xi∈U
    Σk  = cov(U )
  6. выполняем базовую процедуру оптимизации

    EM  (X, k,w, 𝜃,δ)

    где 0 < δ < 1 - максимальное значение изменения параметров(условие остановки EM)

  7. переход на п.2
  8. конец

 
 

4 Реализация

В этом разделе мы рассмотрим реализации базового и последовательного EM-алгоритмов для точек на плоскости (n = 2).

Определим φ(x; 𝜃) как многомерную нормальную (гауссовская) плотность


Рис.1: функция Гаусса
                           (                      )
                       exp--−-12(x-−-μ)TΣ-−1(x-−-μ)--
φ (x;𝜃) = N (x;μ,Σ ) =         ∘ ----n------
                                 (2π)  detΣ
где μ n - мат.ожидание (центр) x X, Σ n×n матрица ковариаций X. Таким образом 𝜃 j это совокупность (μj, Σj) мат.ожидания и матрицы ковариаций для компоненты смеси j.

В этом случае решение задачи (1) по нахождению оптимальных 𝜃j = (μj, Σj) выглядит следующим образом (j = 1,,k)

(
|       ∑m
|||| Nj  =     gij
|||       i=1
|{        1 ∑m
   μj = ---    gijxi
|||       Nj  i=1
|||          ∑m
|||(  Σj = -1-    gij(xi − mj )(xi − mj )T
        Nj  i=1

 
 

4.1 Реализация базового EM-алгоритма

В этом примере будем искать смесь из двух компонент (k = 2).

Начальные значения μ1 - мат.ожидание на X, μ2 генерируется случайным образом, начальные матрицы ковариаций Σ1 = Σ1 = cov(X)

Отметим, что результат зависит от начальных значений параметров, и выбранный способ инициализации является не самым удачным, в данном случае он был выбран из-за его простоты.

На рисунках ниже проиллюстрирована работа ЕМ-алгоритма, который восстанавливает смесь плотностей.

 

Рис.2: начальное состояниеРис.3: результат работы ЕМ

 
Рис.4: изменение параметра g Рис.5: анимация процесса (кликабельно)
 
 

Реализация в системе Octave здесь  ].

 
 

4.2 Реализация последовательного EM-алгоритма

На рисунках ниже проиллюстрирована работа последовательного ЕМ-алгоритма, который восстанавливает смесь плотностей.

 

Рис.7: начальное состояние Рис.8: промежуточный результат работы ЕМ (две компоненты)

 
Рис.9: промежуточный результат работы ЕМ (три компоненты) Рис.10: окончательный результат работы ЕМ (четыре компоненты)
 
  На рисунках ниже проиллюстрировано изменение количества точек с низким значением плотности (отмечены желтым цветом) на разных этапах выполнения алгоритма.
Рис.11: начальное состояние Рис.12: промежуточный результат работы ЕМ (две компоненты)

 
Рис.13: промежуточный результат работы ЕМ (три компоненты) Рис.14: окончательный результат работы ЕМ (четыре компоненты)

 
  На рисунках ниже проиллюстрировано изменение параметра g на разных этапах выполнения алгоритма.

Рис.15: промежуточный результат работы ЕМ (две компоненты) Рис.16: промежуточный результат работы ЕМ (три компоненты)

 
 
Рис.17: промежуточный результат работы ЕМ (четыре компоненты)  

 
 

Реализация в системе Octave здесь  ].

 
 

Список литературы

[1]    Воронцов К.В. Статистические методы классификации – http://shad.yandex.ru/lectures/machine_learning.xml

[2]    GNU Octave – http://www.gnu.org/software/octave/

Яндекс.Метрика
TOP.zp.ua PR-CY.ru
При использовании материалов этого сайта, пожалуйста вставляйте в свой текст ссылку на мою статью.