Цепи Маркова (1121219), страница 87
Текст из файла (страница 87)
Чтобы выполнить итерацию, соответствующую E-шагу,положимf (x; ), x ∈ X , y = Ψ (x) ∈ Y ;(3.9.4)h(x | y; ) =g(y; )здесь h(x | y; ) — условная плотность с.в. X при условии что Ψ (X) == y. Чтобы максимизировать ln g(y; ), запишем логарифм правдоподобияl( 0) (= l(y; 0)) в виде) = ln f(x |0) − ln h(x | y,0ln g(y |0)| ) − H( 0 | ).(3.9.5)0) = Q(0l( 0) = ln g(y |и возьмем условное ожидание E [ · | Ψ (X) = y] . Таким образом, l рассматривается как функция переменной 0 ∈ Θ:Здесь Q( 0 | ) (= Q(y; 0 | )) и H( 0 | ) (= H(y; 0 | )) означают условныеожиданияQ( 0 | ) = E [ln f(X | 0) | Ψ (X) = y] ,(3.9.6)H( 0 | ) = E [ln h(X | y; 0) | Ψ (X) = y] ,(N)):= argmax[Q( |(N+1)причем предполагается, что они существуют для всех пар 0 , ∈ Θ.Теперь определим N-ю итерацию EM-алгоритма как отображение A:(N)∈ Θ 7→ (N+1) = A( (N) ) ∈ Θ, действующее следующим образом.E-шаг.
При заданном (N) находим Q( | (N) ).M-шаг. Выбираем в качестве (N+1) значение, которое максимизируетфункцию ∈ Θ 7→ Q( | (N) ):(ср. с теоремой 3.7.8).Доказательство сходимости алгоритма Баума—Уэлча не использует егоспецифических черт. В частности, структура цепи Маркова используетсятолько при доказательстве монотонности в формулах (3.7.9) и (3.7.30).Это позволяет нам доказать сходимость более широкого класса алгоритмов, а именно т. н. алгоритмов типа EM (expectation-modification, т. е.условное математическое ожидание плюс модификация) и их обобщенийGEM (generalized EM, т.
е. обобщенные EM).Эпиграф к этому параграфу взят из статьи: Meng X.-L. , van Dyk D.The EM algorithm — an old folk-song sung to fast new time // Journ.Roy. Stat. Soc. 1997. V. B 59. P. 511–567. За ним следовал такой текст:«Так же как народную песню, как правило, напевают в течение многихлет, прежде чем ее мелодия становится хорошо узнаваемой, различныеметоды и идеи, относящиеся к EM... можно отыскать в [ранней] литературе... Аналогия с народной песней является точной и в том смысле, чтоона символизирует коллективный вклад в развитие EM-алгоритма...
Бауми другие (1970 г.), пожалуй, внесли наибольшее усовершенствование.»В 1992 г. число публикаций, относящихся к EM, превзошло 1000 (сейчасих, конечно, намного больше). Любопытно, что среди 300 журналов, гдебыли опубликованы вышеупомянутые 1000 статей, лидировал, что неудивительно, Журнал Американского общества статистиков (Journal ofthe American Statistical Association), но Журнал Королевского Статистического общества, Серия B (Journal of the Royal StatisticalSociety), был вытеснен на пятое место Журналом науки о сыроделиии производстве молочных продуктов (Journal of Dairy Science).EM-алгоритм работает следующим образом. У нас есть два выборочных пространства X (недоступные выборки с «полными данными») и Y(наблюдаемые выборки с «неполными данными»), а также многозначноеотображение Ψ из X в Y (это наш механизм наблюдения, он неточен).Это означает, что вместо наблюдения выборочного вектора x ∈ X мынаблюдаем искаженный вектор y = Ψ (x) ∈ Y .
Такая ситуация частовозникает в статистике, когда данные можно группировать или цензурировать, т. е. усекать или отбрасывать. Пусть f(x; ) (= f X (x; )), x ∈ X , —X (y)(3.9.2)b |XT = xT) > L(P | XT = xT) }b ∈ P : L(PD (P) = {Pb∗ =(ср. с теоремой 3.7.12). Аналогично для задачи интерполяции P= Π (P) ∈ D (P), где(3.9.1)b ∈ U : L( ; Z)b > L( ; Z) }D (Z) = {Z551плотность распределения случайного выборочного вектора X, зависящегоот параметра ∈ Θ. Тогда плотность распределения g(y; ) (= g Y (y; ))случайного вектора Y = Ψ (X) задается формулойZf(x; ) dx, y ∈ Y ,(3.9.3)g(y; ) =b ∗ = Φ (Z) лежит в областигде Z§ 3.9.
Обобщения алгоритма Баума—Уэлча. Глобальная сходимость итерацийГлава 3. Статистика цепей Маркова с дискретным временем550∈ Θ] .(3.9.7)1(i= 1) ln ϕ (yi − ).(3.9.12)= 0) ln ϕ (yi) + 1(nXln Lполн ( ) =ii=1Приведенное выше описание EM-алгоритма приводит в этом случаек такой формуле:сящая от и 0 . Очевидно, Q( 0 , ) является квадратичной формой по 0и может быть легко максимизирована.Далее, в силу приведенного выше описания M-шага значение (N+1) ,полученное после (N + 1)-й итерации алгоритма, задается формулой22(3.9.14)(N) )i=1Пример 3.9.2 (двумерные нормальные данные с пропущенными знаXчениями).
Пусть X = X1 — случайный вектор, имеющий двумерное2нормальное распределение N( , Σ), где — действительный двумерный1векторсредних i = EXi , а Σ — положительно определенная дей21112ствительная (2 × 2)-матрица ковариаций видас элементами2122= Var Xi и 21 = 12 = Cov(X1 , X2), i = 1, 2. Предположим, чтомы хотим оценить совокупность параметров = { i , ij , i, j = 1, 2} послучайной выборке размера n, образованной из независимых копий случайного вектора X, причем данные по первой компоненте X 1 пропущены вm1 местах, а данные по второй компоненте X2 пропущены в m2 местахи m1 + m2 6 n; кроме того, позиции, по которым пропущены X 1 , непересекаются с позициями, по которым пропущены X 2 . В этом примереΘ ⊂ R5 ; на самом деле Θ лежит в R × R × R+ × R+ × R.
(Первые двадекартовых множителя R задают область изменения для 1 и 2 , две копиипространства R+ выполняют ту же задачу для 11 и 22 , и последний множитель R — для 12 ; неравенство Коши—Шварца | 212 | 6 11 11 указываетiii=1.wi (которое является смесью с равными весами стандартной нормальной плотности распределения и нормальной плотности распределения со средним∈ R и единичной дисперсией. В данном примере полные данные xсостоят из пар x1 = (y1 , 1), . .
. , xn = (yn , n), где yj ∈ R и j == 0 или 1, j = 1, . . . , n; j определяет, какой плотностью распределения,ϕ (x)или ϕ (x − ), порождена j-я точка наблюдения yj . Функция Ψ (x) = y«уничтожает» все j и оставляет только точки y1 , . . . , yn . То же применимои к случайным выборкам Y и X.В этом примере логарифм правдоподобия ln Lнабл ( ) задается формулойn1X1ln ϕ (yi) + ϕ (yi − ) ,(3.9.11)ln Lнабл ( ) =)(3.9.10)(N)y ∈ R,, y) =yi w i (i=1nP1(ϕ (y) + ϕ (y − )),2(N)Y1Наблюдается выборка из н.о.р.с.в.
Y = ... , имеющая распределениеYng(y; ) =((3.9.9) с плотностью=(N+1)x ∈ R.(N+1)nP1 −x 2 / 2e,2ϕ (x) = √Теперь, чтобы сделать отображение Π взаимно однозначным преобразованием, нам нужно точно определить (N+1) среди точек максимума.Этот выбор может повлиять на различные аспекты реализации EM-алгоритма, как теоретические, так и практические. Очевидно также, что выборначального значения (0) должен быть разумным, критичным.Пример 3.9.1. В этом примере ϕ обозначает стандартную нормальнуюплотность распределения N(0, 1):ϕ (yi − )для i = 1, .
. . , n и c — постоянная, не завиϕ (yi) + ϕ (yi − )(3.9.8)).(N)) > Lнабл ((3.9.13)где wi ( ) =(N+1)Lнабл (i=1(1 − wi ( )) ln ϕ (yi) + wi ( ) ln ϕ (yi − 0) + c,nXQ( 0 , ) =553а логарифм правдоподобия ln Lполн ( ) — формулойЗначение (N+1) зависит от (N) и выборочного вектора y.Мы надеемся (и, вероятно, так и бывало в первых приложениях), чтопоследовательные значения (N) , N = 0, 1, .
. . , полученные с помощьюитераций алгоритма (их часто называют EM-последовательностью), будут «хорошими». В идеале (N) должны сходиться (и достаточно быстро)к ∗ , т. е. к о.м.п., максимизирующей правдоподобие g(y; ) из формулы(3.9.3). К сожалению, это не всегда так, и значительная часть материала этого параграфа направлена на то, чтобы прояснить эту ситуацию.С учетом дальнейших приложений используем в примерах альтернативное обозначение Lполн ( ) = Lполн (x; ) для функции правдоподобия f(x; )и Lнабл ( ) = Lнабл (y; ) для g(y; ), с l( 0) = Lнабл ( 0). При выполненииитераций значение Lнабл ( (N+1) ) всегда не меньше, чем Lнабл ( (N) ):§ 3.9.
Обобщения алгоритма Баума—Уэлча. Глобальная сходимость итерацийГлава 3. Статистика цепей Маркова с дискретным временем552...y1(m)y2(m)∗y2(m+1)...∗(m+m1)y2(m+m1 +1)y1∗...y1(N)∗Полные данные, конечно, представимы массивом векторов y (j)Ti =nX(j)yi , Sil =j=1.(3.9.15) (j) y1=(j) ,.y2(3.9.16)Мы используем для случайных выборок одновременно обозначение Y и Xи отождествляем Y с частью X указанным ранее способом.Логарифм правдоподобия ln Lнабл ( ) = ln Lнабл (y; ) для , вычисленный по наблюдаемым данным y, имеет видnXii−111nXj=m+m1 +1(j)(y1−1)2+mX+m1(j)(y2j=m+1− 2) .(3.9.20)Это замечание наводит на мысль, какой должна быть форма EM-алгоритма в данном примере. Предположим снова, что значение (N) == { i(N) , ij(N) , i, j = 1, 2} было получено при N-й итерации, и обозначимчерез E (N) математическое ожидание, взятое относительно двумерногонормального распределения, с параметрами (N) .
E-шаг (N + 1)-й итерации алгоритма следующий. Нужно вычислить условное математическоеожиданиеhi(N))=E(N)ln Lполн ( ) | yлогарифма правдоподобия, соответствующего полным данным ln L полн (X; ),при условии, что Y = y. Из равенства (3.9.17) видно, что указаннаяпроцедура сводится к вычислению математических ожиданий:2Eи()=(N)Логарифм правдоподобия, соответствующий полным данным, ln L(3.9.19)b i = Ti /n, b il = (Sil − n−1 Ti Tl) /n, i, l = 1, 2.(3.9.17)полн(j)−122(j)yi yl , i, l = 1, 2.Нетрудно видеть, что функция правдоподобия L полн ( ) принадлежитк экспоненциальному семейству с достаточной статистикой, образованнойсовокупностью {T1 , T2 , S11 , S12 , S22 }. Если бы в нашем распоряженииимелся массив полных данных x, то оценка методом максимального правдоподобия (по полным данным) b параметра имела бы видi=11−2j=1(y (j) − ) T Σ−1 (y (j) − ) −21X−mi ln22(3.9.18)j=1Q( ,12ln Lнабл ( ) = −n ln(2 ) − m ln(det Σ) −m1X1 2 12)] ,y1(n)(n)y22 12) +...−2−где22 111 22y1(1)(1)y2+− 2T1 (x=12 S1221 22j = 1, .
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
