Цепи Маркова (1121219), страница 88
Текст из файла (страница 88)
. . , n:+ n(−21 12)11 S22y1(1)y2(1)−+j=12 1122 S11n1 X (j)(y − ) T Σ−1 (y (j) − ) =2+ 2T2 (212) −−221= −n ln(2 ) − n ln( 1121− ( 11 22 − 212) −1 [2∗обозначают m2 наблюдений с пропущенной второй компонентой. Тогда yв развернутом виде соответствует массиву12ln Lполн ( ) = −n ln(2 ) − n ln(det Σ) −y2соответствуют точкам с полностью наблюдаемыми данными, где m = n −∗−m1 − m2 , б)(j) , j = m + 1, .
. . , m + m1 , обозначают m1 наблюденийy2 (j) с пропущенной первой компонентой и y1 , j = m + m1 + 1, . . . , n,555= ln Lполн (x; ), равенна то, что Θ — подпространство в R × R × R+ × R+ × R, не совпадающеес нулем или со всем пространством.)Обозначим массив наблюдаемых данных через y. Для определенности (j) y1пронумеруем данные таким образом, что а) y (j) =(j) , j = 1, .
. . , m,§ 3.9. Обобщения алгоритма Баума—Уэлча. Глобальная сходимость итерацийГлава 3. Статистика цепей Маркова с дискретным временем554E(N)Y1(j) | y(j)Y2 | yи Eи E(N)(N)(Y1(j) ) 2 | y ,j = m + 1, . . . , m + m1 ,(j)(Y2 ) 2 | y , j = m + m1 + 1, . . . , n.(j)y2 ,(j)Y1(j)y1 .=а во второй строке — на=Соответственнозаменить наможно использовать обозначение (j) (j) (j)(j) E (N) Y1 | y2 и E (N) (Y1 ) 2 | y2 , j = m + 1, . . . , m + m1 , (3.9.21)и∗22 (N)=(N)221−(j)y1−(N) 212 )(N) (N)11 22(+(N)12(N)11z (N) =(N)2(j)(Y2(j) ) 2 | Y1(j) = y1(j) = (z (j) (N)) 2 +Здесь(N)E∗22 (N).(N) 1при N-й итерации (ср. с (3.9.23)). Формулы— это значения, вычисленные (j)(j)(j) (j) (3.9.22) для E (N) Y2 | y1 и E (N) (Y2 ) 2 | y1 легко получить перестановкой индексов 1 и 2.M-шаг (N + 1)-й итерации осуществляется простой заменой статистик Ti и Sil на Ti(N) и Sil(N) соответственно, где последние формируются(j)путем подстановки в формулу (3.9.18) вместо пропущенных значений y i(j)и (yi ) 2 , i = 1, 2, их текущих условных ожиданий (3.9.21) и (3.9.22).Вектор-функция C(x) является достаточной статистикой.
Эта конструкциявключает в себя много популярных примеров: многомерное нормальноераспределение, пуассоновское, мультиномиальное, гипергеометрическое(см. том I, § 3.6).Предположим, что функция Ψ: x → y задана и нам доступна выборкаy = Ψ (x). Итерация EM-алгоритма с номером N + 1 осуществляется теперьследующим образом.E-шаг: запишем функциюQ(y;|(N))=E(N)(H(X) | Ψ (X) = y) +e − (grad B( )) T + B( ).+ (grad B( )) T C(y)Здесьe l=EC(y)(3.9.25)[C(X) l | Ψ (X) = y] .(N)(3.9.26)Заметим, что слагаемое в первой строке правой части равенства (3.9.25)не зависит от , а значит, не принимает участия в максимизации по .В отличие от него, третье и четвертое слагаемые зависят от , но неeeгде C(y)определенозависят от (N) .
Второе слагаемое (grad B( )) T C(y),(N).в формуле (3.9.26), зависит и от , и отM-шаг: при заданном значении параметра (N) наша цель — найти максимум функции Q(y; , (N) ) по (при фиксированном y). После отысканияэтого максимума значение, при котором он достигается, присваивается11иj] .Таким образом, приведенные выше математические ожидания (3.9.21)имеют вид(N) (j) (j) (N)(j)(N) y1 − 1E (N) Y2 | y1 = 2 + 12(N).(3.9.24)22B( ) [Cl (x) −j11∂j=1212nX∂1−(grad B( )) T [C(x) − ] =где1)22−=−112 11 (x1+2и дисперсией∗22f(x; ) = exp[(grad B( )) T [C(x) − ] + B( ) + H(x)] ,=nформулойимеет двумерное нормальное распределениеN( , Σ), то распределение Y2 при условии, что Y1 = y1 , является нормальным N( ∗2 , ∗22) со средним∗2(3.9.23)= ...
∈ Rn , задаетсяY1Y2семейство плотностей распределения f(x; ), i, j = 1, 2} задается1(3.9.22)(j)(j) (Y2 ) 2 | y1 , j = m + m1 + 1, . . . , n.Далее, если вектор(N)и E(j) (j)Y 2 | y1(N+1),ij= Sil(N) − n−1 Ti(N) Tl(N) /n.(N+1)il(N)(N+1),iПример 3.9.3 (параметрическое оценивание для экспоненциальныхсемейств). В этом примере Θ = Rn . Напомним, что экспоненциальное E= Ti(N) /n,= {557и его зеркальные отображения(N+1)i(N+1)(j)Y2Соответственно значениеформуламиВ силу независимости выборки, в первой строке условие Y = y можно§ 3.9. Обобщения алгоритма Баума—Уэлча.
Глобальная сходимость итерацийГлава 3. Статистика цепей Маркова с дискретным временем556Эта функция вогнутая и имеет единственный максимум, задающий о.м.п.:b∗ =b∗n1Xxj .ni=1можно также получить как предел различных аппроксимаций,Значениечто дает хорошую практику при реализации EM-алгоритма. Стандартныйметод аппроксимаций — это метод наискорейшего спуска (используемыйдля минимизации выпуклой квадратичной формы).
Но для оценки скоростисходимости нужно следующее неравенство.Пример 3.9.4 (неравенство Канторовича). Пусть Σ — положительно,(3.9.28)При подходящей замене координат матрица Σ становится диагональной, игде= x2ii.Pni=1i=12i xix2i2 Pni=11x2i /i.Pniii=1nP =ϕ( ),():=PnPni/ ii=1x2i . Функция y 7→ 1/y выпуклая при y > 0, точка ϕ ( ) лежитϕ( )>()inf166n(1 +1/yn − )/1n=(411+nn)2на кривой, а точка ( ) является линейной комбинацией точек на кривой.Поэтому минимальное значение частного достигается для некоторого == 1 1 + n n , 1 + n = 1. В этом случае 1 / 1 + 2 / n = ( 1 + n − ) / 1 n ,и мы получаем,так как минимум достигается в точке = ( 1 + n) /2.Пример 3.9.5 (метод наискорейшего спуска для квадратичной формы).Для заданной положительно определенной действительной (n×n) -матрицыΣ и векторов b, x0 ∈ Rn положимj=12n1X(xj − ) T Σ−1 (xj − ).2∈ Rn 7→ −В этом случае логарифм правдоподобия ln Lнабл ( ) = ln Lполн ( ) — этоотрицательная квадратичная форма от (с коэффициентами, зависящимиот выборки x):+)(3.9.27)++gk = Σxk − b и xk+1 = xk −gTk gkgk ,gTk Σgk∈R .−−=4(где − и + являются соответственно наименьшим и наибольшим собственными значениями матрицы Σ.Решение.
Пусть собственные значения i матрицы Σ удовлетворяютнеравенствам0 < − = 1 6 . . . 6 n = +.x∈R ,d>d||x||4(xT Σx) (xT Σ−1 x)i=1j=1xnn1X−(xj − ) T Σ−1 (xj − ) ,2expn1p(2 ) n det Σx1= ... ∈ Rn имеет место следующая оценка:ют многомерное нормальное распределение с известной ковариационной(d × d)-матрицей Σ и неизвестным вектором средних = ∈ Rd .
В этомслучае совместные плотности распределения f X ( · ; ) и gY ( · ; ) совпадаюти задаются выражениемопределенная действительная (n × n)-матрица. Для любого вектора x = X1Y1.. .. X =и неполные данные Y =совпадают и при этом они..XnYn (1) Xj .. образованы н.о.р. векторами Xj = . , j = 1, .
. . , n, которые име(d)Xj559переменной (N+1) (= (N+1) (y)). К сожалению, непосредственная максимизация, которую мы осуществили в примере 3.9.2, бывает возможна оченьредко.Как было сказано ранее, вопрос о сходимости величин (N) , полученныхв процессе итераций, является достаточно тонким и требует дополнительного анализа. Относительно простой случай — это когда полные данные § 3.9. Обобщения алгоритма Баума—Уэлча. Глобальная сходимость итерацийГлава 3. Статистика цепей Маркова с дискретным временем558k = 0, 1, . . .(3.9.29)Тогда для любого x0 ∈ Rn последовательность {xk } сходится при k → ∞к единственной точке минимума x∗ функции12f(x) = xT Σx − xT b.Глава 3. Статистика цепей Маркова с дискретным временем1=и единичной (2 × 2)-матрицей ковари-2аций. В этом примере неизвестный параметрy1.данные совпадают: x = y == , а полные и неполныелами(N)1(N)2}, где(N)== y1 + r (N) cos ϑ (N) ,(N)1(N)2, задается форму-= y2 + r (N) sin ϑ (N) ,где r (0) = 2, ϑ (0) = 0, а(N)|(N)).(3.9.30)0∈ Θ.(3.9.31)∈ Θ это неравенство задает множество∈ Θ : Q( 0 | ) > Q( | ) },(3.9.32)∈ Θ 7→ M( ),и мы вынуждены рассматривать многозначное отображение(N+1)причем гарантировано выполнение условия∈ M((N)),) − ln L((N))=1 (N) 2((r ) − (r (N+1) ) 2) =21= [(r (N) ) 2 − (2 − (r (N) ) −1) 2 ] ,2(3.9.33)поскольку r (N+1) = 2 − (r (N) ) −1 .
Теперь используем элементарную оценку0 < 2 − u−1 6 u для u > 1. Поскольку r (N) > 1 для каждого k, получаем,чтоln L( (N+1) ) − ln L( (N) > 0.Значит, последовательность (N) действительно является GEM-последовательностью.В силу того что r (N) → 1 при N → ∞ последовательность значенийфункции правдоподобия {L( (N) } сходится к величине(2 ) −1 e−1 .n = 0, 1, .
. .Последовательность (N) , имеющая указанное свойство, называетсяGEM-последовательностью.(N+1)Но для последовательности { (N) } любая точка единичного круга с центромв точке y является предельной.0M( ) = {ln L(,| ) > Q( | ),При заданномQ(0неравенство (3.9.30) приводит к свойству монотонности (3.9.8) (что является основополагающей чертой EM-и GEM-алгоритмов).При выполнении GEM-алгоритма мы записываем неравенство∈ Θ,0В терминах аналитической геометрии на плоскости r и ϑ являютсяполярными координатами с центром в наблюдаемом векторе y. Для логарифма правдоподобия ln L( (N) ) = ln L(y | (N) ) видим, что,(i + 1) −1 , k = 1, 2, .
. .i=1| ),0H( | ) > H(Поскольку математическое ожидание H( | ) удовлетворяет, в силу неравенства Гиббса, соотношениюNX) > Q((N)|(N+1)r (N) = 1 + (N + 1) −1 , ϑ (N) =Q(Довольно часто нереально выполнить численно процедуру максимизации на M-шаге. Мы уже столкнулись с этим в примере 3.9.3; в случаемарковских цепей чаще ситуация именно такова, а никак не наоборот. Наэтот случай предусмотрен алгоритм GEM — обобщенное условное ожидание плюс преобразование. При выполнении GEM-алгоритма мы простовыбираем (N+1) таким образом, что(N)GEM-последовательность {y2D(xk) − D(xk+1)(gkT gk) 2= T.D(xk)(gk Σgk) (gkT Σ−1)gkс неизвестным средним−где − и + , как и ранее, являются соответственно минимальным и максимальным собственными значениями матрицы Σ.Указание.
Примените неравенство Канторовича к соотношению+Пример 3.9.6. Приведем пример GEM-последовательности { (N) },для которой L( (N) ) монотонно сходится, в то время как сама последовательность { (N) } не сходится, а имеет в качестве множества предельныхточек единичный круг. Рассмотрим двумерное нормальное распределениеБолее того, если обозначить D(x) = (x − x∗) T Σ (x − x∗), то имеет место2следующая оценка: − 2+−D(xk+1) 6D(xk) ∀ k,+5611§ 3.9. Обобщения алгоритма Баума—Уэлча. Глобальная сходимость итераций560FΦ = {Z ∈ U : Φ (Z) = Z}.(3.9.35)В этом случае можно положитьF (Z) = L( ; Z), Z = ( , P, Q) ∈ U(3.9.36)для любой заданной обучающей последовательности ; см.