鲍姆-韦尔奇算法

在电气工程、计算机科学、统计计算和生物信息学中，鲍姆-韦尔奇算法是用于寻找隐马尔可夫模型未知参数的最大期望算法，它利用前向-后向算法来计算E-Step的统计信息。

 历史 
鲍姆-韦尔奇算法是以其发明者伦纳德·埃绍·鲍姆和劳埃德·理查德·韦尔奇的名字命名的。鲍姆-韦尔奇算法和隐马尔可夫模型在20世纪60年代末和70年代初由鲍姆和他的同事在国防分析研究所的一系列文章中首次描述。HMMs最初主要应用于语音处理领域。20世纪80年代，HMMs开始成为分析生物系统和信息，特别是遗传信息的有用工具。此后，它们成为基因组序列概率建模的重要工具。

 介绍 
隐马尔可夫模型描述了一组“隐含”变量和可观测到的离散随机变量的联合概率。它依赖于假设：第个隐藏变量只与第个隐含变量相关，而与其他先前的隐藏变量无关，而当前观测到的状态仅依赖于当前的隐藏状态。

鲍姆-韦尔奇算法利用最大期望算法，在给定一组观测特征向量的情况下，求出隐马尔可夫模型参数的最大似然估计。

记离散的隐含随机变量为，它总共有种状态(有个不同的值)。设与时间无关，得到与时间无关的随机变量转移矩阵：

初始的状态(即)分布由下式给出：

记观测到的变量为，总共有种取值。同样假设由隐含变量得到的可观测变量与时间无关。在时间，由隐含变量得到的可观察变量的概率是：

由所有可能得和的取值，我们可以得到的矩阵，其中属于所有可能得隐含状态，属于所有的可观测状态。

给出可观测序列：。

我们可以用描述隐马尔科夫链，鲍姆-韦尔奇算法寻找的局部极大值，也就是能够使得观测到的序列出现的概率最大的HMM的参数。

 算法 
初始化参数，可以随机初始化，或者根据先验知识初始化。

 前向过程 
记是参数的条件下，观测的序列是，时刻的状态是的概率。可以通过递归计算：




 后向过程 
记是参数是，在时刻的状态是的条件下，余下部分的观测序列是的概率。




 更新 

 根据贝叶斯公式计算临时变量。
在给定观测序列和参数的情况下，在时间状态是的概率:
在给定观测序列和参数的情况下，在时间状态是，在时间状态是的概率:
和的分母一样，表示给定参数得到观测序列的概率。

 然后更新参数：
，在时间状态是的概率
，等于期望的从状态转换到状态的数量除以从状态开始的转换的总数。
，其中，是期望的从状态得到的观察值等于的数量除以从 状态开始的转换的总数。
重复上面的步骤直到收敛。算法可能过拟合，也不保证收敛到全局最大值。
其中计算和相当于最大期望算法的E-Step，而更新的过程相当于最大期望算法的M-Step。

 例子 
假设我们有一只会下蛋的鸡，每天中午我们都会去拾取鸡蛋。而鸡是否下蛋依赖于一些未知的隐含状态，这里我们简单的假设只有两种隐含状态会决定它是否下蛋。我们不知道这些隐含状态的初始值，不知道他们之间的转换概率，也不知道在每种状态下鸡会下蛋的概率。我们随机初始化他们来开始猜测。

+Transition  State 1  State 2 State 10.5  0.5 State 20.3  0.7+ Emission  No Eggs  Eggs State 10.3  0.7 State 20.8  0.2+ Initial State 10.2 State 20.8
假设我们得到的观测序列是(E=eggs, N=no eggs): N, N, N, N, N, E, E, N, N, N。

这样我们同时也得到了观测状态的转移：NN, NN, NN, NN, NE, EE, EN, NN, NN。

通过上面的信息来重新估计状态转移矩阵。
 Observed sequence  Probability of sequence and state is  then   Highest Probability of observing that sequence NN  0.024  0.3584  , NN  0.024  0.3584  , NN  0.024  0.3584  , NN  0.024  0.3584  , NE  0.006  0.1344  , EE  0.014  0.0490  , EN  0.056  0.0896  , NN  0.024  0.3584  , NN  0.024  0.3584  , Total 0.22  2.4234 

重新估计到的转移概率为(下表中的"Pseudo probabilities")，重新计算所有的转移概率，得到下面的转移矩阵：

+Old Transition Matrix  State 1  State 2 State 10.5  0.5 State 20.3  0.7+New Transition Matrix (Pseudo Probabilities)  State 1  State 2 State 10.0598  0.0908 State 20.2179  0.9705+New Transition Matrix (After Normalization)  State 1  State 2 State 10.3973  0.6027 State 20.1833  0.8167

接下来重新估计Emission Matrix:
 Observed Sequence  Highest probability of observing that sequence if E is assumed to come from  Highest Probability of observing that sequence NE  0.1344  ,  0.1344  , EE  0.0490  ,  0.0490  , EN  0.0560  ,  0.0896 , Total 0.2394   0.2730 

重新估计从隐含状态得到观察结果E的概率是，得到新的Emission Matrix

+Old Emission Matrix  No Eggs  Eggs State 10.3  0.7 State 20.8  0.2+New Emission Matrix (Estimates)  No Eggs  Eggs State 10.0876  0.8769 State 21.0000  0.7385+New Emission Matrix (After Normalization)  No Eggs  Eggs State 10.0908  0.9092 State 20.5752  0.4248

为了估计初始状态的概率，我们分别假设序列的开始状态是和，然后求出最大的概率，再归一化之后更新初始状态的概率。

一直重复上面的步骤，直到收敛。

 代码 
from hmmlearn import hmm
import numpy as np

X = np.array([1, 1, 1, 1, 1, 0, 0, 1, 1, 1]).reshape(-1, 1)
model = hmm.GaussianHMM(n_components=2, covariance_type='full')
model.fit(X)

model.monitor_.history

# pi
model.startprob_

# state transform matrix
model.transmat_

# emission_matrix
np.power(np.e, model._compute_log_likelihood(np.unique(X).reshape(-1, 1)))

Category:估计理论
Category:算法