第十九课.隐马尔科夫模型

隐马尔科夫模型的结构

马尔科夫链与隐马尔科夫模型

隐马尔科夫模型的全称为 Hidden Markov Model（HMM），这是一种统计模型，广泛应用于语音识别，词性自动标注等问题。马尔科夫链（回顾第十四课-马尔科夫链）与HMM的差距体现于Hidden，在这个模型中，首先由一个隐藏的马尔科夫链随机生成一个状态随机序列，再由状态随机序列中的每一个状态对应生成各自的观测，由这些观测构成一个观测随机序列。

实例

为了清晰说明问题，现在通过举例介绍HMM，第一个例子叫做：盒子摸球实验。有3个盒子，编号为1号，2号，3号，每个盒子里都装着个数不等的黑球和白球，具体情况如下描述：

• 1号盒子：黑球2个，白球8个；
• 2号盒子：黑球6个，白球4个；
• 3号盒子：黑球4个，白球6个；

实验过程如下：

每次先随机选择一个盒子，然后从随机选中的盒子中随机摸出一个球，并记录球的颜色，最后把球放回盒子。下一次再随机选择一个盒子，重复以上操作。

在实验过程中，我们只能在每次摸出球之后看到被摸出球的颜色，但不知道盒子的编号。因此，随着实验的进行，会依次出现不同编号的盒子，这个盒子的序列就是隐藏的状态序列$Z$。由于我们只能观测到球的颜色，所以球的颜色序列为观测序列$X$。

现在，我们假设每次盒子随机出现的过程是马尔科夫过程，状态集合为$Q=${盒子1，盒子2，盒子3}，$N=3$。同时，第一次各个盒子出现所满足的概率分布如下：

• 1号盒子出现的概率：0.3；
• 2号盒子出现的概率：0.5；
• 3号盒子出现的概率：0.2；

用$\pi$表示初始状态的概率向量为：$\pi =(0.3,0.5,0.2{)}^T$，这就是状态的初始概率分布。

同时假设各个盒子之间相互转换的概率转移图如下：

从上图可以得到关于盒子的状态转移矩阵，记作$A$：

对于实验中隐藏的状态序列，即三个盒子随机出现的过程如上所述，下面是从盒子中摸球的过程，比如在1号盒子中：黑球2个，白球8个，摸球后并放回。这是最简单的古典概型，从1号盒子中摸出黑球的概率为0.2，摸出白球的概率为0.8，即观测概率，也叫输出概率，它是从特定的隐含状态中生成指定观测的概率。

同样的，我们可以把2号盒子和3号盒子的观测概率求出，得到另一个矩阵（观测概率矩阵$B$）：

观测集合为$V=${黑球，白球}，$M=2$。

通过以上描述，重复7次上述过程，得到两个序列：

• 一个是长度为7的隐藏状态序列：$I=${2号，2号，1号，3号，1号，2号，3号}；注意这个序列是实际存在的，只是我们不能观测到；
• 另一个是长度为7的观测序列：$O=${黑球，黑球，白球，黑球，白球，白球，黑球}，这是我们可以观测到的序列；

在本次盒子摸球实验中，可以得到以下序列关系：

上图直观表达了HMM的两个核心：

• 1.隐：即状态序列，盒子序列是我们看不到的，即隐藏的；
• 2.马尔科夫：指的是整个隐藏状态序列，隐状态之间的转换是一个马尔科夫过程，即隐状态之间的转换遵循指定的概率。

将以上信息反映到概率转移图中，可以看到隐马尔科夫的全部信息：

另一个实例更加贴近生活，一个婴儿有两种状态：{饿了，困了}，但是婴儿不会说话，所以我们无法观测到这两个状态。我们只能从他的行为去推断隐状态，行为有三种：{哭闹，无精打采，爬来爬去}，这就是观测集。

隐状态之间的转移是一个马尔科夫过程，而从状态表现出的行为也会符合一定的概率分布：

这其中的状态集合$Q$，观测集合$V$，状态转移概率矩阵$A$，观测概率矩阵$B$，可以直观从上图得到。

HMM的要素

通过上面的两个实例，用形式化的语言总结HMM的要素：

HMM是一个时序模型，首先由一个隐藏的马尔科夫链按照其设定的状态转移概率，生成一个状态序列，但是这个状态序列是不能观测到的，然后再由每个状态按照观测概率（又称为输出概率）生成对应的观测，由此构成可观测的观测序列。

HMM中所有的隐含状态构成状态集合 Q = { q 1 , q 2 , . . . , q N } Q=\left\{q_{1},q_{2},...,q_{N}\right\} Q={q1​,q2​,...,qN​}，状态个数为$N$，所有的观测构成观测集合 V = { v 1 , v 2 , . . . , v M } V=\left\{v_{1},v_{2},...,v_{M}\right\} V={v1​,v2​,...,vM​}，观测的个数为$M$。经过一段时间$T$之后，生成长度为$T$的状态序列 I = { i 1 , i 2 , . . . , i T } I=\left\{i_{1},i_{2},...,i_{T}\right\} I={i1​,i2​,...,iT​}和对应的观测序列 O = { o 1 , o 2 , . . . , o T } O=\left\{o_{1},o_{2},...,o_{T}\right\} O={o1​,o2​,...,oT​}。

以上是HMM的外在表征，推动HMM运行的是三个要素：状态转移概率矩阵$A$，观测概率矩阵$B$（输出概率矩阵），初始隐含状态的概率向量$\pi$，简写为 λ = { A , B , π } \lambda=\left\{A,B,\pi\right\} λ={A,B,π}。

初始概率向量$\pi =({\pi }_1,{\pi }_2,...,{\pi }_N)$，其中${\pi }_n$表示隐含状态序列中第1个状态为$q_n$的概率，即：$${\pi }_n=P(i_1=q_n)$$
状态转移概率矩阵$A$本质上就是马尔科夫链的转移概率矩阵，$A$是一个$N\times N$的方阵：

其中，$a_{ij}$表示从隐含状态$i$转移到隐含状态$j$的概率，即$a_{ij}=P(i_{t+1}=q_j|i_t=q_i)$；

观测概率矩阵$B$是一个$N\times M$的矩阵：

其中，$b_{ij}$指的是在某时刻$t$，隐含状态为$q_i$的情况下，生成观测$v_j$的概率，即$b_{ij}=P(o_1=v_j|i_t=q_i)$

模型的性质

$t$时刻的隐状态只与前一时刻的隐状态相关：

在隐马尔科夫模型的三要素中，状态转移概率矩阵$A$和初始状态概率向量$\pi$完全确定了隐藏的马尔科夫链，其性质服从马尔科夫性，可以描述为：$$P(i_t|i_{t-1},o_{t-1},i_{t-2},o_{t-2},...,i_1,o_1)=P(i_t|i_{t-1})$$
$t$时刻的观测只与该时刻的隐状态相关：

该性质描述为：$$P(o_t|i_t,i_{t-1},o_{t-1},i_{t-2},o_{t-2},...,i_1,o_1)=P(o_t|i_t)$$

推理问题：HMM的状态解码

隐状态解码问题

前面描述了HMM的概率估计，现在讨论HMM的状态解码。解码，即给定已知的观测序列，求它最有可能对应的状态序列。

即已知模型$\lambda =(A,B,\pi )$，观测序列$O=(o_1,o_2,...,o_T)$，求使得条件概率$P(I|O)$最大的隐状态序列$I=(i_1,i_2,...,i_T)$；

先不考虑输出的观测，只考虑隐状态的转移，目标是找到路径概率最大的一条状态序列。对照下图，图中展现了HMM的状态序列，其中一共包含5种隐状态，状态的序列长度为7，图中横轴是时间轴，纵轴是状态轴：

最大路径概率与维特比算法

从整个隐状态序列的最后面向前看，在第7个时间点，即最后一个时间点，要考虑状态序列的最后一个状态是状态$[1,2,3,4,5]$中的哪一个，即路径以哪个结束可以获得最大概率。

首先关注最后一个时间点（时间点7），问题变成如果状态转移的路径以结束于状态$k$的路径概率最大，概率该如何表示？而时间点6也可以取5种状态的任一种，最终取决于从哪个状态转移到时间点7的状态$k$有最大概率。如下所示：

我们令$X_{6i}$表示到达时间点6时，此时状态为$i$的最大路径概率，状态$i$可以取 { 1 , 2 , 3 , 4 , 5 } \left\{1,2,3,4,5\right\} {1,2,3,4,5}中的任意一个，那么实际上就有$X_{61},X_{62},X_{63},X_{64},X_{65}$五个不同的值。在上图中，对应了虚线框内五种颜色示意的到达时间节点6的五条路径，他们分别是在时间节点6时到达对应状态$i$的最大概率路径。

那么用它乘以对应的状态转移概率，即$X_{6i}{A}_{ik}$，就推进到第7个时间节点。首先固定第7个时间节点结束的状态，比如选择状态1，则可以分别求出从第6个时间节点的状态1，状态2，状态3，状态4，状态5转移到时间点7状态1的概率：$$X_{61}{A}_{11},X_{62}{A}_{21},X_{63}{A}_{31},X_{64}{A}_{41},X_{65}{A}_{51}$$
我们计算出这5个值，取最大值为结束于状态1的最大路径概率。

同理，我们分别假设，令时间点7的结束状态为2,3,4,5，按照上述步骤，取每个结束状态的最大概率路径。然后在这些结束状态的最大路径中选最大的一个，作为最后真正的最大路径。

对于$X_{6i}$，它依赖于时间点5的$X_{5j}$，一直追溯到最早的时间点1为尽头。因此，现在从时间点1出发，描述以上过程：

• 首先在时间点1，我们计算出各种状态出现的概率，由于只有一个节点，因此这个概率值就是此时的最大路径概率$X_{1i}$；
• 然后我们再向前前进到时间点2，对于每一个状态$j$，我们通过利用时间点1的每一个状态的最大路径乘以转移概率，得到5个到达时间点2，状态$j$的路径概率，取最大的一个作为此时的最大路径概率$X_{2j}$，即：$$X_{2j}=max(X_{1i}{A}_{ij})$$其中，$i$遍历1,2,3,4,5每一个状态。同时，我们把最大路径的序列记录下来。
• 以此类推，我们基于时间点$t$中，各状态的最大路径概率，就可以前向计算出$t+1$中各状态的最大路径概率，直到最后一个时间节点$T$，我们得到时间节点$T$的所有状态的最大概率路径$X_{Ti}$，取最大值：$$max(X_{Ti})$$其中，$i$遍历1,2,3,4,5每一个状态。这就是我们要求的最大路径概率，以及结束的状态，然后依据记录的路径得到目标的隐状态序列。

以上即为求最大概率路径的方法，即维特比算法（用于计算最大或最小路径的动态规划算法）。

使用维特比算法解码

下面，我们将维特比算法用于HMM的解码，维特比算法中最大概率路径对应着HMM要求的状态序列，不过现在不光要考虑隐状态的转移概率，还有观测输出概率。

即要寻找一条隐状态$(i_1,i_2,...,i_t)$，用它生成指定的观测序列$(o_1,o_2,...,o_t)$，使得观测序列的存在概率最大。描述为：$${\delta }_t(i)=maxP(i_t=i,i_{t-1},...,i_1,o_t,...,o_1)$$其中，$i=1,2,3,...,N$，表示在时刻$t$，结束于隐状态$i$，同时满足观测序列$(o_1,o_2,...,o_t)$的最大路径概率，这个表达式同时考虑了隐状态和输出观测，是一个联合概率，其中$i$取$1,2,3,...,N$，最终我们取最大的$P$。

以此类推：$${\delta }_{t+1}(i)=maxP(i_{t+1}=i,i_t,...,i_1,o_{t+1},...,o_1)$$其中，$i=1,2,3,...,N$；

对于${\delta }_{t+1}(i)$和${\delta }_t(i)$的递推关系，我们根据HMM进行状态转移，再进行观测输出，所以有：$${\delta }_{t+1}(i)=max({\delta }_t(j)a_{ji}{b}_{i{o}_{t+1}})$$证明如下：$$max({\delta }_t(j)a_{ji}{b}_{i{o}_{t+1}})=max({\delta }_t(j)a_{ji})b_{i{o}_{t+1}}$$$$=max[(maxP(i_t=j,i_{t-1},...,i_1,o_t,...,o_1))P(i_{t+1}=i|i_t=j)P(o_{t+1}|i_{t+1}=i)]$$$$=max[P(i_{t+1}=i,i_t=j,i_{t-1},...,i_1,o_t,...,o_1)]P(o_{t+1}|i_{t+1}=i)$$$$=max[P(i_{t+1}=i,i_t,i_{t-1},...,i_1,o_{t+1},o_t,...,o_1))]={\delta }_{t+1}(i)$$
现在回顾递推公式：${\delta }_{t+1}(i)=max({\delta }_t(j)a_{ji}{b}_{i{o}_{t+1}})$，与上一部分维特比算法求路径概率的例子相比，多了一个观测输出概率$b_{i{o}_{t+1}}$，这个值是已知的，只是比单纯的路径概率多进行了一次观测输出概率的乘积运算。

按照维特比算法从${\delta }_1(i)$起步，一步一步按照上述方法推导到${\delta }_T(i)$，求得最大的概率，同时在递推过程中记录序列。

实例演示

现在用盒子摸球模型演示计算过程，隐状态集合为$Q=${盒子1，盒子2，盒子3}，$N=3$，初始概率$\pi =(0.3,0.5,0.2{)}^T$，状态概率矩阵为$A$：

观测集合$V=${黑球，白球}，$M=2$，观测概率矩阵为$B$：

并且，已知观测序列为 {黑球，白球，黑球}；

首先，在初始化时，对于在$t=1$时刻，隐状态为$i$，观测$o_1$为黑球的概率${\delta }_1(i)$，有：$${\delta }_1(1)={\pi }_1{b}_{1o_1}=0.3\times 0.2=0.06$$$${\delta }_1(2)={\pi }_2{b}_{2o_1}=0.5\times 0.6=0.30$$$${\delta }_1(3)={\pi }_3{b}_{3o_1}=0.2\times 0.4=0.08$$递推到$t=2$时，表示在$t=1$时隐状态为$j$号盒子，输出观测$o_1$为黑球，$t=2$时隐状态为$i$号盒子，输出观测$o_2$为白球的最大概率为：
${\delta }_2(i)=max[{\delta }_1(j)a_{ji}]b_{i{o}_2}$

所以有：

${\delta }_2(1)=max[{\delta }_1(j)a_{j1}]b_{1o_2}=max[{\delta }_1(1)a_{11},{\delta }_1(2)a_{21},{\delta }_1(3)a_{31}]b_{1o_2}$
$=max[0.06\times 0.4,0.3\times 0.3,0.08\times 0.2]0.8=max[0.024,0.09,0.016]0.8=0.072$

同时记录$t=1$时，$j=2$；

${\delta }_2(2)=max[{\delta }_1(j)a_{j2}]b_{2o_2}=max[{\delta }_1(1)a_{12},{\delta }_1(2)a_{22},{\delta }_1(3)a_{32}]b_{2o_2}$
$=max[0.024,0.06,0.048]0.4=0.024$

其中$b_{2o_2}=0.4$，同时记录$t=1$时，$j=2$；

${\delta }_2(3)=max[{\delta }_1(j)a_{j3}]b_{3o_2}=max[{\delta }_1(1)a_{13},{\delta }_1(2)a_{23},{\delta }_1(3)a_{33}]b_{3o_2}$
$=max[0.012,0.15,0.016]0.6=0.09$

其中$b_{3o_2}=0.6$，同时记录$t=1$时，$j=2$；

我们把上述过程的结果可视化：

现在递推到$t=3$：

结果可视化为：

现在，我们从最后一个时间节点3开始，选择状态2，它的路径概率最大；

通过追溯，得到状态序列为：{盒子2，盒子3，盒子2}

基于Python对盒子摸球实验进行状态解码

下面基于Python实现了上述过程的手工计算，HMM使用hmmlearn实现，相关文档为：hmmlearn，其提供了三种模型：

名称	简介	观测类型
hmm.GaussianHMM	观测服从高斯分布的HMM	连续
hmm.GMMHMM	观测服从混合高斯分布的HMM	连续
hmm.MultinomialHMM	观测为离散型的HMM	离散

---

隐马尔科夫模型，要求隐变量为离散型随机变量，观测变量可以是离散型的也可以是连续型的

---

程序如下：

import numpy as np
from hmmlearn import hmm

# 隐状态集合Q
states = ['box1', 'box2', 'box3']
# 观测集合V
observations = ['black', 'white']
# 初始概率pi
start_probability = np.array([0.3, 0.5, 0.2])
# 状态转移矩阵A
transition_probability = np.array([
  [0.4, 0.4, 0.2],
  [0.3, 0.2, 0.5],
  [0.2, 0.6, 0.2]
])

# 观测概率矩阵B
emission_probability = np.array([
  [0.2, 0.8],
  [0.6, 0.4],
  [0.4, 0.6]
])

# 选用MultinomialHMM对离散观测状态建模
model = hmm.MultinomialHMM(n_components=len(states))
model.startprob_ = start_probability
model.transmat_ = transition_probability
model.emissionprob_ = emission_probability

# 观测序列
obervation_list = np.array([0, 1, 0])
# 调用维特比算法对观测序列进行隐含状态解码
logprob, box_list = model.decode(obervation_list.reshape(-1, 1), algorithm='viterbi')

# 对数最大路径概率
print(logprob)

# 输出解码的隐含状态序列
print(box_list)
for i in range(len(obervation_list)):
    print(states[box_list[i]])

输出为：

-3.429596856183853
[1 2 1]
box2
box3
box2

可见，程序计算的结果与手工计算结果一致；

---

对于HMM的学习问题，HMM属于无监督学习，我们需要先设置HMM的隐状态数，再根据观测序列，由EM算法迭代更新模型参数$\lambda =(\pi ,A,B)$

---