Алгоритм Баума – Уэлча

В электротехнике , статистических вычислениях и биоинформатике алгоритм Баума-Уэлча представляет собой частный случай алгоритма ожидания-максимизации, используемого для поиска неизвестных параметров скрытой марковской модели (HMM). Он использует алгоритм вперед-назад для вычисления статистики для шага ожидания.

История

Алгоритм Баума-Уэлча был назван в честь его изобретателей Леонарда Э. Баума и Ллойда Р. Уэлча . Алгоритм и скрытые марковские модели были впервые описаны в серии статей Баума и его коллег в Центре коммуникационных исследований IDA в Принстоне в конце 1960-х - начале 1970-х годов. Одним из первых крупных применений HMM была область обработки речи . В 1980-х годах HMM стали полезным инструментом анализа биологических систем и информации, в частности генетической информации . С тех пор они стали важным инструментом вероятностного моделирования геномных последовательностей.

Описание

Скрытая модель Маркова описывает совместную вероятность набора « скрытых » и наблюдаемых дискретных случайных величин. Он основан на предположении, что i -я скрытая переменная с учетом ( i - 1)-й скрытой переменной не зависит от предыдущих скрытых переменных, а текущие переменные наблюдения зависят только от текущего скрытого состояния.

Алгоритм Баума-Уэлча использует хорошо известный алгоритм EM для нахождения оценки максимального правдоподобия параметров скрытой марковской модели с учетом набора наблюдаемых векторов признаков.

Позволять $X_{t}$ быть дискретной скрытой случайной величиной с $N$ возможные значения (т.е. мы предполагаем, что существуют $N$ штаты в целом). Мы предполагаем, $P(X_{t}\mid X_{t-1})$ не зависит от времени $t$, что приводит к определению независимой от времени матрицы стохастических переходов

$A=\{a_{ij}\}=P(X_{t}=j\mid X_{t-1}=i).$

Распределение начального состояния (т.е. когда $t=1$) определяется

$\pi _{i}=P(X_{1}=i).$

Переменные наблюдения $Y_{t}$ могу взять один из $K$ возможные значения. Мы также предполагаем, что наблюдение в «скрытом» состоянии не зависит от времени. Вероятность определенного наблюдения $y_{i}$ во время $t$ для государства $X_{t}=j$ дается

$b_{j}(y_{i})=P(Y_{t}=y_{i}\mid X_{t}=j).$

Учитывая все возможные значения $Y_{t}$ и $X_{t}$, мы получаем $N\times K$ матрица $B=\{b_{j}(y_{i})\}$ где $b_{j}$ принадлежит всем возможным состояниям и $y_{i}$ принадлежит всем наблюдениям.

Последовательность наблюдений определяется выражением $Y=(Y_{1}=y_{1},Y_{2}=y_{2},\ldots ,Y_{T}=y_{T})$.

Таким образом, мы можем описать скрытую цепь Маркова формулой $\theta =(A,B,\pi )$. Алгоритм Баума – Уэлча находит локальный максимум для $\theta ^{*}=\operatorname {arg\,max} _{\theta }P(Y\mid \theta )$ (т.е. параметры HMM $\theta$ которые максимизируют вероятность наблюдения).

Алгоритм

Набор $\theta =(A,B,\pi )$ со случайными начальными условиями. Их также можно установить, используя предварительную информацию о параметрах, если она доступна; это может ускорить алгоритм, а также направить его к желаемому локальному максимуму.

Процедура пересылки

Позволять $\alpha _{i}(t)=P(Y_{1}=y_{1},\ldots ,Y_{t}=y_{t},X_{t}=i\mid \theta )$, вероятность увидеть наблюдения $y_{1},y_{2},\ldots ,y_{t}$ и находиться в состоянии $i$ во время $t$. Это находится рекурсивно:

1. $\alpha _{i}(1)=\pi _{i}b_{i}(y_{1}),$
2. $\alpha _{i}(t+1)=b_{i}(y_{t+1})\sum _{j=1}^{N}\alpha _{j}(t)a_{ji}.$

Поскольку этот ряд сходится экспоненциально к нулю, алгоритм будет численно неэффективным для более длинных последовательностей. Однако этого можно избежать с помощью слегка модифицированного алгоритма, масштабируя $\alpha$ в переднем и $\beta$ в обратной процедуре ниже.

Обратная процедура

Позволять $\beta _{i}(t)=P(Y_{t+1}=y_{t+1},\ldots ,Y_{T}=y_{T}\mid X_{t}=i,\theta )$ это вероятность окончания частичной последовательности $y_{t+1},\ldots ,y_{T}$ данное начальное состояние $i$ во время $t$. Мы рассчитываем $\beta _{i}(t)$ как,

1. $\beta _{i}(T)=1,$
2. $\beta _{i}(t)=\sum _{j=1}^{N}\beta _{j}(t+1)a_{ij}b_{j}(y_{t+1}).$

Обновлять

Теперь мы можем вычислить временные переменные согласно теореме Байеса:

$\gamma _{i}(t)=P(X_{t}=i\mid Y,\theta )={\frac {P(X_{t}=i,Y\mid \theta )}{P(Y\mid \theta )}}={\frac {\alpha _{i}(t)\beta _{i}(t)}{\sum _{j=1}^{N}\alpha _{j}(t)\beta _{j}(t)}},$

какова вероятность находиться в состоянии $i$ во время $t$ учитывая наблюдаемую последовательность $Y$ и параметры $\theta$

$\xi _{ij}(t)=P(X_{t}=i,X_{t+1}=j\mid Y,\theta )={\frac {P(X_{t}=i,X_{t+1}=j,Y\mid \theta )}{P(Y\mid \theta )}}={\frac {\alpha _{i}(t)a_{ij}\beta _{j}(t+1)b_{j}(y_{t+1})}{\sum _{k=1}^{N}\sum _{w=1}^{N}\alpha _{k}(t)a_{kw}\beta _{w}(t+1)b_{w}(y_{t+1})}},$

какова вероятность находиться в состоянии $i$ и $j$ время от времени $t$ и $t+1$ соответственно, учитывая наблюдаемую последовательность $Y$ и параметры $\theta$.

Знаменатели $\gamma _{i}(t)$ и $\xi _{ij}(t)$ одинаковы; они представляют вероятность наблюдения $Y$ учитывая параметры $\theta$.

Параметры скрытой марковской модели $\theta$ теперь можно обновить:

• $\pi _{i}^{*}=\gamma _{i}(1),$

что является ожидаемой частотой пребывания в состоянии $i$ во время $1$.

• $a_{ij}^{*}={\frac {\sum _{t=1}^{T-1}\xi _{ij}(t)}{\sum _{t=1}^{T-1}\gamma _{i}(t)}},$

которое представляет собой ожидаемое количество переходов из состояния i в состояние j по сравнению с ожидаемым общим количеством переходов из состояния i . Чтобы уточнить, количество переходов из состояния i означает переходы не в другое состояние j , а в любое состояние, включая его самого. Это эквивалентно тому, сколько раз состояние i наблюдается в последовательности от t = 1 до t = T - 1.

• $b_{i}^{*}(v_{k})={\frac {\sum _{t=1}^{T}1_{y_{t}=v_{k}}\gamma _{i}(t)}{\sum _{t=1}^{T}\gamma _{i}(t)}},$

где

$1_{y_{t}=v_{k}}={\begin{cases}1&{\text{if }}y_{t}=v_{k},\\0&{\text{otherwise}}\end{cases}}$

является индикаторной функцией, а $b_{i}^{*}(v_{k})$ - ожидаемое количество раз, когда выходные наблюдения были равны $v_{k}$ находясь в состоянии $i$ больше ожидаемого общего количества раз в состоянии $i$.

Эти шаги теперь повторяются итеративно до достижения желаемого уровня сходимости.

Примечание. Определенный набор данных можно переопределить. То есть, $P(Y\mid \theta _{\text{final}})>P(Y\mid \theta _{\text{true}})$. Алгоритм также не гарантирует глобального максимума.

Несколько последовательностей

Описанный до сих пор алгоритм предполагает единственную наблюдаемую последовательность $Y=y_{1},\ldots ,y_{N}$. Однако во многих ситуациях наблюдается несколько последовательностей: $Y_{1},\ldots ,Y_{R}$. В этом случае информация всех наблюдаемых последовательностей должна использоваться при обновлении параметров. $A$, $\pi$, и $b$. Предполагая, что вы вычислили $\gamma _{ir}(t)$ и $\xi _{ijr}(t)$ для каждой последовательности $y_{1,r},\ldots ,y_{N_{r},r}$, теперь параметры можно обновить:

• $\pi _{i}^{*}={\frac {\sum _{r=1}^{R}\gamma _{ir}(1)}{R}}$
• $a_{ij}^{*}={\frac {\sum _{r=1}^{R}\sum _{t=1}^{T-1}\xi _{ijr}(t)}{\sum _{r=1}^{R}\sum _{t=1}^{T-1}\gamma _{ir}(t)}},$
• $b_{i}^{*}(v_{k})={\frac {\sum _{r=1}^{R}\sum _{t=1}^{T}1_{y_{tr}=v_{k}}\gamma _{ir}(t)}{\sum _{r=1}^{R}\sum _{t=1}^{T}\gamma _{ir}(t)}},$

где

$1_{y_{tr}=v_{k}}={\begin{cases}1&{\text{if }}y_{t,r}=v_{k},\\0&{\text{otherwise}}\end{cases}}$

является индикаторной функцией

Пример

Предположим, у нас есть курица, от которой мы собираем яйца каждый день в полдень. Теперь, отложила ли курица яйца для сбора, зависит от некоторых неизвестных факторов, которые скрыты. Однако мы можем (для простоты) предположить, что курица всегда находится в одном из двух состояний, влияющих на то, несет ли она яйца, и что это состояние зависит только от состояния в предыдущий день. Теперь мы не знаем состояние в начальной отправной точке, мы не знаем вероятности перехода между двумя состояниями и не знаем вероятность того, что курица отложит яйцо в определенном состоянии. Для начала мы сначала угадаем матрицы перехода и выбросов.

Переход Штат 1 Штат 2 Штат 1 0.5 0.5 Штат 2 0.3 0.7

Эмиссия Нет яиц Яйца Штат 1 0.3 0.7 Штат 2 0.8 0.2

Исходный Штат 1 0.2 Штат 2 0.8

Затем мы берем набор наблюдений (E = яйца, N = нет яиц): N, N, N, N, N, E, E, N, N, N.

Это дает нам набор наблюдаемых переходов между днями: NN, NN, NN, NN, NE, EE, EN, NN, NN.

Следующим шагом является оценка новой матрицы перехода. Например, вероятность последовательности NN и состояния ⁠ $S_{1}$⁠ тогда ⁠ $S_{2}$⁠ определяется следующим образом: $P(S_{1})\cdot P(N|S_{1})\cdot P(S_{1}\rightarrow S_{2})\cdot P(N|S_{2}).$

Наблюдаемая последовательность	Наивысшая вероятность наблюдения этой последовательности, если состояние ⁠ $S_{1}$⁠ тогда ⁠ $S_{2}$⁠	Наивысшая вероятность наблюдения этой последовательности
НН	0.024 = 0.2 × 0.3 × 0.5 × 0.8	0.3584	⁠ $S_{2}$⁠ , ⁠ $S_{2}$⁠
НН	0.024 = 0.2 × 0.3 × 0.5 × 0.8	0.3584	⁠ $S_{2}$⁠ , ⁠ $S_{2}$⁠
НН	0.024 = 0.2 × 0.3 × 0.5 × 0.8	0.3584	⁠ $S_{2}$⁠ , ⁠ $S_{2}$⁠
НН	0.024 = 0.2 × 0.3 × 0.5 × 0.8	0.3584	⁠ $S_{2}$⁠ , ⁠ $S_{2}$⁠
NE	0.006 = 0.2 × 0.3 × 0.5 × 0.2	0.1344	⁠ $S_{2}$⁠ , ⁠ $S_{1}$⁠
ЭЭ	0.014 = 0.2 × 0.7 × 0.5 × 0.2	0.0490	⁠ $S_{1}$⁠ , ⁠ $S_{1}$⁠
В	0.056 = 0.2 × 0.7 × 0.5 × 0.8	0.0896	⁠ $S_{2}$⁠ , ⁠ $S_{2}$⁠
НН	0.024 = 0.2 × 0.3 × 0.5 × 0.8	0.3584	⁠ $S_{2}$⁠ , ⁠ $S_{2}$⁠
НН	0.024 = 0.2 × 0.3 × 0.5 × 0.8	0.3584	⁠ $S_{2}$⁠ , ⁠ $S_{2}$⁠
Общий	0.22	2.4234

Таким образом, новая оценка для ⁠ $S_{1}$от ⁠ до ⁠ $S_{2}$⁠ переход сейчас ${\frac {0.22}{2.4234}}=0.0908$ (в следующих таблицах они называются «псевдовероятностями»). Затем мы вычисляем ⁠ $S_{2}$от ⁠ до ⁠ $S_{1}$⁠ , ⁠ $S_{2}$от ⁠ до ⁠ $S_{2}$⁠ и ⁠ $S_{1}$от ⁠ до ⁠ $S_{1}$⁠ вероятности перехода и нормализуем их так, чтобы их сумма составляла 1. Это дает нам обновленную матрицу перехода:

Старая матрица перехода Штат 1 Штат 2 Штат 1 0.5 0.5 Штат 2 0.3 0.7

Новая матрица перехода (псевдовероятности) Штат 1 Штат 2 Штат 1 0.0598 0.0908 Штат 2 0.2179 0.9705

Новая матрица перехода (после нормализации) Штат 1 Штат 2 Штат 1 0.3973 0.6027 Штат 2 0.1833 0.8167

Далее мы хотим оценить новую матрицу выбросов,

Наблюдаемая последовательность	Наивысшая вероятность наблюдения этой последовательности если предполагается, что E происходит от ⁠ $S_{1}$⁠		Наивысшая вероятность наблюдения этой последовательности
NE	0.1344	⁠ $S_{2}$⁠ , ⁠ $S_{1}$⁠	0.1344	⁠ $S_{2}$⁠ , ⁠ $S_{1}$⁠
ЭЭ	0.0490	⁠ $S_{1}$⁠ , ⁠ $S_{1}$⁠	0.0490	⁠ $S_{1}$⁠ , ⁠ $S_{1}$⁠
В	0.0560	⁠ $S_{1}$⁠ , ⁠ $S_{2}$⁠	0.0896	⁠ $S_{2}$⁠ , ⁠ $S_{2}$⁠
Общий	0.2394		0.2730

Новая оценка E, полученная от ⁠ $S_{1}$⁠ выброс сейчас ${\frac {0.2394}{0.2730}}=0.8769$.

Это позволяет нам рассчитать матрицу выбросов, как описано выше в алгоритме, путем сложения вероятностей соответствующих наблюдаемых последовательностей. Затем мы повторяем, если N произошло от ⁠ $S_{1}$⁠ и если бы N и E произошли от ⁠ $S_{2}$⁠ и нормализовать.

Старая матрица выбросов Нет яиц Яйца Штат 1 0.3 0.7 Штат 2 0.8 0.2

Новая матрица выбросов (оценки) Нет яиц Яйца Штат 1 0.0404 0.8769 Штат 2 1.0000 0.7385

Новая матрица выбросов (после нормализации) Нет яиц Яйца Штат 1 0.0441 0.9559 Штат 2 0.5752 0.4248

Чтобы оценить начальные вероятности, мы предполагаем, что все последовательности начинаются со скрытого состояния ⁠ $S_{1}$⁠ и вычислите наибольшую вероятность, а затем повторите для ⁠ $S_{2}$⁠ . Затем мы снова нормализуем, чтобы получить обновленный исходный вектор.

Наконец, мы повторяем эти шаги до тех пор, пока полученные вероятности не сойдутся удовлетворительно.

Приложения

Распознавание речи

Скрытые марковские модели были впервые применены для распознавания речи Джеймсом К. Бейкером в 1975 году. Распознавание непрерывной речи происходит с помощью следующих шагов, смоделированных HMM. Анализ характеристик сначала проводится на временных и/или спектральных характеристиках речевого сигнала. Это создает вектор наблюдения. Затем этот признак сравнивается со всеми последовательностями блоков распознавания речи. Этими единицами могут быть фонемы , слоги или целые слова. Система декодирования лексикона применяется для ограничения исследуемых путей, поэтому исследуются только слова в лексиконе системы (словаре слов). Подобно декодированию словаря, системный путь дополнительно ограничен правилами грамматики и синтаксиса. Наконец, применяется семантический анализ, и система выводит распознанное высказывание. Ограничением многих приложений HMM для распознавания речи является то, что текущее состояние зависит только от состояния на предыдущем временном шаге, что нереально для речи, поскольку зависимости часто имеют длительность в несколько временных шагов. Алгоритм Баума-Уэлча также имеет обширные применения при решении HMM, используемых в области синтеза речи.

Криптоанализ

Алгоритм Баума-Уэлча часто используется для оценки параметров HMM при расшифровке скрытой или зашумленной информации и, следовательно, часто используется в криптоанализе . В сфере безопасности данных наблюдатель хотел бы извлечь информацию из потока данных, не зная всех параметров передачи. Это может включать в себя реверс-инжиниринг кодера канала . HMM и, как следствие, алгоритм Баума-Уэлча также использовались для идентификации произнесенных фраз в зашифрованных вызовах VoIP. Кроме того, криптоанализ HMM является важным инструментом для автоматического исследования данных синхронизации кэша. Это позволяет автоматически обнаруживать критическое состояние алгоритма, например значения ключей.

Приложения в биоинформатике

Поиск генов

Прокариотический

( Программное обеспечение GLIMMER Gene Locator and Interpolated Markov ModelER) было ранней программой поиска генов , используемой для идентификации кодирующих областей в ДНК прокариот . GLIMMER использует интерполированные марковские модели (IMM) для идентификации кодирующих областей и отличия их от некодирующей ДНК . Было показано, что последняя версия (GLIMMER3) обладает повышенной специфичностью и точностью по сравнению со своими предшественниками в отношении прогнозирования сайтов инициации трансляции, демонстрируя среднюю точность 99% в определении местоположения 3'-локаций по сравнению с подтвержденными генами у прокариот.

Эукариотический

Веб -сервер GENSCAN представляет собой локатор генов, способный анализировать эукариотические последовательности длиной до одного миллиона пар оснований (1 Мбит). GENSCAN использует общую неоднородную трехпериодическую марковскую модель кодирующих областей ДНК пятого порядка. Кроме того, эта модель учитывает различия в плотности и структуре генов (например, в длине интронов), которые встречаются в разных изохорах . В то время как большинство интегрированных программ для поиска генов (на момент выпуска GENSCAN) предполагали, что входные последовательности содержат ровно один ген, GENSCAN решает общий случай, когда присутствуют частичные, полные или множественные гены (или даже отсутствие гена вообще). Было показано, что GENSCAN точно предсказывает расположение экзона с точностью 90% и специфичностью 80% по сравнению с аннотированной базой данных.

Обнаружение изменений количества копий

Вариации числа копий (CNV) представляют собой распространенную форму вариаций структуры генома у человека. Двумерный HMM с дискретными значениями (dbHMM) использовался для присвоения хромосомным областям семи различных состояний: незатронутых областей, делеций, дупликаций и четырех переходных состояний. Решение этой модели с использованием Баума-Велча продемонстрировало способность предсказывать местоположение точки разрыва CNV примерно до 300 п.н. на основе экспериментов на микрочипах . Такая величина разрешения обеспечивает более точные корреляции между различными CNV и популяциями , чем это было возможно ранее, что позволяет изучать частоты популяций CNV. Он также продемонстрировал шаблон прямого наследования для конкретного CNV .

Реализации

• Accord.NET на C#
• ghmm C с привязками Python , поддерживающими как дискретную, так и непрерывную эмиссию.
• Библиотека Jajapy Python , которая реализует Баума-Уэлча на различных типах марковских моделей ( HMM , MC , MDP , CTMC ).
• HiddenMarkovModels.jl Пакет для Юлии .
• Функция HMMFit в RHmm для R. пакете
• хммпоезд в MATLAB
• ржавчина био в Rust

См. также

• Алгоритм Витерби
• Скрытая модель Маркова
• ЭМ-алгоритм
• Максимальная вероятность
• Распознавание речи
• Биоинформатика
• Криптоанализ

Внешние ссылки

• Комплексный обзор методов и программного обеспечения HMM в биоинформатике - Профиль скрытых марковских моделей
• Ранние публикации HMM Баума:
  • Метод максимизации, используемый при статистическом анализе вероятностных функций цепей Маркова
  • Неравенство с приложениями к статистическому оцениванию вероятностных функций марковских процессов и к модели экологии
  • Статистический вывод для вероятностных функций цепей Маркова с конечным состоянием
• Лекция Шеннона Уэлча, в которой рассказывается о том, как можно эффективно реализовать алгоритм:
  • Скрытые марковские модели и алгоритм Баума-Уэлча , Информационный бюллетень Общества теории информации IEEE, декабрь 2003 г.
• Альтернатива алгоритму Баума – Уэлча, алгоритм подсчета путей Витерби:
  • Дэвис, Ричард И.А.; Ловелл, Брайан С.; «Сравнение и оценка алгоритмов обучения ансамбля HMM с использованием критериев обучения, тестирования и числа условий» , Pattern Analysis and Applications, vol. 6, нет. 4, стр. 327–336, 2003.
• Интерактивная электронная таблица для обучения алгоритму прямого и обратного движения (таблица и статья с пошаговым руководством)
• Формальный вывод алгоритма Баума – Уэлча. Архивировано 28 февраля 2012 г. на Wayback Machine.
• Реализация алгоритма Баума – Уэлча.