Алгоритм Баума-Велша

Алгоритм Баума-Велша (англ. Baum–Welch algorithm) — алгоритм для нахождения неизвестных параметров скрытой Марковской модели. Использует алгоритм прямого-обратного хода.

История

Скрытые Марковские модели (HMMs) и алгоритм Баума-Велша впервые были описаны в заметках Леонарда Баума и его сверстников в конце [math]1960[/math]. Одно из первых основных приложений на основе HMMs было использовано в области обработки речи. В [math]1980[/math] HMMs стали эффективным инструментом в анализе биологических систем и информации, особенно в генном анализе.

Описание алгоритма

Пусть [math]Q_t[/math] — это дискретная случайная переменная, принимающая одно из [math]N[/math] значений [math](1..N)[/math]. Будем полагать, что данная модель Маркова, определенная как [math]P(Q_t \mid Q_{t - 1})[/math] однородна по времени, то есть независима от [math]t[/math]. Тогда можно задать [math]P(Q_t \mid Q_{t - 1}) [/math] как независящую от времени стохастическую матрицу перемещений [math]A = \{a_{ij}\} = p(Q_t = j \mid Q_{t - 1} = i)[/math]. Особый случай для времени [math]t = 1[/math] определяется начальным распределением [math]\pi_i = P(Q_1 = i)[/math].

Будем считать, что мы в состоянии [math]j[/math] в момент времени [math]t[/math], если [math]Q_t = j[/math]. Последовательность заданных состояний определяется как [math]q = \{q_1 \dots q_T \}[/math], где [math]q_t \in \{ 1\ldots N\}[/math] является состоянием в момент времени [math]t[/math].

Наблюдение может иметь одно из [math]L[/math] возможных значений, [math]Q_t \in \{o_1 \dots o_L\}[/math]. Вероятность заданного вектора наблюдений в момент времени [math]t[/math] для состояния [math]j[/math] определяется как [math]b_j(o_t) = P(O_t = o_t \mid Q_t = j)\ ( B = \{ b_{ij}\}[/math] — это матрица [math]L[/math] на [math]N)[/math]. Заданная последовательность наблюдений [math]O[/math] выражается как [math] O = (O_1 = o_1, \dots ,O_T = o_T)[/math].

Следовательно, мы можем описать скрытую модель Маркова с помощью [math] \lambda = (A, B, \pi)[/math]. При заданном векторе наблюдений [math]O[/math] алгоритм Баума-Велша находит [math] \lambda^*=\max_\lambda P(O\mid\lambda)[/math]. [math]\lambda[/math] максимизирует вероятность наблюдений [math]O[/math].

Исходные данные: [math] \lambda = (A, B, \pi)[/math] со случайными начальными условиями. Алгоритм итеративно обновляет параметр [math]\lambda[/math] до схождения в одной точке.

Прямая процедура

[math]a_i(t) = p(O_1 = o_1 \dots O_t = o_t, Q_t = \lambda_i )[/math], что является вероятностью получения заданной последовательности [math]\{ o_1 \dots o_t \}[/math] для состояния [math]i[/math] в момент времени [math]t[/math].

[math]a_i(t)[/math] можно вычислить рекурсивно:

[math]1.\,[/math] [math]a_i(1) = \pi_i \cdot b_i(O_1) [/math];

[math]2.\,[/math] [math]a_j(t + 1) = b_j(O_{t + 1})\displaystyle\sum^N_{i=1}a_i(t) \cdot a_{ij}[/math].

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

Данная процедура позволяет вычислить вероятность конечной заданной последовательности [math]\{ o_{t + 1} \dots o_T \}[/math] при условии, что мы начали из исходного состояния [math]i[/math], в момент времени [math]t[/math].

[math]\beta_i(t)[/math] можно вычислить рекурсивно:

[math]1.\,[/math] [math]\beta_i(T) = 1[/math];

[math]2.\,[/math] [math]\beta_i(t) = \displaystyle\sum^N_{j = 1}\beta_j(t + 1)a_{ij}b_j(O_{t + 1})[/math].

Обновление переменных

Определим временные переменные:

[math]\gamma_i(t) = p(Q_t=i\mid O,\;\lambda)=[/math] [math]\dfrac{\alpha_i(t)\beta_i(t)}{\displaystyle\sum^N_{j=1}\alpha_j(t)\beta_j(t)}[/math]

[math]\xi_{ij}(t) = p(Q_t=i,\;Q_{t+1}=j\mid O,\;\lambda)=[/math] [math]\dfrac{\alpha_i(t)a_{ij}\beta_j(t+1)b_j(o_{t+1})}{\displaystyle\sum^N_{i=1}\displaystyle\sum^N_{j=1}\alpha_i(t)a_{ij}\beta_j(t+1)b_j(O_{t+1})}[/math].

Имея [math]\gamma[/math] и [math]\xi[/math], можно определить:

[math]\bar\pi_i=\gamma_i(1)[/math],

[math]\bar{a}_{ij}=\dfrac{\displaystyle\sum^{T-1}_{t=1}\xi_{ij}(t)}{\displaystyle\sum^{T-1}_{t=1}\gamma_i(t)}[/math],

[math]\bar{b}_i(k)=\dfrac{\displaystyle\sum^T_{t=1}\delta_{O_t,\;o_k}\gamma_i(t)}{\displaystyle\sum^T_{t=1}\gamma_i(t)}[/math].

Используя новые переменные [math] A, B, \pi[/math] итерации продолжаются до схождения.

Пример

Предположим, у нас есть курица, с которой мы собираем яйца. Снесла ли курица яйца — зависит от некоторых неизвестных факторов. Для простоты предположим, что существуют лишь два состояния, которые определяют есть ли яйца. В начальный момент нам неизвестно текущее состояние, также нам неизвестна вероятность перехода из одного состояния в другое. Для начала возьмем произвольные матрицы переходов и состояний.

Переходы Состояние 1 Состояние 2 Состояние 1 [math]0.5[/math] [math]0.5[/math] Состояние 2 [math]0.3[/math] [math]0.7[/math]

Состояния Яйца не отложены Яйца отложены Состояние 1 [math]0.3[/math] [math]0.7[/math] Состояние 2 [math]0.8[/math] [math]0.2[/math]

Начальное состояние Состояние 1 [math]0.2[/math] Состояние 2 [math]0.8[/math]

Рассмотрим набор наблюдений ([math]E[/math] — яйца отложены, [math]N[/math] — яйца не отложены): [math]NN, NN, NN, NN, NE, EE, EN, NN, NN[/math].

Следующим шагом оценим новую матрицу переходов:

Последовательность	Вероятность последовательности и состояний	Наибольшая вероятность наблюдения
[math]NN[/math]	[math]0.024[/math]	[math]0.3584\,[/math] [math]S_2,[/math] [math]S_2[/math]
[math]NN[/math]	[math]0.024[/math]	[math]0.3584\,[/math] [math]S_2,[/math] [math]S_2[/math]
[math]NN[/math]	[math]0.024[/math]	[math]0.3584\,[/math] [math]S_2,[/math] [math]S_2[/math]
[math]NN[/math]	[math]0.024[/math]	[math]0.3584\,[/math] [math]S_2,[/math] [math]S_2[/math]
[math]NE[/math]	[math]0.006[/math]	[math]0.1344\,[/math] [math]S_2,[/math] [math]S_1[/math]
[math]EE[/math]	[math]0.014[/math]	[math]0.0490\,[/math] [math]S_1,[/math] [math]S_1[/math]
[math]EN[/math]	[math]0.056[/math]	[math]0.0896\,[/math] [math]S_2,[/math] [math]S_2[/math]
[math]NN[/math]	[math]0.024[/math]	[math]0.3584\,[/math] [math]S_2,[/math] [math]S_2[/math]
[math]NN[/math]	[math]0.024[/math]	[math]0.3584\,[/math] [math]S_2,[/math] [math]S_2[/math]
Итог	[math]0.22[/math]	[math]2.4234[/math]

Таким образом получаем новую оценку перехода из [math]S_1[/math] в [math]S_2[/math], которая составляет [math]\dfrac{0.22}{2.4234}[/math][math] = 0.0908[/math]. После этого можно подсчитать вероятность переходов из [math]S_2[/math] в [math]S_1[/math], [math]S_2[/math] в [math]S_2[/math], [math]S_1[/math] в [math]S_1[/math] и изменим их так, чтобы в суммы вероятностей давали [math]1[/math]. В итоге получаем новую матрицу переходов:

Старая матрица Состояние 1 Состояние 2 Состояние 1 [math]0.5[/math] [math] 0.5[/math] Состояние 2 [math]0.3 [/math] [math]0.7[/math]

Новая матрица (Псевдовероятности) Состояние 1 Состояние 2 Состояние 1 [math]0.0598[/math] [math]0.0908[/math] Состояние 2 [math]0.2179[/math] [math] 0.9705[/math]

Новая матрица (После изменения) Состояние 1 Состояние 2 Состояние 1 [math]0.3973[/math] [math] 0.6027[/math] Состояние 2 [math]0.1833[/math] [math] 0.8167[/math]

Далее оценим новую матрицу состояний:

Последовательности	Наибольшая вероятность наблюдения Если допустимо, что E получено из [math]S1[/math]	Наибольшая вероятность наблюдения
[math]NE[/math]	[math]0.1344\,[/math] [math] S_2, [/math] [math]S_1[/math]	[math]0.1344\,[/math] [math]S_2,[/math][math] S_1[/math]
[math]EE[/math]	[math]0.0490\,[/math] [math] S_1, [/math] [math]S_1[/math]	[math]0.0490\,[/math] [math]S_1,[/math][math] S_1[/math]
[math]EN[/math]	[math]0.0560\,[/math] [math] S_1, [/math] [math]S_2[/math]	[math]0.0896\,[/math] [math]S_1,[/math] [math]S_2[/math]
Итог	[math]0.2394[/math]	[math]0.2730[/math]

Новая оценка для [math]E[/math], полученная из [math]S_1[/math], составляет [math]\dfrac{0.2394}{0.2730}[/math] [math] = 0.8769[/math].

Благодаря этому, возможно рассчитать матрицу состояний:

Старая матрица Яйца не отложены Яйца отложены Состояние 1 [math]0.3[/math] [math]0.7[/math] Состояние 2 [math]0.8[/math] [math]0.2[/math]

Новая матрица (Оценка) Яйца не отложены Яйца отложены Состояние 1 [math]0.0876[/math] [math]0.8769[/math] Состояние 2 [math]1.0000[/math] [math] 0.7385[/math]

Новая матрица (После изменения) Яйца не отложены Яйца отложены Состояние 1 [math]0.0908[/math] [math]0.9092[/math] Состояние 2 [math]0.5752[/math] [math]0.4248[/math]

Для оценки начальной вероятности, мы предполагаем, что все последовательности начаты со скрытого состояния [math]S_1[/math] и рассчитаны с высокой вероятностью, а затем повторяем для [math]S_2[/math]. После нормализации получаем обновленный исходный вектор.

Повторяем эти шаги до тех пор, пока вероятности не сойдутся.

Псевдокод

 // T — конечный момент времени
 int[] DynamicOptionalStateSequance([math]\lambda[/math], d):
     double [math]\gamma[/math][1, i] = [math]\pi[/math][i] * b[i, d[1]]
     int [math]\psi[/math][1, i] = []
     int ans[]
     for t = 2 to T
         for i = 1 to n
             if [math]\gamma[/math][t, j] < [math]\gamma[/math][t - 1, i] * a[i, j] * b[j, d[t]]
                 [math]\gamma[/math][t, j] = [math]\gamma[/math][t - 1, i] * a[i, j] * b[j, d[t]]
                 [math]\psi[/math][t, j] = i
     ans[T] = 1 
     for i = 2 to n
         if [math]\gamma[/math][T, i] > [math]\gamma[/math][T, i - 1]
             ans[T] = i    
     for t = T - 1 downto 1
         ans[t] = [math]\psi[/math][t + 1, ans[t + 1]]
 return ans

См. также

• Скрытые Марковские модели
• Алгоритм "Вперед-Назад"
• Алгоритм Витерби

Источники информации

• Википедия — Алгоритм Баума-Велша

• Лекция "Скрытые Марковские Модели" Сергея Николенко

• Wikipedia — Baum–Welch algorithm