机器学习算法笔记：HMM隐马尔可夫模型

Hidden Markov Model

隐马尔可夫模型是一种概率图模型。机器学习模型可以从频率派和贝叶斯派两个方向考虑，在频率派的方法中的核心是优化问题，而在贝叶斯派的方法中，核心是积分问题，发展出一系列的积分方法：如变分推断，MCMC 等。

概率图模型最基本的模型可分为有向图（贝叶斯网络）和无向图（马尔可夫随机场），如果样本之间存在关联，可认为样本中附带时序信息，使得样本间不独立全同分布，这种模型就是动态模型，隐变量随着时间发生变化，观测变量也随之变化：

根据状态变量的特点，可以分为：

• HMM，状态变量（隐变量）是离散的，观测变量没有要求
• Kalman 滤波，状态变量是连续的，线性的
• 粒子滤波，状态变量是连续，非线性的

模型定义

HMM 用概率图表示为：

上图表示了三个时刻的隐变量变化，用 $\lambda =(\pi ,A,B)$ 表示，其中 $\pi$ 为初始概率分布，$A$ 为状态转移矩阵，$B$ 为发射矩阵。

• $o_t$ 表示观测变量，$O$ 为观测序列，$V=-v_1,v_2,\cdots ,v_M$ 表示观测值域（能取的值），
• $i_t$ 表示状态变量，$I$ 为状态序列，$Q=q_1,q_2,\cdots ,q_N$ 表示状态变量值域，
• π = { π ( 1 ) , π ( 2 ) , ⋯   , π ( N ) } {\color{blue}\pi}=\{\pi(1),\pi(2),\cdots,\pi(N)\} π={π(1),π(2),⋯,π(N)}表示初始状态，${\sum }_{i=1}^N\pi (i)=1$
• $A=(a_{ij}=p(i_{t+1}=q_j|i_t=q_i))$ 表示状态转移矩阵
• $B=(b_j(k)=p(o_t=v_k|i_t=q_j))$ 表示发射矩阵。
• ${\lambda }^{(t)}=({\pi }^{(t)},A^{(t)},B^{(t)})$ 表示 $t$ 时刻参数

两个基本假设

• 齐次 Markov 假设： $$p(i_{t+1}|i_t,i_{t-1},\cdots ,i_1,o_t,o_{t-1},\cdots ,o_1)=p(i_{t+1}|i_t)$$

• 观测独立假设： $$p(o_t|i_t,i_{t-1},\cdots ,i_1,o_{t-1},\cdots ,o_1)=p(o_t|i_t)$$

三个任务

HMM的三个任务：

• 1、Evaluation：$p(O|\lambda )$，Forward-Backward 算法
• 2、Learning： $\lambda =\underset{\lambda }{\mathrm{argmax}}p(O|\lambda )$，EM 算法（Baum-Welch）
• 3、Decoding：$I=\underset{I}{\mathrm{argmax}}p(I|O,\lambda )$，Vierbi 算法$$\begin{array}{rl}& \text{预测问题：}p(i_{t+1}|o_1,o_2,\cdots ,o_t)\\ & \text{滤波问题：}p(i_t|o_1,o_2,\cdots ,o_t)\end{array}$$

Evaluation

Evaluation 问题可表示为：
$$p(O|\lambda )=\sum _{I}p(I,O|\lambda )=\sum _{I}p(O|I,\lambda )p(I|\lambda )$$

上式中：
$$p(I|\lambda )=p(i_1,i_2,\cdots ,i_t|\lambda )=p(i_t|i_1,i_2,\cdots ,i_{t-1},\lambda )p(i_1,i_2,\cdots ,i_{t-1}|\lambda )$$

根据齐次 Markov 假设：
$$p(i_t|i_1,i_2,\cdots ,i_{t-1},\lambda )=p(i_t|i_{t-1})=a_{i_{t-1}{i}_t}$$

所以：
$$p(I|\lambda )={\pi }_{i_1}\prod _{t=2}^T{a}_{i_{t-1},i_t}$$ 又由于： $$p(O|I,\lambda )=\prod _{t=1}^T{b}_{i_t}(o_t)$$

于是：
$$\begin{array}{rl}p(O|\lambda )& =\sum _I{\pi }_{i_1}\prod _{t=2}^T{a}_{i_{t-1},i_t}\prod _{t=1}^T{b}_{i_t}(o_t)\\ & =\underset{\text{复杂度}O(N^T)}{\underbrace{\sum _i\cdots \sum _T}}{\pi }_{i_1}\prod _{t=2}^T{a}_{i_{t-1},i_t}\prod _{t=1}^T{b}_{i_t}(o_t)\end{array}$$

上式 $I$ 为状态序列，序列长度为 $T$，每个状态变量 $i$ 有 $N$ 种状态，所以其复杂度为 $O(N^T)$，计算量太大，要用更加简单的算法计算 $p(O|\lambda )$ 。

前向算法

记 ${\alpha }_t(i)=p(o_1,o_2,\cdots ,o_t,i_t=q_i|\lambda )$，所以 ${\alpha }_T(i)=p(O,i_T=q_i|\lambda )$则：
$$p(O|\lambda )=\sum _{i=1}^{N}p(O,i_T=q_i|\lambda )=\sum _{i=1}^N{\alpha }_T(i)$$

上式中，对状态变量 $i_T$ 在其状态空间上求和 ∑ i = 1 N { i T = q i } = 1 \sum^N\limits_{i=1}\{i_T=q_i\}=1 i=1∑N​{iT​=qi​}=1。

推导 ${\alpha }_{t+1}(j)$与${\alpha }_t(i)$间的递推关系：
$$\begin{array}{rl}{\alpha }_{t+1}(j)& =p(o_1,o_2,\cdots ,o_{t+1},i_{t+1}=q_j|\lambda )\\ & =\sum _{i=1}^{N}p(o_1,o_2,\cdots ,o_{t+1},i_{t+1}=q_j,i_t=q_i|\lambda )\\ & =\sum _{i=1}^{N}p(o_{t+1}|o_1,o_2,\cdots ,i_{t+1}=q_j,i_t=q_i|\lambda )p(o_1,\cdots ,o_t,i_t=q_i,i_{t+1}=q_j|\lambda )\\ & =\sum _{i=1}^{N}p(o_{t+1}|i_{t+1}=q_j)p(o_1,\cdots ,o_t,i_t=q_i,i_{t+1}=q_j|\lambda )\text{观测独立假设}\\ & =\sum _{i=1}^N\underset{b_j(o_t)}{\underbrace{p(o_{t+1}|i_{t+1}=q_j)}}\underset{p(i_{t+1}=q_j|i_t=q_i,\lambda )=a_{ij}}{\underbrace{p(i_{t+1}=q_j|o_1,\cdots ,o_t,i_t=q_i,\lambda )}}\underset{{\alpha }_t(i)}{\underbrace{p(o_1,\cdots ,o_t,i_t=q_i|\lambda )}}\\ & =\sum _{i=1}^N{b}_j(o_t)a_{ij}{\alpha }_t(i)\end{array}$$

利用齐次 Markov 假设得到递推公式，也称为前向算法。

后向算法

定义：
$$\begin{gathered} {\beta }_t(i)=p(o_{t+1},\cdots ，o_T|i_t=q_i,\lambda ) \\ ⋮ \\ {\beta }_1(i)=p(o_2,\cdots ，o_T|i_1=q_i,\lambda ) \end{gathered}$$

则有：
$$\begin{array}{rl}p(O|\lambda )& =p(o_1,\cdots ,o_T|\lambda )\\ & =\sum _{i=1}^{N}p(o_1,o_2,\cdots ,o_T,i_1=q_i|\lambda )\text{引入i1}\\ & =\sum _{i=1}^{N}p(o_1,o_2,\cdots ,o_T|i_1=q_i,\lambda )\underset{\text{初始概率分布}{\pi }_i}{\underbrace{p(i_1=q_i)}}\\ & =\sum _{i=1}^{N}p(o_1,o_2,\cdots ,o_T|i_1=q_i,\lambda ){\pi }_i\\ & =\sum _{i=1}^N\underset{\text{观测独立假设}}{\underbrace{p(o_1|o_2,\cdots ,o_T,i_1=q_i,\lambda )}}\underset{{\beta }_1(i)}{\underbrace{p(o_2,\cdots ,o_T|i_1=q_i,\lambda )}}{\pi }_i\\ & =\sum _{i=1}^N{b}_i(o_1){\pi }_i{\beta }_1(i)\end{array}$$

${\beta }_1(i)$ 同样可以通过递推公式得到：
$$\begin{array}{rl}{\beta }_t(i)& =p(o_{t+1},\cdots ,o_T|i_t=q_i)\\ & =\sum _{j=1}^{N}p(o_{t+1},o_{t+2},\cdots ,o_T,i_{t+1}=q_j|i_t=q_i)\\ & =\sum _{j=1}^{N}p(o_{t+1},\cdots ,o_T|i_{t+1}=q_j,i_t=q_i)p(i_{t+1}=q_j|i_t=q_i)\\ & =\sum _{j=1}^{N}p(o_{t+1},\cdots ,o_T|i_{t+1}=q_j)a_{ij}\\ & =\sum _{j=1}^{N}p(o_{t+1}|o_{t+2},\cdots ,o_T,i_{t+1}=q_j)p(o_{t+2},\cdots ,o_T|i_{t+1}=q_j)a_{ij}\\ & =\sum _{j=1}^N{b}_j(o_{t+1})a_{ij}{\beta }_{t+1}(j)\end{array}$$

上述两种算法的复杂度均为 $O(T{N}^2)$，使得计算量大大降低。

Learning

为了学习得到参数的最优值，要用到 MLE ： $${\lambda }_{MLE}=\underset{\lambda }{\mathrm{argmax}}p(O|\lambda )$$

上式难以直接求解，需要采用 EM 算法（这里也叫 Baum Welch 算法）进行迭代求解，EM 算法的迭代公式为： $${\theta }^{t+1}=\underset{\theta }{\mathrm{argmax}}{\int }_z\log p(X,Z|\theta )p(Z|X,{\theta }^t)dz$$

上式 $X$ 是观测变量，$Z$ 是隐变量序列，$\theta$ 为模型参数，分别与这里的 $O,I,\lambda$ 对应，于是：
$$\begin{array}{rl}{\lambda }^{t+1}& =\underset{\lambda }{\mathrm{argmax}}\sum _I\log p(O,I|\lambda )p(I|O,{\lambda }^t)\\ & =\underset{\lambda }{\mathrm{argmax}}\sum _I\log p(O,I|\lambda )\frac{p(O,I|{\lambda }^t)}{\boxed{p(O,{\lambda }^t)}}{}_{\text{框中O为λ无关项}}\\ & =\underset{\lambda }{\mathrm{argmax}}\sum _I\log p(O,I|\lambda )p(O,I|{\lambda }^t)\end{array}$$

上式中 $p(O|{\lambda }^t)$ 和$\lambda$ 无关，去掉该项不影响参数 $\lambda$ 的优化。

由 Evaluation 中的推导可知：
$$\begin{array}{rl}Q(\lambda ,{\lambda }^t)& =\sum _I\log p(O,I|\lambda )p(O,I|{\lambda }^t)\\ & =\sum _I[\log {\pi }_{i_1}+\sum _{t=2}^T\log a_{i_{t-1},i_t}+\sum _{t=1}^T\log b_{i_t}(o_t)]p(O,I|{\lambda }^t)\end{array}$$

提取上式中 $\boxed{\pi }$ 相关项：
$$\begin{array}{rl}{\pi }^{t+1}& =\underset{\pi }{\mathrm{argmax}}\sum _I[\log {\pi }_{i_1}p(O,I|{\lambda }^t)]\\ & =\underset{\pi }{\mathrm{argmax}}\sum _I[\log {\pi }_{i_1}\cdot p(O,i_1,i_2,\cdots ,i_T|{\lambda }^t)]\\ & =\underset{\pi }{\mathrm{argmax}}\sum _{i_1}\underset{\text{相当于对i2,i2,⋯,iT 求边缘分布}}{\underbrace{\sum _{i_2}\cdots \sum _{i_T}[\log {\pi }_{i_1}\cdot p(O,i_1,i_2,\cdots ,i_T|{\lambda }^t)]}}\\ & =\underset{\pi }{\mathrm{argmax}}\sum _{i_1}[\log {\pi }_{i_1}\cdot p(O,i_1|{\lambda }^t)]\end{array}$$

上式中：$i_1$ 有$N$个状态：$i_1=q_i$，$\pi$ 的约束条件 $st:\sum _i{\pi }_i=1$，可定义 Lagrange 函数： $$L(\pi ,\eta )=\sum _{i=1}^N\underset{l_i}{\underbrace{\log {\pi }_i\cdot p(O,i_1=q_i|{\lambda }^t)+\eta (\sum _{i=1}^N{\pi }_i-1)}}$$

对求和符号中的被求和项求偏导：
$$\frac{\partial l_i}{\partial {\pi }_i}=\frac{1}{{\pi }_i}p(O,i_1=q_i|{\lambda }^t)+\eta =0$$

对上式求和：
$$\begin{gathered} \frac{\partial L}{\partial {\pi }_i}=\sum _{i=1}^{N}p(O,i_1=q_i|{\lambda }^t)+{\pi }_i\eta =0 \\ \Rightarrow \eta =-p(O|{\lambda }^t) \end{gathered}$$

于是可得到：
$${\pi }_i^{t+1}=\frac{p(O,i_1=q_i|{\lambda }^t)}{p(O|{\lambda }^t)}$$

Decoding

Decoding 问题表示为：
$$I=\underset{I}{\mathrm{argmax}}p(I|O,\lambda )$$

就是找到一个由 $q_i$ 组成的序列，使得概率最大，可采用动态规划的思想求解

定义：
$${\delta }_t(i)=\underset{i_1,\cdots ,i_{t-1}}{\max }p(o_1,\cdots ,o_t,i_1,\cdots ,i_{t-1},i_t=q_i)$$

则： $${\delta }_{t+1}(j)=\underset{1\le i\le N}{\max }{\delta }_t(i)a_{ij}{b}_j(o_{t+1})$$

从上一步到下一步的概率再求最大值，记录序列的路径(每个时刻对应状态$q_i$的下标)：
$${\psi }_{t+1}(j)=\underset{1\le i\le N}{\mathrm{argmax}}{\delta }_t(i)a_{ij}$$

小结

HMM 是一种动态模型，是由混合树形模型和时序结合起来的一种模型（类似 GMM + Time）。对于类似 HMM 的这种状态空间模型，普遍的除了学习任务（采用 EM ）外，还有推断任务，推断任务包括：

• 解码 Decoding：$p(z_1,z_2,\cdots ,z_t|x_1,x_2,\cdots ,x_t)$

• 似然概率：$p(X|\theta )$

• 滤波：$p(z_t|x_1,\cdots ,x_t)$，Online $$p(z_t|x_{1:t})=\frac{p(x_{1:t},z_t)}{p(x_{1:t})}=C{\alpha }_t(z_t)$$

• 平滑：$p(z_t|x_1,\cdots ,x_T)$，Offline $$p(z_t|x_{1:T})=\frac{p(x_{1:T},z_t)}{p(x_{1:T})}=\frac{{\alpha }_t(z_t)p(x_{t+1:T}|x_{1:t},z_t)}{p(x_{1:T})}$$ 根据概率图的条件独立性，有： $$p(z_t|x_{1:T})=\frac{{\alpha }_t(z_t)p(x_{t+1:T}|z_t)}{p(x_{1:T})}=C{\alpha }_t(z_t){\beta }_t(z_t)$$ 这个算法叫做前向后向算法。

• 预测：$p(z_{t+1},z_{t+2}|x_1,\cdots ,x_t),p(x_{t+1},x_{t+2}|x_1,\cdots ,x_t)$ $$p(z_{t+1}|x_{1:t})=\sum _{z_t}p(z_{t+1},z_t|x_{1:t})=\sum _{z_t}p(z_{t+1}|z_t)p(z_t|x_{1:t})$$$$p(x_{t+1}|x_{1:t})=\sum _{z_{t+1}}p(x_{t+1},z_{t+1}|x_{1:t})=\sum _{z_{t+1}}p(x_{t+1}|z_{t+1})p(z_{t+1}|x_{1:t})$$

参考文献

【1】隐马尔可夫模型
【2】NLP —— 图模型（一）隐马尔可夫模型（Hidden Markov model，HMM）