Цепи Маркова (1121219), страница 86
Текст из файла (страница 86)
е.sXpij =иqjk =Xsm=1sXm=1Xκm=1∂m∂L( ; Z)m −1∂pimL( ; Z)∂ pim∂qjmL( ; Z)∂ qjmj −1 −1∂L( ; Z),∂ j∂L( ; Z)∂ pij(3.8.24)∂L( ; Z).∂ qjk(3.8.25)pijqjk(3.8.23)Эти уравнения возникали и раньше, например в определении (3.7.5),а также в примере 3.7.3, и в общем контексте в (3.7.22) – (3.7.25) (см. также(3.7.32) и (3.7.38 б)).
Эту систему нелинейных уравнений в общем случаенельзя решить аналитически.Поскольку максимизацию можно осуществить по каждой переменнойотдельно, рассмотрим максимизацию только по переменным p ij .Лемма 3.8.4. Имеет место следующее уравнение:n−1X∂L( ; Z) =∂ pijm (i)qjm+1m+1 (j).m=1Д о к а з а т е л ь с т в о. ЗапишемL( ; Z) = P Z (b(X) = ; Z) =sXj=1n (j) n (j);(3.8.20)j> 0,b ∗ = (b∗ , Pb∗, Qb ∗),Φ (Z) = Z545где , i и lj — множители Лагранжа. Стационарная точка из внутренностиобласти удовлетворяет соотношениямФормулы (3.8.17) – (3.8.19) образуют базу современных вычислительных приемов, интенсивно используемых во многих приложениях; см.,например, статью L.R.
Rabiner. A tutorial on hidden Markov models and selected applications in speech recognition. Proc. IEEE, 77:2 (1989), 257 –288,а также [RJ] . Доступное введение в с.м.м. в биологии изложено в статье A.Krogh. An introduction to hidden Markov models for biological sequences. In:S.L, Salzberg, D.B. Searls and S. Kasif. Computational methods in molecular biology. Amsterdam: Elsevier, 1999, pp.
45–63; см. также [DRKM] .Еще раз повторим, что равенства (3.8.17) – (3.8.19) задают точку§ 3.8. Скрытые марковские модели, II. Обучающий алгоритм Баума—УэлчаГлава 3. Статистика цепей Маркова с дискретным временем544(3.8.26)m+1 (j)=Xsпоскольку n (j) = 1, этот множитель можно опустить. Далее, воспользуемся соотношением (3.8.12) при m = n − 1:m (i)piji=1qjm+1.547что в силу соотношения (3.8.14) равноn−2 (i)qjn−1n−1 (j).(3.8.30)Далее, объединяя слагаемые II и III и меняя порядок суммирования, получаемXs hsiX∂(l)pqpq(3.8.31)II + III =rm m n .n−2lr r n−1Отсюда следует, что§ 3.8. Скрытые марковские модели, II.
Обучающий алгоритм Баума—УэлчаГлава 3. Статистика цепей Маркова с дискретным временем546r,l=1n.(3.8.27)Рекуррентное соотношение (3.8.14) приводит нас к равенствуsXДалее, непосредственным дифференцированием получаемn−1 (i)qjn+s hX∂∂L( ; Z) =∂ pij∂ pijl,m=1i(l)plm qm n .n−1(3.8.28)n=n−1 (r).Поэтому правая часть соотношения (3.8.31) принимает видs hX∂∂ pijr,l=1i(l)pir qrn−2n−1n−1 (r).Итак, уравнение (3.8.28) принимает вид+ n−2 (i)qjs hX∂+∂ pijm=1r,l=1Подставляя частную производную ∂ n−1 (m) /∂ pij во второе слагаемоев правой части равенства (3.8.28), получимn (j)∂ pijn−1n−1 (j) +n−2 (l)inn−1 (i)qj∂L( ; Z) =∂ pijЗначит, нам нужно вычислить только двойную сумму в правой части.Имеемis h ∂P(m)pmj qj n−1 , если i = j, n−2 (i)qj n−1 +n−2∂m=1 ∂ pij(l)=hisP∂∂ pij n−1если i 6= j.n−2 (m) pmj qm n−1 ,prm qmm=1n−1 (i)pij qjm=1i,j=1∂ pijL( ; Z) =sXplr qrn−1n−1 (r).(3.8.32)+r=1sXr,m,l=11(m 6= i)pjr qrn+r=1h ∂∂ pijl=1n−2 (l)plm qm∂ pijn−1in−2 (l)plj qj n−1 pjr qr n +pmr qrn= I + II + III.(3.8.29)r=1n−2 (i)qjn−1pjr qrn=I=n−2 (i)qjn−1Xsr=1pjr qrnm (i)qjm+1,m+1 (j),m=0что и завершает доказательство.Лемма 3.8.4 приводит уже к другому виду уравнений (3.7.22) – (3.7.24),(3.8.4), (3.8.7) – (3.8.8) и (3.8.17) – (3.8.19).
Действительно, подставим соотношение (3.8.26) в (3.8.24) и поменяем порядок интегрирования в знаменателе. ПолучимXn−1sX(i)pq(j),mij j m+1 m+1sXТеперь найдем значение каждого из трех слагаемых в правой частиравенства (3.8.29). Во-первых,n−1X∂L( ; Z) =∂ pijn−1s XsX∂n−2 (i)qj==sXn∂ pijr,l=1i(r)prl qln− 1s hX∂Двойную сумму в правой части уравнения (3.8.32) можно подвергнуть аналогичной процедуре, подставляя производную ∂ n−2 (l) ∂ pij и преобразуяполученные суммы так же, как это делалось выше. Поскольку r (0) не содержит ни одной вероятности pij , процесс дифференцирования обрываетсяпри l = 0.
Итак, мы получимm=1j=1548Глава 3. Статистика цепей Маркова с дискретным временемpij qjm+1m+1 (j)sXи обратная рекурсия (3.8.14) приведет к формулам=m (i).j=1Значит, имеет место равенствоpijj=1m (i) m (i)n−1X∂L( ; Z) =∂ pijsXm=1что немедленно дает нам правую часть соотношения (3.8.17).Основываясь на уравнениях (3.8.23) – (3.8.25), перепишем преобразование Баума—Уэлча Φ для задачи фильтрации с.м.м. в видеb∗, Qb∗ ,Φ : ( , P, Q) 7→ b∗ , Pгдеb∗ =jиbp∗ij =bq∗jk =Xs∂mm=1sXpimm=1Xκm=1qjm∂L( ; Z)m −1∂L( ; Z)∂ pim∂L( ; Z)∂ qjm∂L( ; Z),∂ j(3.8.34)pij∂L( ; Z)∂ pij(3.8.35)qjk∂L( ; Z).∂ qjk(3.8.36) −1 −1(3.8.33)jВажно подчеркнуть, что в силу леммы 3.8.2, преобразование (3.8.33) имеетсущественное преимущество с вычислительной точки зрения: оно сводитсяк сумме суперпозиций локальных преобразований.
Это приводит к огромной экономии в вычислениях, требуя ns2 операций, в то время как прямойподход требует nsn операций.Анализ сходимости итераций ΦN отображения Φ опирается на основныегеометрические идеи, восходящие к первой половине XX в. Он проведенв § 3.9, а результат сформулирован в теореме 3.9.9 (для задачи фильтрациис.м.м.). Аналогичный результат имеет место для задачи интерполяции.Значительная часть данного параграфа служила иллюстрацией важности выбора математических методов для проведения компьютерных вычислений.
Завершим его историей о профессиональных вычислителях (т. н.«computors»).§ 3.9. Обобщения алгоритма Баума—Уэлча. Глобальная сходимость итераций549Недавно авторы этой книги натолкнулись на воспоминания А. А. Самарского, выдающегося русского прикладного математика и действительного члена Российской академиинаук, о раннем этапе развития параллельных вычислений в СССР.
(См. Губарев В. С. Белыйархипелаг Сталина. М.: Молодая гвардия, 2004.) В конце 1940-х гг. в Советском Союзеусиленными темпами создавалась собственная атомная бомба, что требовало проведенияогромного объема вычислений. В частности, ежедневно необходимо было численно решатьсистемы, состоящие из сотен линейных уравнений. Советская компьютерная промышленностьтого периода выпускала только механические арифмометры. Однако вычисления были выполнены быстро и надежно. Советское решение проблемы было весьма элегантным: Самарский,возглавлявший вычислительную группу, взял себе в подчинение около 30 юных девушек«компьютеров», только что окончивших Московский институт геодезии и картографии.
Каждая из девушек должна была решить на своем персональном арифмометре дюжину уравненийи передать свои результаты другой девушке для сравнения и дальнейшего использованияв соответствии со специально разработанным алгоритмом параллельных вычислений. Вовремя первого испытательного взрыва (август 1949 г.), ученые сумели численно предсказатьрезультаты испытания с точностью до 30 %, что, согласно Самарскому, превосходило уровень точности, достигнутый американцами (которые уже пользовались первыми прототипамиЭВМ).На эффективность советской вычислительной системы того времени могло повлиять то,что неспособность верно провести вычисления рассматривалась как акт саботажа и моглаиметь тяжелые последствия.
Блестящий физик или инженер, работавший с радиоактивнымиматериалами, с легкостью мог превратиться в заключенного и отправиться в шахту добыватьтот же радиоактивный материал, но уже без всякой защиты.В это время специально натренированные вычислители (иногда их называли computors,в отличие от вычислительных устройств computers) использовались во многих странах. Так,в некоторых британских университетах была специальная должность («computor»), закрепленная за профессорами математики.
Обязанностью ассистента было проводить вычисления,комбинируя наборы специализированных (механических или электрических) калькуляторов,подходящие для решения данной задачи.§ 3.9. Обобщения алгоритма Баума—Уэлча.Глобальная сходимость итерацийЯ чувствовал себя подобно старому менестрелю, который пел своюпесню в течение 18 лет, а теперь с огромным удовлетворениемобнаружил, что его фольклор стал темой могучей симфонии.Х. O. Хартли (1912–1980), американский статистикКак было установлено выше, преобразование Баума —Уэлча для задачи фильтрации с.м.м. задается с помощью (эквивалентных) формул(3.7.28), (3.8.20), (3.8.33), а для задачи интерполяции с помощью формулы(3.7.33). Для определенности далее будет рассматриваться только задачабез ограничений.
В случае задачи фильтрации преобразование переводитb∗,изначальную модель Z в модифицированную, или уточненную модель Zгде X (y) — это прообраз {x : Ψ (x) = y}. Параметр неизвестен и будетоцениваться методом максимального правдоподобия, т. е. путем максимизации функции g(y; ), или, что эквивалентно, ln g(y; ), на множестве∈ Θ. Поскольку выборка x недоступна наблюдениям, заменим логарифм правдоподобия ln f(x; ) его условным ожиданием при заданном y.В описании алгоритма это проделывается для произвольного ∈ Θ, но припоследовательных итерациях в качестве будем выбирать значение (N) ,полученное после N-й итерации; см. ниже.Равенство (3.9.3) — это отправная точка так называемого E-шагаEM-алгоритма.