机器学习笔记补充——隐马尔可夫模型

0. 马尔可夫模型介绍

1. 隐马尔可夫模型的基本概念

  隐马尔可夫模型(hidden Markov model,HMM)是可用于标注问题的统计学习模型，描述隐藏的马尔可夫链随机生成观测序列的过程，属于生成模型。在语音识别、自然语言处理、生物信息、模式识别等领域有广泛应用。

1.1 隐马尔可夫模型的定义

  **定义：**隐马尔可夫模型是关于时序的概率模型，描述由一个隐藏的马尔可夫链随机生成不可观测的状态随机序列，再由各个状态生成一个观测从而产生观测随机序列的过程。隐藏的马尔可夫链随机生成的状态的序列，称为状态序列(state sequence)；每个状态生成一个观测，而由此产生的观测的随机序列，称为观测序列(observation sequence)。序列的每个位置又可以看作是一个时刻。
  定义很概括很抽象，可以结合之后的例子进行理解。

1.2 符号说明：

• Q = { q 1 , q 2 , . . . , q N } Q=\{q_1,q_2,...,q_N\} Q={q1​,q2​,...,qN​}是所有可能的状态的集合；
• V = { v 1 , v 2 , . . . , v M } V=\{v_1,v_2,...,v_M\} V={v1​,v2​,...,vM​}是所有可能的观测的集合；
• $I=(i_1,i_2,...,i_T)$是长度为T的状态序列；
• $O=(o_1,o_2,...,o_T)$是对应I的观测序列；
• $A=[a_{ij}{]}_{N\times N}$是状态转移概率矩阵，其中$a_{ij}=P(i_{t+1}=q_j|i_t=q_i),i=1,2,...,N;j=1,2,...,N$表示在时刻t时，从状态$q_i$转移到时刻t+1的状态$q_j$的概率。
• $B=[b_j(k){]}_{N\times M}$是观测概率矩阵，其中$b_j(k)=P(o_t=v_k|i_t=q_j),k=1,2,...,M;j=1,2,...,N$表示在t时刻处于状态$q_j$的条件下，生成观测量$v_k$的概率。
• $\pi$是初始状态概率向量：$\pi =({\pi }_i)$，其中${\pi }_i=P(i_1=q_i),i=1,2,...,N$是时刻t=1处于状态$q_i$的概率。

   隐马尔可夫模型可以由$\pi ,A$决定状态序列，B决定观测序列，三者称为隐马尔可夫三要素。所以，隐马尔可夫模型的参数可以用如下方式表示：
$$\lambda =(A,B,\pi )$$

1.3 隐马尔可夫模型的前提假设

  1).齐次马尔可夫性假设：隐藏的马尔可夫链在任意时刻t的状态只依赖于其前一时刻的状态，与其他时刻的状态及观测无关，也与时刻t无关：
$$P(i_t|i_{t-1},o_{t-1},...,i_1,o_1)=P(i_t|i_{t-1}),t=1,2,...,T$$
  2).观测独立性:假设任意时刻的观测只依赖于该时刻的马尔可夫链的状态，与其他观测及状态无关：
$$P(o_t|i_T,o_T,i_{T-1},o_{T-1},...,i_{t+1},o_{t+1},i_t,i_{t-1},o_{t-1},...,i_1,o_1=P(o_t|i_t)$$

1.4 隐马尔科夫模型的示例

  已知：

盒子编号	1	2	3	4
红球数	5	3	6	8
白球数	5	7	4	2

  按照下述方式抽球，产生一个描述每次抽球颜色的观测序列：
  1).开始，从4个盒子中以等概率随机选取1个盒子，从这个盒子里随机抽出1个球，记录其颜色然后放回；
  2).然后，从当前盒子随机转移到下一个盒子，规则是：如果当前盒子是盒1，那么下一个盒子一定是盒2；如果当前盒子是盒2或3，那么分别以0.4和0.6的概率转移到左边或右边的盒子；如果当前是盒4，那么各以0.5的概率留在盒4或转移到盒3；
  3).确定转移的盒子后，再从该盒子里抽球，记录颜色后放回；
  4).重复5次上述过程，得到观测序列：
$$O=(红，红，白，白，红)$$
  在上述问题中，隐马尔可夫模型的状态集合为：
Q = { 盒 1 ， 盒 2 ， 盒 3 ， 盒 4 } ， N = 4 Q=\{盒1，盒2，盒3，盒4\}，N=4 Q={盒1，盒2，盒3，盒4}，N=4
  观测集合是：
V = { 红 ， 白 } ， M = 2 V=\{红，白\}，M=2 V={红，白}，M=2
  状态序列盒观测序列的长度均为T=5；
  初始概率分布为：$\pi =(0.25,0.25,0.25,0.25{)}^T$
  状态转移概率分布为：
$$A=\left[\begin{array}{cccc}0& 1& 0& 0\\ 0.4& 0& 0.6& 0\\ 0& 0.4& 0& 0.6\\ 0& 0& 0.5& 0.5\end{array}\right]$$
  观测概率分布为：
$$B=\left[\begin{array}{cc}0.5& 0.5\\ 0.3& 0.7\\ 0.6& 0.4\\ 0.8& 0.2\end{array}\right]$$
  所以，根据上述分析可以得到，该案例中观测序列的生成过程为：
  1).按照初始状态分布$\pi$产生状态$i_1$；
  2).令t=1；
  3).按照状态$i_t$的观测概率分布$b_{i_t}(k)$生成$o_t$；
  4).按照状态$i_t$的状态转移概率分布 { a i t i t + 1 } \{a_{i_ti_{t+1}}\} {ait​it+1​​}产生状态$i_{t+1},i_{t+1}=1,2,...,N$；
  5).令t=t+1，如果t<T，重复上述步骤3,4，否则终止。

1.5 隐马尔可夫模型的3个基本问题

  研究隐马尔可夫模型通常有三个基本问题：
  1).概率计算问题：给定模型$\lambda =(A,B,\pi )$和观测序列$O=(o_1,o_2,...,o_T)$，计算在模型$\lambda$下观测序列O出现的概率$P(O|\lambda )$
  2).学习问题：已知观测序列$O=(o_1,o_2,...,o_T)$，估计模型$\lambda =(A,B,\pi )$参数，使得该模型产生观测序列概率$P(O|\lambda )$最大。用极大似然估计的方法估计参数。
  3).预测问题：也称为解码问题。已知模型$\lambda =(A,B,\pi )$和观测序列$O=(o_1,o_2,...,o_T)$，求对给定观测序列条件概率$P(I|O)$最大的状态序列$I=(i_1,i_2,...,i_T)$。即给定观测序列，求最有可能的对应状态序列。

2.概率计算问题算法

  理论上，要解决隐马尔可夫的概率计算问题，直接根据已知条件$\lambda =(A,B,\pi )$和$O=(o_1,o_2,...,o_T)$直接计算就可以得到。但是，这种做法实际上并不可行，因为算法时间复杂度是$O(T{N}^T)$，在数据量较大的情况下是不行的。所以，通常使用下述算法进行计算：前向-后向算法。

2.1 前向算法

  前向概率的定义：给定隐马尔可夫模型$\lambda$，定义到时刻t的部分观测序列为$o_1,o_2,...,o_t$且状态为$q_i$的概率为前向概率，记作：
$${\alpha }_t(i)=P(o_1,o_2,...,o_t,i_t=q_i|\lambda )$$
  可以用递推法求前向概率${\alpha }_t(i)$和观测序列概率$P(O|\lambda )$：

• 输入：隐马尔可夫模型参数$\lambda$和观测序列O；
• 输出：观测序列概率$P(O|\lambda )$；
  1). 计算前向概率的初值：
  $${\alpha }_1(i)={\pi }_i{b}_i(o_1),i=1,2,...,N$$
  2). 递推公式进行递推，对$t=1,2,…,T-1:
  $${\alpha }_{t+1}(i)=[\sum _{j=1}^N{\alpha }_t(j)a_{ji}]b_i(o_{t+1}),i=1,2,...,N$$
  3). 递推结束后，计算观测序列概率：
  $$P(O|\lambda )=\sum _{i=1}^N{\alpha }_T(i)$$
  说明：
  1).前向概率的初值就是模型分别处于各个状态$q_i$下，观测到$o_1$的概率；
  2).递推公式表达的就是，对应状态$q_i$的第t+1刻的前向概率是，第t刻的所有前向概率分别从自己的状态$q_j$转移到状态$q_i$的概率之和，再乘以由第t+1刻的状态$q_i$观测到$o_{t+1}$的概率，如下图：
  3).第三步是因为${\alpha }_T(i)=P(o_1,o_2,...,o_T,i_T=q_i|\lambda$，表示观测到序列$o_1,o_2,...,o_T$的最后一项$o_T$是由状态$i_T=q_i$状态观测得到的概率。所以将每种状态的概率相加就有$P(O|\lambda )={\sum }_{i=1}^N{\alpha }_T(i)$。
  4).易知，上述算法的时间复杂度为$O(N^{2}T)$，相比于直接计算的指数级复杂度，已经属于可接受范围内。

  案例可以参见李航老师的《统计学习方法》教材，由于上述递推过程理解难度不大，故不附上案例。

2.2 后向算法

  后向概率的定义：给定隐马尔可夫模型$\lambda$，定义在时刻t且状态为$q_i$的条件下，从t+1到T的部分观测序列为$o_{t+1},o_{t+2},...,o_T$的概率为后向概率，记作：
$${\beta }_t(i)=P(o_{t+1},o_{t+2},...,o_T|i_t=q_i,\lambda )$$
  可以用递推法求后向概率${\beta }_t(i)$和观测序列概率$P(O|\lambda )$：

• 输入：隐马尔可夫模型参数$\lambda$和观测序列O；
• 输出：观测序列概率$P(O|\lambda )$；
  1). 计算前向概率的初值：
  $${\beta }_T(i)=1,i=1,2,...,N$$
  2). 递推公式进行递推，对$t=T-1,T-2,…,1:
  $${\beta }_t(i)=\sum _{j=1}^N{a}_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j),i=1,2,...,N$$
  3). 递推结束后，计算观测序列概率：
  $$P(O|\lambda )=\sum _{i=1}^N{\pi }_i{b}_i(o_1){\beta }_1(i)$$
  说明：
  1).后向概率的初值是人为规定的；
  2).递推公式表达的就是，从后向前递推，${\beta }_t(i)$其实就是在t+1时刻的后向概率情况下，增加上从第t时刻的$q_i$状态分别转移到第t+1时刻的$q_j$状态，并观测到$o_{t+1}$的概率，如下图：
  3).第三步是因为$P(O|\lambda )={\beta }_0(i)$。所以将每种状态的0时刻后向概率相加就有$P(O|\lambda )={\sum }_{i=1}^N{\pi }_i{b}_i(o_1){\beta }_1(i)$。
  4).易知，上述算法的时间复杂度为$O(N^{2}T)$，相比于直接计算的指数级复杂度，已经属于可接受范围内。

2.3 一些概率与期望值的计算

  1).给定模型$\lambda$和观测O，在t时刻处于状态$q_i$的概率记为：
$${\gamma }_t(i)=P(i_t=q_i|O,\lambda )$$
  该概率可以通过前向后向概率计算：
$${\gamma }_t(i)=P(i_t=q_i|O,\lambda )=\frac{P(i_t=q_i,)|\lambda }{P(O|\lambda )}=\frac{{\alpha }_t(i){\beta }_t(i)}{P(O|\lambda )}=\frac{{\alpha }_t(i){\beta }_t(i)}{{\sum }_{j=1}^N{\alpha }_t(j){\beta }_t(j)}$$
  其中，根据前向概率和后向概率的定义可知，${\alpha }_t(i){\beta }_t(i)=P(i_t=q_i,O|\lambda )$，且$P(O|\lambda )={\sum }_{i=1}^{N}P(i_t=q_i,O|\lambda )$，从而可以进行上述推导。

  2).给定模型$\lambda$和观测O，在t时刻处于状态$q_i$且在t+1时刻处于状态$q_j$的概率记为：
$${ϵ}_t(i,j)=P(i_t=q_i,i_{t+1}=q_j|O,\lambda )$$
  该概率可以通过前向后向概率计算：
$${ϵ}_t(i,j)=\frac{P(i_t=q_i,i_{t+1}=q_j,O|\lambda )}{P(O|\lambda )}=\frac{P(i_t=q_i,i_{t+1}=q_j,O|\lambda )}{{\sum }_{i=1}^N{\sum }_{j=1}^{N}P(i_t=q_i,i_{t+1}=q_j,O|\lambda )}=\frac{{\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)}{{\sum }_{i=1}^N{\sum }_{j=1}^N{\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)}$$
  其中，$P(i_t=q_i,i_{t+1}=q_j,O|\lambda )={\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)$

  3).一些有用的期望值：
  i.在观测O下状态i出现的期望值：${\sum }_{t=1}^T{\gamma }_t(i)$；
  ii.在观测O下由状态i转移的期望值：${\sum }_{t=1}^{T-1}{\gamma }_t(i)$；
  iii.在观测O下由状态i转移到状态j的期望值：${\sum }_{t=1}^{T-1}{ϵ}_t(i,j)$

3.学习问题算法

  隐马尔可夫模型参数的学习，可以根据已知条件是只包含观测序列还是也包含对应状态序列划分为监督学习和无监督学习。

3.1 监督学习隐马尔可夫模型

  监督学习方法很简单，就是已知数据集 { ( O 1 , I 1 ) , ( O 2 , I 2 ) , . . . , ( O S , I S ) } \{(O_1,I_1),(O_2,I_2),...,(O_S,I_S)\} {(O1​,I1​),(O2​,I2​),...,(OS​,IS​)}。
  计算样本中t时刻处于状态i而t+1时刻处于状态j的样本频数为$A_{ij}$，那么状态转移概率的估计是：
$${\hat{a}}_{ij}=\frac{A_{ij}}{{\sum }_{j=1}^N{A}_{ij}},i=1,2,...,N;j=1,2,...,N$$
  计算样本中状态为j并观测为k的样本频数为$B_{jk}$，那么状态为j观测为k的概率估计是：
$${\hat{b}}_j(k)=\frac{B_{jk}}{{\sum }_{k=1}^M{B}_{jk}},j=1,2,...,N;k=1,2,...,M$$
  初始状态概率${\pi }_i$的估计值就是样本中初始状态为$q_i$的频数。
  有监督学习学习隐马尔可夫模型比较简单，因为其中的隐藏信息其实已经显示出来了。但是有监督学习需要标注样本，这个人工成本非常高，在实际情况中基本不会使用。

3.2 无监督学习隐马尔可夫模型——Baum-Welch算法(EM算法)

  由EM算法的特性可知，EM算法可以用来估计含有隐变量的马尔可夫模型，且很明显此处可观测变量为观测序列O，不可观测的隐变量为状态序列I，则隐马尔科夫模型可以写成如下形式：
$$P(O|\lambda )=\sum _{I}P(O|I,\lambda )P(I|\lambda )$$
  具体算法步骤如下：

• 输入：观测数据$O=(o_1,o_2,...,o_T)$；
• 输出：隐马尔可夫模型参数；
  1). 初始化：对n=0,取$a_{ij}^{(0)},b_j(k{)}^{(0)},{\pi }_i^{(0)}$，得到初始模型参数${\lambda }^{(0)}=(A^{(0)},B^{(0)},{\pi }^{(0)})$；
  2). 递推：对n=1,2,…,求解：
  $$a_{ij}^{(n+1)}=\frac{{\sum }_{t=1}^{T-1}{ϵ}_t(i,j)}{{\sum }_{t=1}^{T-1}{\gamma }_t(i)}$$
  $$b_j(k{)}^{(n+1)}=\frac{{\sum }_{t=1,o_t=v_k}^T{\gamma }_t(j)}{{\sum }_{t=1}^T{\gamma }_t(i)}$$
  $${\pi }_i^{(n+1)}={\gamma }_1(i)$$
  上式右端各值按已知观测序列O和当前模型参数${\lambda }^{(n+1)}=(A^{(n+1)},B^{(n+1)},{\pi }^{(n+1)})$
  计算，${\gamma }_t(j),{ϵ}_t(i,j)$按2.3节的公式计算。
  3).算法结束，得到模型参数：${\lambda }^{(n+1)}=(A^{(n+1)},B^{(n+1)},{\pi }^{(n+1)})$

4.预测问题算法

4.1 近似算法

  近似算法的思路很简单，就是在每个时刻t都选择当前时刻概率最大的状态$i_t^{\ast }$作为预测的状态，以此类推，可以得到预测状态序列$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$。
  假设已知隐马尔可夫模型$\lambda$和观测序列O，在时刻t处处于状态$q_i$的概率是：
$${\gamma }_t(i)=\frac{{\alpha }_t(i){\beta }_t(i)}{P(O|\lambda )}=\frac{{\alpha }_t(i){\beta }_t(i)}{{\sum }_{j=1}^N{\alpha }_t(j){\beta }_t(j)}$$
  选择当前时刻t概率最大的状态作为预测的t时刻状态：
$$i_t^{\ast }=\arg \underset{1\le i\le N}{\max }[{\gamma }_t(i)],t=1,2,\dots ,T$$
  从而得到预测的状态序列$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$。
  近似算法的优点是非常直观，实现也简单，缺点是只能保证某时刻的状态是当前时刻下的最大概率状态，而不能保证最终得到的预测状态序列整体是概率最大的序列。另外，这种算法前后预测值之间没有联系和约束，很可能会产生两个相邻状态之间状态转移概率为0的情况($a_{ij}=0$)。尽管如此，这种算法在不少场合仍是非常有用的。

4.2 维特比算法(Viterbi algorithm)

  维特比算法实质上是用 动态规划 解决隐马尔可夫模型的预测问题，即用动态规划求解概率最大路径，是个最优路径问题。
  熟悉动态规划算法的读者可能知道，动态规划的特点就是：

• “无后效性”：后续发生的状态不影响前面已经求得的状态序列最优解；
• “全局最优解的局部解也是局部最优解”：也就是说，假设求出的最优序列是$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$，那么它的局部序列$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_k^{\ast })，k<T$ 肯定是状态1~k的最优序列，因为如果不是的话，完全可以将这部分序列换成最优子序列，从而可以得到更优的全局序列，这与最开始假设的最优序列是矛盾的。

  维特比算法的思想就是如此，只是它的最优子序列是从T时刻开始往前取的，实质是一样的。先说明一下两个符号：
  1).${\delta }_t(i)$表示定义在时刻 t 状态为集合{i|i=1,2,…,N}中任意一种的所有可能路径中，概率最大的那一条路径的概率值。
$${\delta }_t(i)=\underset{i_1,i_2,....,i_{t-1}}{\max }P(i_t=i,i_{t-1},...,i_1,o_t,...,o_1|\lambda ),i=1,2,...,N$$
$$递推公式：{\delta }_{t+1}(i)=\underset{i_1,i_2,....,i_t}{\max }P(i_t=i,i_t,...,i_1,o_{t+1},...,o_1|\lambda )$$
$$=\underset{1\le j\le N}{\max }[{\delta }_t(j)a_{ji}]b_i(o_{t+1}),i=1,2,...,N；t=1,2,...,T-1$$
  其中$a_{ji}$是状态i转移到状态j的转移概率矩阵元素；$b_i(o_{t+1})$是状态i的观测概率矩阵元素。

  2).${\Psi }_t(i)$是在时刻 t 状态为 i 的所有单个路径$(i_1,i_2,\dots ,i_{t-1},i)$中概率最大路径的第t-1个结点的状态，即可以简单理解为：${\Psi }_t(i)=i_{t-1}$，但是这个$i_{t-1}$是有上述前提条件的。
$${\Psi }_t(i)=\arg \underset{1\le j\le N}{\max }[{\delta }_{t-1}(j)a_{ji}],i=1,2,\dots ,N$$

  维特比算法的具体步骤：

• 输入：马尔可夫模型$\lambda =(A,B,\pi )$和观测$O=(o_1,o_2,...,o_T)$；
• 输出：最优预测序列(最优路径)：$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$；
• 1).初始化
  $${\delta }_1(i)={\pi }_i{b}_i(o_1),{\Psi }_1(i)=0,i=1,2,...,N$$
• 2).对$t=2,3,...,T$进行如下递推：
  $${\delta }_t(i)=\underset{1\le j\le N}{\max }[{\delta }_{t-1}(j)a_{ji}]b_i(o_t),i=1,2,...,N$$
  $${\Psi }_t(i)=\arg \underset{1\le j\le N}{\max }[{\delta }_{t-1}(j)a_{ji}],i=1,2,...,N$$
• 3).得到最优预测序列概率$P^{\ast }$和T时刻最优状态$i_T^{\ast }$：
  $$P^{\ast }=\underset{1\le i\le N}{\max }{\delta }_T(i),i_T^{\ast }=\arg \underset{1\le i\le N}{\max }[{\delta }_T(i)]$$
• 4).最优路径回溯：
  $$i_t^{\ast }={\Psi }_{t+1}(i_{t+1}^{\ast }),t=T-1,T-2,...,1$$
  从而得到最优预测序列：$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$。

  维特比算法的案例也很简单很直观，读者可以自行参阅李航老师《统计学习方法》第10章，以巩固加深对维特比算法的理解。