【NLP】隐马尔可夫模型详解

前言

  本篇我们将介绍隐马尔可夫模型，隐马尔可夫模型$（HiddenMarkovModel,HMM）$是关于时序的概率模型，描述一个由一个隐藏的马尔科夫链随机生成不可测的状态随机序列，再由各个状态生成一个观测而产生观测随机序列的过程。（其中的一些名词在后面会给出定义）

  $HMM$属于生成模型，广泛应用于语音识别、自然语言处理、模式识别等领域。

  本篇将从如下几个方面进行介绍：

• 1、生成模型与判别模型
• 2、马尔可夫链
• 3、隐马尔可夫模型
• 4、概率计算算法
• 5、学习算法
• 6、预测算法

一、生成模型与判别模型

  在概率模型中，利用已知变量推测出未知变量分布的核心就是如何基于可观测变量推测出未知变量的条件分布

假定所关心的变量集合为 $Y$ ，可观测变量集合为$O$，其他变量集合为 $R$：

• 生成模型：考虑联合分布P(Y,R,O)P(Y,R,O)
• 判别模型：考虑条件分布P(Y,R|O)P(Y,R|O)

给定一组观测变量值，推断就要由 $P(Y,R,O)$或$P(Y,R|O)$得到条件概率分布$P(Y|O)$

二、马尔可夫链

我们先来看一个例子：

  假设每天的天气有三种情况：晴天、阴天、多云，根据我们的经验我们知道，今天是晴天，明天是晴天的概率比较大，是阴天或多云的概率较小，今天是阴天明天是阴天的概率较大，是多云的概率次之，是晴天的概率较小。

• 我们记：晴天（状态1）、阴天（状态2）、多云（状态3）
• 记：$P(1|2)$ 为今天是阴天，明天是晴天的概率

P(1|1)=0.8  P(1|2)=0.1  P(1|3)=0.2
P(2|1)=0.1  P(2|2)=0.7  P(2|3)=0.3
P(3|1)=0.1  P(3|2)=0.2  P(3|3)=0.5

画出状态转换图如下：

• 上图中箭头表示从一个状态到另一个状态的概率，晴天后阴天的概率为0.1，其他依次类推。
• 从上图可以看出，一个状态的转移只依赖于之前的n个状态，当n取1时就是马尔可夫假设。

由此得出马尔可夫链的定义：

  马尔可夫链是随机变量 $S_1,\dots ,S_t$ 的一个数列（状态集），这些变量的范围，即他们所有可能取值的集合，被称为“状态空间”，而 $S_t$ 的值则是在时间 $t$ 的状态。如果 $S_{t+1}$ 对于过去状态的条件概率分布仅是 $S_t$ 的一个函数，则：

$$P(S_{t+1}=x|S_1,...,S_t)=P(S_{t+1}=x|S_t)$$

其中 $x$ 为过程中的某个状态

  一个含有 $N$ 个状态的马尔可夫链有 $N^2$ 个状态转移。每一个转移的概率叫做 状态转移概率 (state transition probability)，就是从一个状态转移到另一个状态的概率。这所有的 $N^2$ 个概率可以用一个状态转移矩阵 来表示：

$$A=\left[\begin{array}{ccc}0.8& 0.1& 0.1\\ 0.1& 0.7& 0.2\\ 0.2& 0.3& 0.5\end{array}\right]$$

• 矩阵$A$的行表示当前状态，列表示当前状态下转移到其他的状态的概率
• 例如第一行表示：当前状态为1（晴天）下，下一状态为1的概率为0.8，为2的概率为0.1，为3的概率为0.1
• 且每一行的和都为1

三、隐马尔可夫模型

  在很多时候，马尔可夫过程不足以描述我们发现的问题，假设我们并不知道当前的天气情况，但是可以通过是否有风、热等进行推测，如果有风，可能是阴天；如果热，可能是晴天。由此我们将产生两个状态集，一个是可观测的状态集$O$和一个隐藏状态集$S$，我们的目的之一是借由可观测状态预测隐藏状态，这里为了简化描述，将多云这个状态去掉。重写画出下图：（例子主要是为了表达想法，所以不是那么确切）

由此得到 O = { o 1 , o 2 } ， S = { s 1 , s 2 } O = \{o_1,o_2\}，S = \{s_1,s_2\} O={o1​,o2​}，S={s1​,s2​}，$o_1$表示有风，$o_2$表示热，$s_1$表示晴天，$s_2$表示阴天。

例如，在晴天的状态下表现出有风、热的概率分别为$(0.3,0.7)$

  上面的例子中，可以观察到的 状态序列和隐藏的状态序列是概率相关的。于是我们可以将这种类型的过程建模为有一个隐藏的马尔科夫过程和一个与这个隐藏马尔科夫过程概率相关的并且可以观察到的状态集合。这就是隐马尔可夫模型。

  隐马尔可夫模型 (Hidden Markov Model，HMM) 是一种统计模型，用来描述一个含有隐含未知参数的马尔可夫过程。

  通过转移矩阵，由前面的式子：$P(S_{t+1}=x|S_1,...,S_t)=P(S_{t+1}=x|S_t)$ 我们可以知道 $P(S_t|S_{t-1})$，那么怎么表示 $P(O_t|S_t)$（观测到的状态相当于对隐藏的真实状态的一种估计）？

在$HMM$中使用另一个矩阵（观测概率矩阵）：

$$B=\left[\begin{array}{cc}0.3& 0.7\\ 0.9& 0.1\end{array}\right]$$

• 矩阵行代表隐藏状态，列代表可观测的状态，矩阵每一行概率值的和为1。
• 其中第1行第1列，$P(有风|1)=0.3$ 表示在晴天（1状态）下，有风的概率为0.3

  观测概率矩阵可视为马尔可夫模型的另一个假设，独立性假设：假设任意时刻的观测只依赖于该时刻的马尔可夫链的状态，与其它观测状态无关。

$$P(o_t|s_1,s_2,...,s_t,o_1,o_2,...,o_t)=P(o_t|s_t)式（2）$$

隐马尔可夫模型定义：

  隐藏的马尔可夫链随机生成的状态的序列称为“状态序列”；每个状态生成一个观测，而由此产生的观测的随机序列称为“观测序列”。序列的每一个位置又可以看作是一个时刻。

  隐马尔可夫模型由 初始概率分布、状态转移概率分布以及观测概率分布 确定。定义如下：

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

  隐马尔可夫模型由 初始概率向量$\pi$ 、状态转移概率矩阵$A$ 和 观测概率矩阵 $B$决定。$\pi$和$A$决定状态序列，$B$决定观测序列。因此，隐马尔可夫模型$\lambda$可以用三元符号表示，即：
$$\lambda =(A,B,\pi )$$

$A,B,\pi$称为隐马尔可夫模型的三要素。

另外，隐马尔可夫模型涉及三个基本问题：

（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(S|O)$最大的状态序列$S=(s_1,s_2,...,s_T)$，即给定观测序列，求最有可能的对应的状态序列。

四、概率计算问题

  根据前面的介绍，我们可以得到如下图模型：

状态转移概率：
$$p(q_t|q_1,q_2,...,q_{t-1},y_1,y_2,...,y_{t-1})=p(q_t|q_{t-1})$$
观测概率：
$$p(y_t|q_1,q_2,...,q_t,y_1,y_2,...,y_{t-1})=p(y_t|q_t)$$
状态转移概率矩阵：
$$A_{K\times K}=\left[\begin{array}{cc}p(q_t=1|q_{t-1}=1)& ...\\ ...& \\ ...& p(q_t=k|q_{t-1}=k)\end{array}\right]$$
观测概率矩阵：（离散状态下）
$$B_{K\times L}=\left[\begin{array}{ccc}p(y_1=1|q_t=1)& ...& p(y_1=l|q_t=1)\\ ...& & \\ ...& ...& p(q_t=l|q_{t-1}=k)\end{array}\right]$$

初始状态概率向量：

$$\pi =p(q_i)$$

马尔科夫模型$\lambda$可表示为：
$$\lambda =(A,B,\pi )$$

上面例子要求的是：
$$\begin{gathered} p(Y|\lambda )=\sum _{q_1=1}^K\sum _{q_2=1}^K\sum _{q_3=1}^{K}p(y_1,y_2,y_3,q_1,q_2,q_3) \\ =\sum _{q_1=1}^K\sum _{q_2=1}^K\sum _{q_3=1}^{K}p(y_3|y_1,y_2,q_1,q_2,q_3)p(y_1,y_2,q_1,q_2,q_3)贝叶斯公式 \\ =\sum _{q_1=1}^K\sum _{q_2=1}^K\sum _{q_3=1}^{K}p(y_3|q_3)p(q_3|y_1,y_2,q_1,q_2)p(y_1,y_2,q_1,q_2)马尔可夫过程 \\ =\sum _{q_1=1}^K\sum _{q_2=1}^K\sum _{q_3=1}^{K}p(y_3|q_3)p(q_3|q_2)p(y_1,y_2,q_1,q_2)马尔可夫过程 \end{gathered}$$
这时发现，最后一步$p(y_1,y_2,q_1,q_2)$出现了递归，因此上式:
$$=\sum _{q_1=1}^K\sum _{q_2=1}^K\sum _{q_3=1}^{K}p(y_3|q_3)p(q_3|q_2)p(y_2|q_2)p(q_2|q_1)p(q_1)$$

其中$p(q_3|q_2)$可以用$A$表示，$p(y_3|q_3)$可以用$B$表示，$p(q_1)$可以用$\pi$表示

通过上面的式子我们便可以解决三个状态的概率计算问题，如果是多个状态呢？

马尔科夫模型参数$\lambda =(A,B,\pi )$，求$p(y_1,y_2,...,y_T|\lambda )$:

$$p(y_1,y_2,...,y_T|\lambda )=\sum _{q_1=1}^K\sum _{q_2=1}^K...\sum _{q_T=1}^{K}p(y_1,y_2,...,y_T,q_1,q_2,...,q_T)$$
同样利用贝叶斯公式和马尔可夫过程得到：
$$p(Y|\lambda )=\sum _{q_1=1}^K\sum _{q_2=1}^K...\sum _{q_T=1}^K\pi (q_1)\prod _{t=2}^T{a}_{q_{t-1},q_t}{b}_{q_t}(y_t)$$
其中$a_{q_{t-1},q_t}，b_{q_t}(y_t)$为$A,B中的元素$，我们发现这个计算量特别大（$T{K}^T$）无法像上面三个状态一样通过直接计算法计算出概率。

这时候我们就可以使用 前向-后向算法 来计算观测序列概率 $P(Y|\lambda )$

4.1 前向算法

  如下图，观测序列为 $y_1,y_2,...,y_{t-1},y_t,y_{t+1},...,y_T$，状态序列为 $q_1,q_2,...,q_i,q_{i+1},...,q_T$，其中每个状态 $q_i$ 表示在此时刻状态为 $i,(i=1,2,...,N)$ ，有 $N$ 个可选状态（比如前面我们所说的晴天、阴天、多云，也就是说每个状态序列都有三个可选的状态），每个观测值 $y_t$ 表示在时刻 $t$ 的观测值，且都有 $M$ 个观测值（比如前面所说的有风、热，也就是说每个观测到的值都有两种可能），图中，红色框起来的为前向算法，后向算法与之对应。

给定隐马尔可夫模型 $\lambda$ ，定义到时刻 $t$ 部分观测序列为 $y_1,y_2,...,y_t$，且状态为 $q_i$ 的概率的前向概率记为：
$${\alpha }_t(i)=p(y_1,y_2,...,y_t,q_i=i|\lambda )$$
可以递归地求得前向概率 ${\alpha }_t(i)$ 和观测序列概率 $P(Y|\lambda )$

推导过程：

$${\alpha }_1(i)=p(y_1,q_1=i)=p(y_1|q_1=i)p(q_1)=b_i(y_1)\pi (q_1)$$

上式表示：在开始时刻的状态为 $i$， 观测值为 $y_1$的联合概率，也就等于在开始时刻的状态为 $q_1=i$的条件下观测值为 $y_1$的概率 $b_i(y_1)$（也就是前面我们提到的观测概率矩阵）乘以开始时刻每种状态的出现的概率 $\pi (q_1)$（也就是前面的初始状态概率）

$${\alpha }_2(j)=p(y_1,y_2,q_2=i)=\sum _{i=1}^{K}p(y_1,y_2,q_1=i,q_2=j)$$

这里加入了 $q_1$ ，积分出去即可，离散的积分就是求和

$$=\sum _{i=1}^{K}p(y_2|q_2=j)p(q_2=j|q_1=i)p(y_1,q_1=i)$$
$$=p(y_2|q_2=j)\sum _{i=1}^{K}p(q_2=j|q_1=i){\alpha }_1(i)$$
$$=b_j(y_2)[\sum _{i=1}^K{a}_{i,j}{\alpha }_1(i)]$$

归纳得到：
$$\begin{gathered} {\alpha }_T(j)=b_j(y_T)[\sum _{i=1}^K{a}_{i,j}{\alpha }_{T-1}(i)] \\ =p(y_1,y_2,...,y_T,q_T=j) \end{gathered}$$

此时计算复杂度为$（K^{2}T）$

$$P(Y|\lambda )=p(y_1,y_2,...,y_T)=\sum _{j=1}^K{\alpha }_T(j)$$

前向算法描述如下：

输入：隐马尔可夫模型$\lambda$，观测序列 $Y$

输出：观测序列概率$P(Y|\lambda )$

（1）初值
$${\alpha }_1(i)={\pi }_i{b}_i(y_1)，i=1,2,..,N$$
（2）递推 \quad 对$t=1,2,...,T-1$
$${\alpha }_{t+1}(i)=[\sum _{j=1}^N{\alpha }_t(j)a_{ij}]b_i(y_{t+1})，i=1,2,...,N$$
（3）终止
$$p(Y|\lambda )=\sum _{i=1}^N{\alpha }_T(i)$$

• 步骤（1）：初始化前向概率，是初始时刻的状态 $q_1=i$ 和观测 $y_1$ 的联合概率
• 步骤（2）：前向概率的递推公式，计算到时刻 $t+1$ 部分观测序列为 $y_1,y_2,....,y_t,y_{t+1}$且在时刻 $t+1$ 处于状态 $q_i$ 的前向概率

其中：

1. ${\alpha }_t(j)$是到时刻 $t$ 观测到 $y_1,y_2,...,y_t$并且在时刻 $t$ 处于状态 $q_j$ 的前向概率，那么 ${\alpha }_t(j)a_{ij}$就是到时刻 $t$ 观测到 $y_1,y_2,...,y_t$并且在时刻 $t$ 处于状态 $q_j$ 而在时刻 $t+1$ 到达状态 $q_i$ 的联合概率，对这个乘积在时刻 $t$ 的所有可能的 $N$ 个状态 $q_j$ 求和，其结果就是到时刻 $t$ 观测为 $y_1,y_2,...,y_t$并且在时刻 $t+1$ 处于状态 $q_i$ 的联合概率；
2. $[{\sum }_{j=1}^N{\alpha }_t(j)a_{ij}]b_i(y_{t+1})$恰好是到时刻 $t+1$ 观测到 $y_1,y_2,...,y_t,y_{t+1}$并且在时刻 $t+1$ 处于状态 $q_i$ 的前向概率 ${\alpha }_{t+1}(i)$

• 步骤（3）：给出 $P(Y|\lambda )$的计算公式

实例计算：

考虑盒子和球模型$\lambda =(A，B,\pi )$，状态集合 Q = { 1 , 2 , 3 } Q=\{1,2,3\} Q={1,2,3}，观测集合 V = { 红 ， 白 } V= \{红，白\} V={红，白} ，
$$A=\left[\begin{array}{ccc}0.5& 0.2& 0.3\\ 0.3& 0.5& 0.2\\ 0.2& 0.3& 0.5\end{array}\right]$$
$$B=\left[\begin{array}{cc}0.5& 0.5\\ 0.4& 0.6\\ 0.7& 0.3\end{array}\right]$$
$$\pi =(0.2,0.4,0.4{)}^T$$
设$T=3，Y=(红,白,红)$，用前向算法计算$P(Y|\lambda )$

解：
（1）计算初值
$$\begin{gathered} {\alpha }_1(1)={\pi }_1{b}_1(o_1)=0.2\ast 0.5=0.1 \\ {\alpha }_1(2)={\pi }_2{b}_2(o_1)=0.4\ast 0.4=0.16 \\ {\alpha }_1(3)={\pi }_3{b}_3(o_1)=0.4\ast 0.7=0.28 \end{gathered}$$
（2）递推计算
$$\begin{gathered} {\alpha }_2(1)=[\sum _{i=1}^3{\alpha }_1(i)a_{i1}]b_1(o_2)=[0.1\ast 0.5+0.16\ast 0.3+0.28\ast 0.2]\ast 0.5=0.077 \\ {\alpha }_2(2)=[\sum _{i=1}^3{\alpha }_1(i)a_{i2}]b_2(o_2)=0.184\ast 0.6=0.1104 \\ {\alpha }_2(3)=[\sum _{i=1}^3{\alpha }_1(i)a_{i3}]b_3(o_2)=0.202\ast 0.3=0.0606 \\ {\alpha }_3(1)=[\sum _{i=1}^3{\alpha }_2(i)a_{i1}]b_1(o_2)=[0.077\ast 0.5+0.1104\ast 0.3+0.0606\ast 0.2]\ast 0.5=0.04187 \\ {\alpha }_3(2)=[\sum _{i=1}^3{\alpha }_2(i)a_{i2}]b_2(o_2)=0.03551 \\ {\alpha }_3(3)=[\sum _{i=1}^3{\alpha }_2(i)a_{i3}]b_3(o_2)=0.05284 \end{gathered}$$
（3）终止
$$P(Y|\lambda )=\sum _{i=1}^3{\alpha }_3(i)=0.13022$$

前向算法优点：

• 前向算法高效的关键是其局部计算前向概率，然后利用路径结构将前向概率“递推”到全局，得到$P(Y|\lambda )$
• 减少计算量的原因在于每一个计算直接引用前一个时刻的计算结果，避免重复计算

4.2 后向算法

  后向算法与前向算法类似，只是一个计算前向概率一个计算后向概率。

给定隐马尔可夫模型 $\lambda$，定义在时刻 $t$ 状态为 $q_t$ 的条件下，从 $t+1$ 到 $T$ 的 部分观测序列 为 $y_{t+1},y_{t+2},...,y_T$ 的概率为后向概率，记为：
$${\beta }_t(i)=P(y_{t+1},y_{t+2},...,y_T|q_t=i,\lambda )$$
同样可以用递归的方法求得后验概率 ${\beta }_t(i)$ 和观测序列概率 $P(Y|\lambda )$

后向算法描述如下：

输入：隐马尔可夫模型$\lambda$，观测序列 $Y$

输出：观测序列概率$P(Y|\lambda )$

（1）初值
$${\beta }_T(i)=1，i=1,2,..,N$$
（2）递推 \quad 对$t=T-1,T-2,...,1$
$${\beta }_t(i)=\sum _{j=1}^N{a}_{ij}{b}_j(y_{t+1}){\beta }_{t+1}(j)，i=1,2,...,N$$
（3）终止
$$p(Y|\lambda )=\sum _{i=1}^N{\pi }_i{b}_i(y_1){\beta }_1(i)$$

• 步骤（1）：初始化后向概率，对最终时刻的所有状态 $q_i$ 规定 ${\beta }_T(i)=1$；
• 步骤（2）：后向概率的递推公式，计算在时刻 $t$ 状态为 $q_i$ 条件下时刻 $t+1$ 之后的部分观测序列为 $y_{t+1},y_{t+2},...,y_T$ 的后向概率 ${\beta }_t(i)$，只需要考虑在时刻 $t+1$ 所有可能的 $N$ 个状态 $q_j$ 的转移概率（即 $a_{ij}$），以及在此状态下的观测 $y_{t+1}$ 的观测概率（即 $b_j(y_{t+1})$），然后考虑状态 $q_j$ 之后的观测序列的后向概率（即 ${\beta }_{t+1}(j)$）
• 步骤（3）：给出 $P(Y|\lambda )$的计算公式

利用前向概率和后向概率的定义可以将观测序列概率 $P(Y|\lambda )$ 统一写为：

$$P(Y|\lambda )=\sum _{i=1}^N\sum _{j=1}^N{\alpha }_t(i)b_j(y_{t+1}){\beta }_{t+1}(j),t=1,2,...,T-1$$

五、学习问题

  上一节我们解决了在已知隐马尔可夫模型 $\lambda$的条件下 $P(Y|\lambda )$ 概率计算问题，那么学习问题就是根据观测序列来估计模型 $\lambda$ 参数，使得在新的模型参数 $\lambda$ 下观测序列概率 $P(Y|\lambda )$ 最大，即用极大似然估计的方法进行参数估计。

隐马尔可夫模型定义：$\lambda =(A,B,\pi )$

其中：

$A$是状态转移概率矩阵：
$$A=[a_{ij}{]}_{N\ast N}$$
$B$是观测概率矩阵
$$B=[b_j(k){]}_{N\ast M}$$
$\pi$是初始状态概率向量：
$$\pi =({\pi }_i)$$

5.1 监督学习方法

\quad\quad 假设已给定训练数据包含$S$个长度相同的观测序列和对应的状态序列 { ( 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​)}，那么可以利用极大似然估计法来估计隐马尔可夫模型的参数：

（1）转移概率$a_{ij}$的估计：
$$\widehat{a_{ij}}=\frac{A_{ij}}{{\sum }_{j=1}^N{A}_{ij}}，i=1,2,...,N；j=1,2,...,N$$
其中$A_{ij}$表示样本中时刻$t$处于状态$i$时刻$t+1$转移到状态$j$的频数；
（2）观测概率$b_j(k)$的估计：
$$\widehat{b_j(k)}=\frac{B_{jk}}{{\sum }_{k=1}^M{B}_{jk}}，j=1,2,...,N；k=1,2,...,M$$
其中$B_{jk}$表示样本中状态为$j$并且观测为$k$的频数；
（3）初始状态概率${\pi }_i$的估计：
$$\widehat{{\pi }_i}为S个样本中初始状态为q_i的频数$$

然而，监督学习所需的训练数据集需要人工标注，导致代价很高，我们得到的大部分数据是只有观测序列的，我们就会利用非监督学习的方法

5.2 非监督学习方法-Baum-Welch算法（由EM算法实现）

\quad\quad 假定给定训练数据集只包含$S$个长度为$T$的观测序列 { O 1 , O 2 , . . . , O S } \{O_1,O_2,...,O_S\} {O1​,O2​,...,OS​}而没用对应的状态序列，目标是学习隐马尔可夫模型$\lambda =(A,B,\pi )$。

\quad\quad 我们将观测序列数据看着观测数据$O$，状态序列数据看作不可观测的隐数据$I$，那么隐马尔可夫模型事实上就是一个含有隐变量的概率模型：
$$P(O|\lambda )=\sum _{I}P(O|I,\lambda )P(I|\lambda )$$
EM算法一文中，我们提到$EM$算法一般用于含有隐变量的概率模型参数的极大似然估计，或极大后验概率估计。

$EM$算法的迭代公式：

$${\lambda }^{(i+1)}=arg\underset{\lambda }{\max }Q(\lambda ,{\lambda }^{(i)})$$
其中：$$Q(\lambda ,{\lambda }^{(i)})=\sum _{I}P(O,I|\lambda )P(I|O,{\lambda }^{(i)})$$
又由贝叶斯公式：$$P(O,I|{\lambda }^{(i)})=P(I|O,{\lambda }^{(i)})\ast P(O|{\lambda }^{(i)})$$
$$P(I|O,{\lambda }^{(i)})=\frac{P(O,I|{\lambda }^{(i)})}{P(O|{\lambda }^{(i)})}$$
分母项为常数，求$Q$ 函数极大化可以不考虑，所以：
$$Q(\lambda ,{\lambda }^{(i)})\propto \sum _{I}P(O,I|\lambda )P(O,I|{\lambda }^{(i)})$$
$O$为观测变量数据，$I$为隐变量数据，联合分布$P(O,I|\lambda )$，条件分布$P(I|O,\lambda )$

（1）确定完全数据的对数似然函数

观测数据$O=(o_1,o_2,...,o_T)$，隐数据$I=(i_1,i_2,...,i_T)$，

完全数据是$(O,I)=(o_1,o_2,...,o_T,i_1,i_2,...,i_T)$，完全数据的对数似然函数为：
$$logP(O,I|\lambda )$$
（2）EM算法的E步：求Q函数$Q(\lambda ,\overline{\lambda })$
$$Q(\lambda ,\overline{\lambda })\propto \sum _{I}log(P(O,I|\lambda ))P(O,I|\overline{\lambda })$$
其中，$\overline{\lambda }$是隐马尔可夫模型参数的当前估计值，$\lambda$是使$Q$函数极大化的模型参数。

且：
$$P(O,I|\lambda )={\pi }_{i_1}{b}_{i_1}(o_1)a_{i_1{i}_2}{b}_{i_2}(o_2)\cdot \cdot \cdot a_{i_{T-1}{i}_T}{b}_{i_T}(o_T)$$
于是函数$Q(\lambda ,\overline{\lambda })$可以写成：
$$Q(\lambda ,\overline{\lambda })=\sum _{I}log{\pi }_{i_1}P(O,I|\overline{\lambda })+\sum _I(\sum _{t=1}^{T-1}log{a}_{i_t{i}_{t+1}})P(O,I|\overline{\lambda })+\sum _I(\sum _{t=1}^{T}log{b}_{i_t}(o_t))P(O,I|\overline{\lambda })$$
式中求和都是对所有训练数据的序列总长度$T$进行的

（3）EM算法的M步：极大化Q函数求模型参数

由于Q函数中有3项相加，所以可以对每一项分别求极大化：

• 第一项：
  $$\sum _{I}log{\pi }_{i_1}P(O,I|\overline{\lambda })=\sum _{i=1}^N{\pi }_{i}P(O,i_1=i|\overline{\lambda })$$
  其中，${\pi }_i$满足约束条件${\sum }_{i=1}^N{\pi }_i=1$，可以用拉格朗日乘子法

写出拉格朗日函数：
$$L{M}_1=\sum _{i=1}^N{\pi }_{i}P(O,i_1=i|\overline{\lambda })+\gamma (\sum _{i=1}^N{\pi }_i-1)$$
N是可能的状态数
对其求偏导数，并且令为0
$$\frac{\partial L{M}_1}{\partial {\pi }_i}=\frac{P(O,i_1=i|\overline{\lambda })}{{\pi }_i}+\gamma =0$$
$$\frac{\partial L{M}_1}{\partial \gamma }=\sum _{i=1}^N{\pi }_i-1=0$$
所以：
$$P(O,i_1=i|\overline{\lambda })=-\gamma {\pi }_i$$
$$\sum _{i=1}^{N}P(O,i_1=i|\overline{\lambda })=-\gamma \sum _{i=1}^N{\pi }_i=-\gamma$$
$${\pi }_i=\frac{P(O,i_1=i|\overline{\lambda })}{-\gamma }$$
$${\pi }_i=\frac{P(O,i_1=i|\overline{\lambda })}{{\sum }_{i=1}^{N}P(O,i_1=i|\overline{\lambda })}$$

• 第二项：
  $$\sum _I\left(\sum _{t=1}^{T-1}log{a}_{i_t{i}_{t+1}}\right)P(O,I|\overline{\lambda })=\sum _{i=1}^N\sum _{j=1}^N\sum _{t=1}^{T-1}log{a}_{ij}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })$$

其中，${\sum }_{j=1}^N{a}_{ij}=1$，同样利用拉格朗日乘子法

写出拉格朗日函数：
$$L{M}_2=\sum _{i=1}^N\sum _{j=1}^N\sum _{t=1}^{T-1}log{a}_{ij}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })+\sum _{i=1}^N{\gamma }_i(\sum _{j=1}^N{a}_{ij}-1)$$
N是可能的状态数
对其求偏导数，且令为0
$$\frac{\partial L{M}_2}{\partial a_{ij}}=\frac{{\sum }_{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })}{a_{ij}}+\sum _{i=1}^N{\gamma }_i=0$$
$$\frac{\partial L{M}_2}{\partial {\gamma }_i}=\sum _{j=1}^N{a}_{ij}-1=0$$
$$\sum _{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })=-a_{ij}\sum _{i=1}^N{\gamma }_i$$
$$a_{ij}=\frac{{\sum }_{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })}{-{\sum }_{i=1}^N{\gamma }_i}$$
$$\sum _{j=1}^N\sum _{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })=\sum _{j=1}^N-a_{ij}\sum _{i=1}^N{\gamma }_i=-\sum _{i=1}^N{\gamma }_i\sum _{j=1}^N{a}_{ij}=-\sum _{i=1}^N{\gamma }_i$$
$$a_{ij}=\frac{{\sum }_{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })}{{\sum }_{j=1}^N({\sum }_{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda }))}$$
$$=\frac{{\sum }_{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })}{{\sum }_{t=1}^{T-1}P(O,i_t=i,|\overline{\lambda })}$$

• 第三项：
  $$\sum _I(\sum _{t=1}^{T}log{b}_{i_t}(o_t))P(O,I|\overline{\lambda })=\sum _{j=1}^N\sum _{t=1}^{T}log{b}_j(o_t)P(O,i_t=j|\overline{\lambda })$$

其中，${\sum }_{k=1}^M{b}_j(k)=1$，同样使用拉格朗日乘子法

写出拉格朗日函数：
$$L{M}_3=\sum _{j=1}^N\sum _{t=1}^{T}log{b}_j(o_t)P(O,i_t=j|\overline{\lambda })+\gamma (\sum _{k=1}^M{b}_j(k)-1)$$
N是可能的状态数，M是可能的观测数
对其求偏导数，且令为0：
注意：只有在$o_t=v_k$时$b_j(k)$对$b_j(o_t)$的偏导数才不为0,以$I(o_t=v_k)$表示
$$\frac{\partial L{M}_3}{\partial b_j(o_t)}=\frac{{\sum }_{t=1}^{T}P(O,i_t=j|\overline{\lambda })}{b_j(o_t)}+\gamma =0$$
$$\frac{\partial L{M}_3}{\partial \gamma }=\sum _{j=1}^M{b}_j(o_t)-1=0$$
$$\sum _{t=1}^{T}P(O,i_t=j|\overline{\lambda })=-b_j(o_t)\gamma$$
$$b_j(o_t)=\frac{{\sum }_{t=1}^{T}P(O,i_t=j|\overline{\lambda })}{-\gamma }$$
$$\sum _{j=1}^N\sum _{t=1}^{T}P(O,i_t=j|\overline{\lambda })=-\gamma \sum _{j=1}^N{b}_j(k)=-\gamma$$
$$b_j(k)=\frac{{\sum }_{t=1}^{T}P(O,i_t=j|\overline{\lambda })I(o_t=v_k)}{{\sum }_{t=1}^{T}P(O,i_t=j|\overline{\lambda })}$$

Baum-Welch算法描述如下：

输入：观测数据$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}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })}{{\sum }_{t=1}^{T-1}P(O,i_t=i,|\overline{\lambda })}$$

$$b_j(k{)}^{(n+1)}=\frac{{\sum }_{t=1}^{T}P(O,i_t=j|\overline{\lambda })I(o_t=v_k)}{{\sum }_{t=1}^{T}P(O,i_t=j|\overline{\lambda })}$$

$${\pi }_i^{(n+1)}=\frac{P(O,i_1=i|\overline{\lambda })}{{\sum }_{i=1}^{N}P(O,i_1=i|\overline{\lambda })}$$

右端各值按模型${\lambda }^{(n)}=(A^{(n)},B^{(n)},{\pi }^{(n)})$计算

（3）终止。得到模型参数${\lambda }^{(n+1)}=(A^{(n+1)},B^{(n+1)},{\pi }^{(n+1)})$

六、预测问题

  预测问题：已知模型$\lambda =(A,B,\pi )$和观测序列$O=(o_1,o_2,...,o_T)$，求对给定观测序列条件概率$P(S|O)$最大的状态序列$S=(s_1,s_2,...,s_T)$，即给定观测序列，求最有可能的对应的状态序列。

6.1 近似算法

思想描述：

  在每个时刻 $t$ 选择在该时刻最有可能出现的状态 $i_t^{\ast }$ ，从而得到一个状态序列 $I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$，这便是我们预测的结果。

• 在给定的隐马尔可夫模型 $\lambda$ 和观测序列 $O$下，在时刻 $t$ 处于状态 $q_i$ 的概率 ${\gamma }_t(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$ 最有可能的状态 $i_t^{\ast }$，表示为：

$$i_t^{\ast }=\arg \underset{1\le i\le N}{\max }[{\gamma }_t(i)],t=1,2,...,T$$

从而得到状态序列 $I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$。

近似算法计算简单，但是不能保证预测的状态序列整体上是最优的，因为预测的状态序列可能有实际不发生的部分。事实上，该方法得到的状态序列中有可能存在转移概率为0的相邻状态，即对某些 $i,j,a_{ij}=0$。

6.2 维特比算法

  维特比算法实际上是用 动态规划 解隐马尔可夫模型预测问题，即用 动态规划 求概率最大路径（这是最优路径），一条路径对应一个预测的状态序列。

动态规划思想描述如下：

  如果最优路径在时刻 $t$ 通过结点 $i_t^{\ast }$，那么这一路径从结点 $i_t^{\ast }$ 到终点 $i_T^{\ast }$的部分路径，对于从 $i_t^{\ast }$ 到$i_T^{\ast }$ 的所有可能的部分路径来说，也必须是最优的。因为如果不是最优，那么必定存在另一条比原来的路径更优的路径。

  根据上面的思想，我们只需要从时刻 $t=1$ 开始，递归地计算在时刻 $t$ 状态为 $i$ 的各条部分路径的最大概率，直到时刻 $t=T$ 状态为 $i$ 的各条路径的最大概率。时刻 $t=T$ 的最大概率即为最优路径的概率 $P^{\ast }$，最优路径的终结点 $i_T^{\ast }$ 也就得到了。之后，为了找出最优路径的各个节点，从终结点 $i_T^{\ast }$ 开始，由后向前逐步求得结点 $i_{T-1}^{\ast },i_{T-2}^{\ast },...,i_1^{\ast }$，最终得到最优路径 $I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$，这便是维特比算法。

定义变量 $\delta ,\psi$

• ${\delta }_t(i)$：表示在时刻 $t$ 状态为 $i$ 的所有单个路径 $(i_1,i_2,...,i_t)$ 中的概率最大值

$${\delta }_t(i)=\underset{i_1,i_2,...,i_{t-1}}{\max }P(i_t=i,i_{t-1},...,t_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+1}=i,i_t,...,t_1,o_{t+1},...,o_1|\lambda )$$
$$=\underset{1\le j\le N}{\max }[{\delta }_t(j)a_{ji}]b_t(o_{t+1}),i=1,2,...,N;t=1,2,..,T-1$$

• ${\psi }_t(i)$：表示在时刻 $t$ 状态为 $i$ 的所有单个路径 $(i_1,i_2,...,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,...,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),i=1,2,...,N$$
$${\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 }=\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）最优路径回溯，对 $t=T-1,T-2,...,1$

$$i_t^{\ast }={\psi }_{t+1}(i_{t+1}^{\ast })$$

求得最优路径 $I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$

实例计算：

定义隐马尔可夫模型 $\lambda =(A,B,\pi )$
$$A=\left[\begin{array}{ccc}0.5& 0.2& 0.3\\ 0.3& 0.5& 0.2\\ 0.2& 0.3& 0.5\end{array}\right]$$
$$B=\left[\begin{array}{cc}0.5& 0.5\\ 0.4& 0.6\\ 0.7& 0.3\end{array}\right]$$
$$\pi =(0.2,0.4,0.4{)}^T$$
已知观测序列 $O=(红，白，红)$，试求最优状态序列，即最优路径 $I^{\ast }=(i_1^{\ast },i_2^{\ast },i_3^{\ast })$

（1）初始化

在 $t=1$ 时，对每一个状态 $i$，$i=1,2,3$，求状态为 $i$ 观测 $o_1$ 为红的概率，记此概率为 ${\delta }_1(i)$，则：

$${\delta }_1(i)={\pi }_i{b}_i(o_1)={\pi }_i{b}_i(红),1=1,2,3$$

代入实际数据：

$${\delta }_1(1)=0.2\times 0.5=0.10,{\delta }_1(2)=0.4\times 0.4=0.16,{\delta }_1(3)=0.4\times 0.7=0.28$$

记 ${\psi }_1(i)=0,i=1,2,3$

（2）递推，对 $t=2,3,...,T$

在 $t=2$时，对每个状态 $i$，$i=1,2,3$，求在 $t=1$时状态为 $j$ 观测为红并在 $t=2$ 时状态为 $i$ 观测 $o_2$ 为白的路径的最大概率，记为 ${\delta }_2(i)$，则：

$${\delta }_2(i)=\underset{1\le j\le 3}{\max }[{\delta }_1(j)a_{ji}]b_i(o_2)$$

同时，对每个状态 $i$，$i=1,2,3$，记录概率最大路径的前一个状态 $j$：

$${\psi }_2(i)=arg\underset{1\le j\le 3}{\max }[{\delta }_1(j)a_{ji}],i=1,2,3$$

计算：

$${\delta }_2(1)=\underset{1\le j\le 3}{\max }[{\delta }_1(j)a_{j1}]b_1(o_2)$$
= max ⁡ j { 0.10 × 0.5 , 0.16 × 0.3 , 0.28 × 0.2 } × 0.5 = 0.028 =\max_{j}\{0.10 \times0.5,0.16\times0.3,0.28\times0.2\}\times 0.5=0.028 =jmax​{0.10×0.5,0.16×0.3,0.28×0.2}×0.5=0.028

$${\psi }_2(1)=3$$
$${\delta }_2(2)=0.0504,{\psi }_2(2)=3$$
$${\delta }_2(2)=0.042,{\psi }_2(2)=3$$

同样，在 $t=3$ 时，

$${\delta }_3(i)=\underset{1\le j\le 3}{\max }[{\delta }_2(j)a_{ji}]b_i(o_3)$$
$${\psi }_3(i)=arg\underset{1\le j\le 3}{\max }[{\delta }_2(j)a_{ji}]$$
$${\delta }_3(1)=0.00756,{\psi }_3(1)=3$$
$${\delta }_3(2)=0.01008,{\psi }_3(2)=3$$
$${\delta }_3(3)=0.0147,{\psi }_3(3)=3$$

（3）终止

以 $P^{\ast }$表示最优路径的概率，则

$$P^{\ast }=\underset{1\le i\le 3}{\max }{\delta }_3(i)=0.0147$$

最优路径的终点是 $i_3^{\ast }$：

$$i_3^{\ast }=arg\underset{i}{\max }[{\delta }_3(i)]=3$$

（4）最优路径回溯，对 $t=T-1,T-2,...,1$

由最优路径的终点 $i_3^{\ast }$，逆向找到 $i_2^{\ast },i_1^{\ast }$：

• 在 $t=2$ 时，$i_2^{\ast }={\psi }_3(i_3^{\ast })={\psi }_3(3)=3$
• 在 $t=1$ 时，$i_1^{\ast }={\psi }_2(i_2^{\ast })={\psi }_2(3)=3$

于是求得最优路径，即最优状态序列 $I^{\ast }=(i_1^{\ast },i_2^{\ast },i_3^{\ast })=(3,3,3)$

七、总结

  至此，关于隐马尔可夫模型的理论就介绍完了，那么隐马尔可夫模型到底有什么应用呢？隐含马尔可夫模型最早成功的使用场景是语音识别，后来陆续成功的应用在机器翻译、词性标注、拼写错误、手写体识别、图像处理、基因序列分析等很多it领域，目前也被用于股票预测和投资。总而言之，隐马尔可夫模型可以解决那些根据已知一个数据序列求另一个数据序列的问题。

  隐马尔可夫模型主要就是三个基本问题，而概率计算问题是为了学习过程中的概率计算并使其最大化，最终的目标便是根据学习问题得到最优模型去预测最优的状态序列。

• 前向-后向算法用来解决概率计算问题
• Baum-Welch算法用来解决学习问题
• 维特比算法用来解决预测问题

本篇主要内容参考《统计学习方法》
下篇博客我们将用python来实现，三个基本问题的算法实现，链接