本文已参与「新人创作礼」活动，一起开启掘金创作之路。

0. 背景介绍

判别式模型,对条件分布进行建模；生成式模型，对联合分布进行建模。概率图模型是一类利用图来表达变量间相关关系的概率模型，既包含判别式模型，也包含生成式模型。

1. 隐马尔科夫模型（Hidden Markov Model, HMM）

1.1. HMM的数学定义

HMM是关于时序的概率模型，描述由一个隐藏的马尔科夫链随机生成不可观测的状态随机序列，再由各个状态生成一个可观测的随机序列的过程。

它是一种生成式模型，因为它需要先生成状态变量，然后在状态的条件下生成观测变量，根据条件概率的公式，等价于生成状态变量和观测变量的联合概率。所以是生成式模型。

HMM模型常用于序列标注问题建模，假设时间长度为n，对于 $i=1,2,\cdots,n$ 时刻，HMM中的变量可以分为两组，用 $X=\{x_1,x_2,\cdots,x_n\},x_i\in O=\{o_1,o_2,\cdots,o_M\}$ 表示观测变量序列（每个观测变量有 $M$ 个可能的取值）； $Y=\{y_1,y_2,\cdots,y_n\},y_i\in S=\{s_1,s_2,\cdots,s_N\}$ （每个状态变量有 $N$ 个可能的取值）表示状态变量序列，也叫隐变量。并且有如下假设：在任一时刻，观测变量的取值仅仅依赖于该时刻的状态变量，即 $x_t$ 仅由 $y_t$ 决定；

P(x_t|y_T,x_T,y_{T-1},x_{T-1},\cdots,y_1,x_1)=P(x_t|y_t))\\ \tag{1-1}

同时， $y_t$ 仅依赖于 $y_{t-1}$ ，即状态变量序列是一个马尔科夫链。由此可以得到所有变量的联合概率分布为（注意根据假设，条件概率的条件省去了无关变量）：

P(x_1,y_1,\cdots,x_n,y_n)=P(y_1)P(x_1|y_1)\prod_{i=2}^{n}P(y_i|y_{i-1})P(x_i|y_i)\\ \tag{1-2}

除了前面提到的状态空间 $Y$ 和观测空间 $X$ ，要想确定一个HMM模型还需要如下三组参数：

状态转移概率矩阵（隐藏状态到隐藏状态之间的概率转换）： $A=[a_{ij}]_{N\times N}$ ，其中 $a_{ij}=P(y_{t+1}=s_j|y_t=s_i),1\leq i,j\leq N$ .为了便于后文公式推导中的符号说明，我们引入 $a_{s_is_j}=P(y_{t+1}=s_j|y_t=s_i),1\leq i,j\leq N$ 。有时候也用省略的写法 $a_{ij}$
观测概率矩阵（也叫发射概率，隐藏状态到观测状态之间的概率转换）： $B=[b_{ij}]_{N\times M}$ ，其中 $b_{ij}$ 表示 $P(x_t=o_j|y_t=s_i),1\leq i\leq N,1\leq j\leq M$ .为了便于后文公式推导中的符号说明，我们引入 $b_{s_io_j}$ 表示 $P(x_t=o_j|y_t=s_i),1\leq i\leq N,1\leq j\leq M$ 。有时候也用省略的写法 $b_{ij}$
初始状态概率向量： $\pi=(\pi_1,\pi_2,\cdots,\pi_N)$ ，其中 $\pi_i=P(y_1=s_i)$

记一个HMM为 $\lambda=[A,B,\pi]$ ，则可按照如下过程产生观测序列：

设置t=1，根据初始状态概率 $\pi$ 选择初始状态 $y_1$ （比如根据概率最大的原则选取）;
根据状态 $y_t$ 和输出观测概率矩阵 $B$ 确定观测变量 $x_t$ ；
根据状态 $y_t$ 和状态转移矩阵 $A$ 确定 $y_{t+1}$ ；
若 $t< n$ ，设置 $t=t+1$ ，并跳转到2；否则停止。

1.2. HMM的实际应用问题描述

在实际应用中，我们常常关注HMM的下面三个问题：

观测序列概率计算问题：给定模型 $\lambda=[A,B,\pi]$ 和观测序列 $X=\{x_1,x_2,\cdots,x_n\}$ ，如何有效计算出观测序列出现的概率 $P(X|\lambda)$ 。换言之，如何评估模型与观测序列之间的匹配程度（因为计算出观测序列的时候，可能会有多个选择要选择最匹配的那个结构）。例如，根据时间序列的历史值预测未来值（已知n时刻及之前的x，对n+1时刻可能的x都拿来计算观测序列的概率，挑使得观测序列概率最大的状态作为n+1时刻的x）；这个过程称为前向算法
状态序列预测问题算法：给定模型 $\lambda=[A,B,\pi]$ 和观测序列 $X=\{x_1,x_2,\cdots,x_n\}$ ，如何找到与此序观测序列最匹配的状态序列 $Y=\{y_1,y_2,\cdots,y_n\},y_i\in\{s_1,s_2,\cdots,s_N\}$ 。换言之，如何根据观测序列推断出隐藏状态。例如，在语音识别ASR中，如何根据语音信号（观测值）确定文字序列（隐藏状态）；主流的算法是维特比算法
参数学习问题：给定观测序列 $X=\{x_1,x_2,\cdots,x_n\}$ ，如何调整模型参数 $\lambda=[A,B,\pi]$ 使得该序列出现的概率 $P(X|\lambda)$ 最大。换言之，如何训练模型使其能最好滴描述观测数据。例如，根据梯度下降的方式在序列样本上学习得到最优的模型参数；EM算法前后向算法

1.2.1. 观测序列概率计算问题

给定模型 $\lambda=[A,B,\pi]$ 和观测序列 $X=\{x_1,x_2,\cdots,x_n\}$ ，如何有效计算出观测序列出现的概率 $P(X|\lambda)$ 。

1.2.1.1. 直接暴力计算法

此方法是最简单粗暴的算法，但是复杂度过高，实际中一般不采用。其他算法一般都是在此算法上优化得到的，因此有必要了解一下直接计算法的原理。直接计算法的思路是，通过列举所有可能的长度为T的状态序列 $Y=\{y_1,y_2,\cdots,y_T\}$ ，求各个状态序列Y与观测序列 $X=\{x_1,x_2,\cdots,x_T\}$ 的联合概率 $P(X,Y|\lambda)$ ，然后利用全概率对所有可能的状态序列求和得到 $P(X|\lambda)$ 。对于某个具体的状态序列 $Y=\{y_1,y_2,\cdots,y_T\}$ ，其概率为：

P(Y|\lambda)=\pi_{y_1}a_{y_1y_2}a_{y_2y_3}\cdots a_{y_{T-1}y_T}\tag{1-3}

给定上述状态序列，观测序列 $X=\{x_1,x_2,\cdots,x_T\}$ 的概率为：

P(X|Y,\lambda)=b_{y_1,x_1}b_{y_2,x_2}\cdots b_{y_Tx_T}\tag{1-4}

则状态序列和观测序列的联合概率可以表示为：

\begin{aligned} P(X,Y|\lambda)&=P(X|Y,\lambda)P(Y|\lambda)\\ &=\pi_{y_1}b_{y_1x_1}a_{y_1y_2}b_{y_2x_2}a_{y_2y_3}\cdots b_{y_{T-1}x_{T-1}}a_{y_{T-1}y_T}b_{y_Tx_T}\tag{1-5} \end{aligned}

则观测序列的概率可以有全概率公式求得:

\begin{aligned} P(X|\lambda)&=\sum_YP(X|Y,\lambda)P(Y|\lambda)\\ &=\sum_{y_1\in S,y_2\in S,\cdots,y_T\in S}\pi_{y_1}b_{y_1x_1}a_{y_1y_2}b_{y_2x_2}a_{y_2y_3}\cdots b_{y_{T-1}x_{T-1}}a_{y_{T-1}y_T}b_{y_Tx_T}\tag{1-6} \end{aligned}

上式的计算复杂度高达 $O(TN^T)$

1.2.1.2. 基于动态规划的前向算法

直接计算法复杂度过高的原因在于存在大量冗余计算。例如，当t时刻之前的状态序列固定时，穷举t时刻之后的状态序列的各种取值情况时，都会重复计算t时刻之前的各种取值情况。因此可以采用动态规划的思想，求得计算过程中的递推式，减少重复计算。定义到t时刻时，观测序列 $\{x_1,x_2,\cdots,x_t\}$ 且状态为 $s_i$ 的概率为前向概率，记作：

\alpha_t(s_i)=P(x_1,x_2,\cdots,x_t,y_t=s_i|\lambda)\tag{1-7}

进行如下递推计算过程：初值t=1：

\alpha_1(s_i)=\pi_{s_i}b_{s_ix_i},i=1,2,\cdots,N\tag{1-8}

对于 $t=1,2,\cdots,T-1$ :

\alpha_{t+1}(s_i)=\left[\sum_{j=1}^N\alpha_t(s_j)a_{s_js_i}\right]b_{s_ix_{t+1}},i=1,2,\cdots,N\tag{1-9}

终止值：

P(O|\lambda)=\sum_{i=1}^N\alpha_T(s_i)\tag{1-10}

1.2.1.3. 《统计学习方法》例10.2python代码

项目地址：github.com/lankuohsing…

# 状态 1 2 3之间的转移概率矩阵
A = [[0.5,0.2,0.3],
     [0.3,0.5,0.2],
     [0.2,0.3,0.5]]
# 初始状态概率
pi = [0.2,0.4,0.4]
# 每个袋子里，红白求的概率
# red white
B = [[0.5,0.5],
     [0.4,0.6],
     [0.7,0.3]]

# In[]
#前向算法
def hmm_forward(A,B,pi,O):
    T = len(O)# 观测序列长度
    N = len(A[0])# 状态个数
    #step1 初始化
    alpha = [[0]*T for _ in range(N)]# 每行代表不同的状态，每列代表不同的观测时刻
    for i in range(N):
        alpha[i][0] = pi[i]*B[i][O[0]]

    #step2 计算alpha(t)
    for t in range(1,T):
        for i in range(N):
            temp = 0
            for j in range(N):
                temp += alpha[j][t-1]*A[j][i]
            alpha[i][t] = temp*B[i][O[t]]

    #step3
    proba = 0
    for i in range(N):
        proba += alpha[i][-1]
    return proba,alpha

A = [[0.5,0.2,0.3],
     [0.3,0.5,0.2],
     [0.2,0.3,0.5]]
B = [[0.5,0.5],
     [0.4,0.6],
     [0.7,0.3]]
pi = [0.2,0.4,0.4]
O = [0,1,0]
proba,alpha=hmm_forward(A,B,pi,O)  #结果为 0.130218
print(proba)

1.2.2. 状态序列预测问题算法

给定模型 $\lambda=[A,B,\pi]$ 和观测序列 $X=\{x_1,x_2,\cdots,x_n\}$ ，如何找到与此序观测序列最匹配的状态序列 $Y=\{y_1,y_2,\cdots,y_n\},y_i\in\{s_1,s_2,\cdots,s_N\}$ 。

1.2.2.1. 直接暴力计算法

穷举所有的状态序列的可能情况，分别计算它们的概率并求最大。此方法时间复杂度过高，基本没有实用价值

1.2.2.2. Viterbi算法

（注：下文的起点，都是默认指0时刻从起点出发的意思）该算法是基于动态规划的最短路径算法，在其他最短路径问题中也具有广泛应用。它的思想是：记末尾时刻为 $T$ ，如果最优路径在时刻t通过结点 $i_t^*$ ，那么这条路径从结点 $i_t^*$ 到终点 $i_T^*$ 的部分路径，对于从 $i_t^*$ 到 $i_T^*$ 的所有子路径来说必须是最优的，同理，这条路径中从起点到结点 $i_t^*$ 的部分路径也是从起点到结点 $i_t^*$ 的所有子路径中最优的。根据这一原理，我们只需从时刻 $t=1$ 开始，递推地计算从起点到在时刻t状态为i的结点的各条子路径中的最优路径，并记录下每个结点的上一时刻的最优结点，当计算到终点时就完成了整个最优路径的计算。

定义从起点到时刻t状态为i的所有子路径 $(i_1,i_2,\cdots,i_t)$ 中概率最大值为：

\delta_t(i)=\mathop {\max }\limits_{i_1,i_2,\cdots,i_{t-1}}P(i_t=i,i_{t-1},\cdots,i_1,x_t,\cdots,x_1|\lambda),i=1,2,\cdots,N\tag{1-11}

那么可以得到 $\delta_{t+1}(i)$ 的递推式：

\begin{aligned} \delta_{t+1}(i)&=\mathop {\max }\limits_{i_1,i_2,\cdots,i_{t}}P(i_t=i,i_{t},\cdots,i_1,x_{t+1},\cdots,x_1|\lambda)\\ &=\mathop {\max }\limits_{1\leq j\leq M}[\delta_t(j)a_{ji}]b_{ix_{t+1}},i=1,2,\cdots,N;t=1,2,\cdots,T-1\\ \tag{1-12} \end{aligned}

定义从起点到在时刻t状态为i的结点所有单个路径 $\{(i_1,i_2,\cdots,i_{t-1},t)\}$ 中概率最大的路径的第 $t-1$ 时刻的结点为：

\psi_t(i)=arg \mathop{\max }\limits_{1\leq j\leq N}[\delta_{t-1}(j)a_{ji}]\tag{1-13}

那么viterbi算法的流程就呼之欲出了： (1) 赋予初值：

\begin{aligned} \delta_1(i)&=\pi_ib_{ix_1},i=1,2,\cdots,N\\ \psi_i(i)&=0,i=1,2,\cdots,N\\ \tag{1-14} \end{aligned}

(2) 递推式, $t=2,3,\cdots,T$ ：

\begin{aligned} \delta_t(i)&=\mathop {\max }\limits_{1\leq j\leq N}[\delta_{t-1}(j)a_{ji}]b_{ix_{t+1}},i=1,2,\cdots,N\\ \psi_t(i)&=arg \mathop{\max }\limits_{1\leq j\leq N}[\delta_{t-1}(j)a_{ji}],i=1,2,\cdots,N\\ \tag{1-15} \end{aligned}

(3) 终止，计算完整的最优路径和概率：

\begin{aligned} P^*&=\mathop {\max }\limits_{1\leq i\leq N}\delta_T(i)\\ i_T^*&=arg \mathop{\max }\limits_{1\leq j\leq N}[\delta_{t}(i)]\\ i_t^*&=\psi_{t+1}(i_{t+1}^*),t=T-1,T-2,\cdots,1\\ \tag{1-16} \end{aligned}

1.3. HMM实际例子

在这里插入图片描述

1.4. 马尔科夫随机场（Markov Random Field, MRF）

MRF是典型的马尔科夫网，是一种著名的无向图模型，图中每个结点表示一个或一组变量，结点之间的边表示两个变量之间的依赖关系。定义势函数（亦称“因子”）为变量子集上的非负实函数，主要用于定义概率分布函数。对于结点的一个子集，若其中任意两点结点之间都有边连接，则称该节点子集为一个“团”。若在一个团中加入任何一个节点都不再形成团，则称该团为“极大团”。对于n各变量 $X=\{x_1,x_2,\cdots,x_n\}$ ，记其所有团构成的集合为 $C$ ，与团 $Q\in C$ 对应的变量集合记为 $X_Q$ ，则联合概率 $P(X)$ 定义为

P(X)=\frac{1}{Z}\prod_{Q\in C}\psi(X_Q)\tag{1-2}

其中 $\psi(X_Q)$ 为与团 $Q$ 对应的是函数， $Z=\sum_X\prod_{Q\in C}\psi(X_Q)$ 为规范化因子，用以确保上式的值域为0~1. 为了减少需要计算的团的数目，有如下推论：若 $Q$ 不是极大团，则它必被一个极大团 $Q^*$ ，于是 $P(X)$ 可以基于极大团来定义。假设所有极大团构成的集合为 $C^*$ ，则有

P(X)=\frac{1}{Z^*}\prod_{Q\in C^*}\psi(X_Q)\tag{1-3}

假设结点集合A中的结点到结点集B中的点都必须经过结点集C中的结点，则称结点集A和B被结点集C分离，C称为“分离集”。对于马尔科夫随机场，有“全局马尔科夫性”：给定两个变量子集的分离集，则这两个变量子集条件独立，也即有

P(x_A,x_B|x_C)=P(x_A|x_C)P(x_B|x_C)\tag{1-4}

具体证明见《机器学习》（周志华）324. 由全局马尔科夫性可以得到两个推论：

局部马尔科夫性：给定某变量的邻接变量，则该变量条件独立与其他变量
成对马尔科夫性：给定所有其他变量，两个非邻接变量条件独立。

2. 条件随机场（Conditional Random Field， CRF）

2.1. CRF的数学定义

令 $X=\{x_1,x_2,\cdots,x_n\},x_i\in \{o_1,o_2,\cdots,o_M\}$ 为观测序列， $Y=\{y_1,y_2,\cdots,y_n\},y_i\in\{s_1,s_2,\cdots,s_N\}$ 为状态序列（标记或者叫标签），CRF的目的是构建条件概率模型 $P(Y|X)$ 。注意， $Y$ 的各分量之间可以是相关的，例如nlp的词性标注任务中，观测数据为句子的单词序列，标记数据为相应的词性序列。令 $G=\langle V,E\rangle$ 表示结点与标记变量 $Y$ 中元素一一对应的无向图， $y_v$ 表示与结点 $v$ 对应的标记变量， $n(v)$ 表示结点 $v$ 的邻接节点，若图 $G$ 中每个变量 $y_v$ 都满足马尔科夫性，即

P\left(y_v|X,Y_{V\backslash\{v\}}\right)=P(y_v|X,Y_n(v))\tag{2-1}

则 $(Y,X)$ 构成一个条件随机场。

以链式条件随机场为例，CRF的条件概率被定义为

P(Y|X)=\frac{1}{Z}exp\left(\sum_j\sum_{i=1}^{n-1}\lambda_j t_j(y_{i+1},y_i,X,i)+\sum_k\sum_{i=1}^{n}\mu_k s_k(y_i,X,i)\right)\tag{2-2}

其中 $t_j(y_{i+1},y_i,X,i)$ 是定义在观测序列的两个相邻标记位置上的转移特征函数，用于刻画向量标记变量之间的相关关系以及观测序列对它们的影响； $s_k(y_i,X,i)$ 是定义在观测序列的标记位置 $i$ 上的状态特征函数，用于刻画观测序列对标记变量的影响； $\lambda_j$ 和 $\mu_k$ 为参数， $Z$ 为规范化因子。

2.2. 特征函数的例子

以词性标注为例， $t_j(y_{i+1},y_i,X,i)=1,if\ y_{i+1}=[P],y_i=[V]\ and\ x_i="knock";0\ otherwise$

$s_k(y_i,X,i)=1,if\ y_i=[V]\ and\ x_i="knock";0\ otherwise$

参考资料：

《机器学习》——周志华
《统计学习方法》——李航
blog.csdn.net/continueoo/…
blog.csdn.net/continueoo/…

机器学习之概率图模型