隐马尔科夫模型HMM之前后向算法Python代码实现，包括2个优化版本

☕️ 本文系列文章汇总：

（1）HMM开篇：基本概念和几个要素

（2）HMM计算问题：前后向算法

（3）HMM学习问题：Baum-Welch算法

（4） HMM预测问题：维特比算法

本篇算法原理分析及公式推导请参考： HMM计算问题：前后向算法

之前发布的四篇隐马尔科夫模型系列学习博文，将隐马尔科夫模型包含的各个算法的原理及公式都详细介绍了，所谓无代码实现无以成体系，那么接下来我会带大家把算法用代码实现一遍，语言主要为python。本篇先来对前后向算法进行coding~

1. 参数初始化：

        """初始化模型参数:param pi: 初始状态概率向量:param A: 已学习得到的状态转移概率矩阵，这里直接引用课本中的例子10.2:param B: 已学习得到的概率矩阵，这里直接引用课本中的例子10.2:param V: 已知的观测集合，同样使用例子中的值"""self.pi = piself.A = Aself.B = Bself.V = V

2. 根据原理公式，定义前向算法计算过程：

    def forward(self, O):"""前向算法:param O: 已知的观测序列:return: P(O|λ)"""row, col = len(O), self.A.shape[0]alpha_t_plus_1 = np.zeros((row, col), dtype=float)for t, o in enumerate(O):if t == 0:# 初值α 公式10.15for i, p in enumerate(self.pi):obj_index = self.V.index(o)alpha_t_plus_1[t][i] = p * self.B[i][obj_index]else:# 递推 公式10.16for i in range(self.A.shape[0]):alpha_ji = 0.# 公式10.16里中括号的内容for j, a in enumerate(alpha_t_plus_1[t-1]):alpha_ji += (a * self.A[j][i])obj_index = self.V.index(o)# 公式10.16alpha_t_plus_1[t][i] = alpha_ji * self.B[i][obj_index]return alpha_t_plus_1

该方法主要输出的是前向算法中alpha的中间结果，即原理中提到的。

3. 定义后向算法计算过程：

    def backward(self, O):"""后向算法:param O: 已知的观测序列:return: P(O|λ)"""row, col = len(O), self.A.shape[0]betaT = np.zeros((row+1, col), dtype=float)for t, o in enumerate(O[::-1]):if t == 0:# 初值β 公式10.19betaT[t][:] = [1] * self.A.shape[0]continueelse:# 反向递推 公式10.20for i in range(self.A.shape[0]):beta_t = 0.obj_index = self.V.index(O[t - 1])for j, b in enumerate(betaT[t-1]):beta_t += (self.A[i][j] * self.B[j][obj_index] * b)betaT[t][i] = beta_tbetaT[-1][:] = [self.pi[i] * self.B[i][self.V.index(O[0])] * betaT[-2][i] for i in range(self.A.shape[0])]return betaT

该方法主要输出的是前向算法中beta的中间结果，即原理中提到的。

4. 定义概率计算函数：

    def cal_prob(self, O, opt):if opt == 'f':metrix = self.forward(O)# 计算前向算法P(O|λ) 公式10.17return sum(metrix[-1])elif opt == 'b':# 计算后向算法P(O|λ) 公式10.21metrix = self.backward(O)return sum(metrix[-1])

将上述方法封装成一个类：

import numpy as npclass FB:def __init__(self, pi, A, B, V):"""初始化模型参数:param pi: 初始状态概率向量:param A: 已学习得到的状态转移概率矩阵，这里直接引用课本中的例子10.2:param B: 已学习得到的概率矩阵，这里直接引用课本中的例子10.2:param V: 已知的观测集合，同样使用例子中的值"""self.pi = piself.A = Aself.B = Bself.V = Vdef cal_prob(self, O, opt):if opt == 'f':metrix = self.forward(O)# 计算P(O|λ) 公式10.17return sum(metrix[-1])elif opt == 'b':# 计算P(O|λ) 公式10.21metrix = self.backward(O)return sum(metrix[-1])def forward(self, O):"""前向算法:param O: 已知的观测序列:return: P(O|λ)"""row, col = len(O), self.A.shape[0]alpha_t_plus_1 = np.zeros((row, col), dtype=float)for t, o in enumerate(O):if t == 0:# 初值α 公式10.15for i, p in enumerate(self.pi):obj_index = self.V.index(o)alpha_t_plus_1[t][i] = p * self.B[i][obj_index]else:# 递推 公式10.16for i in range(self.A.shape[0]):alpha_ji = 0.# 公式10.16里中括号的内容for j, a in enumerate(alpha_t_plus_1[t-1]):alpha_ji += (a * self.A[j][i])obj_index = self.V.index(o)# 公式10.16alpha_t_plus_1[t][i] = alpha_ji * self.B[i][obj_index]return alpha_t_plus_1def backward(self, O):"""后向算法:param O: 已知的观测序列:return: P(O|λ)"""row, col = len(O), self.A.shape[0]betaT = np.zeros((row+1, col), dtype=float)for t, o in enumerate(O[::-1]):if t == 0:# 初值β 公式10.19betaT[t][:] = [1] * self.A.shape[0]continueelse:# 反向递推 公式10.20for i in range(self.A.shape[0]):beta_t = 0.obj_index = self.V.index(O[t - 1])for j, b in enumerate(betaT[t-1]):beta_t += (self.A[i][j] * self.B[j][obj_index] * b)betaT[t][i] = beta_tbetaT[-1][:] = [self.pi[i] * self.B[i][self.V.index(O[0])] * betaT[-2][i] for i in range(self.A.shape[0])]return betaT

用书上的例子运行一下：

if __name__ == '__main__':from time import time# 课本例子10.2pi = [0.2, 0.4, 0.4]a = np.array([[0.5, 0.2, 0.3],[0.3, 0.5, 0.2],[0.2, 0.3, 0.5]])b = np.array([[0.5, 0.5],[0.4, 0.6],[0.7, 0.3]])O = ['红', '白', '红']f = FB(pi=pi, A=a, B=b, V=['红', '白'])# start = time()resf = f.forward(O)resb = f.backward(O)# print(time()-start)print('α：{}\n前向算法的概率计算结果：{}'.format(resf, f.cal_prob(O, opt='f')))print('β：{}\n后向算法的概率计算结果：{}：'.format(resb, f.cal_prob(O, opt='b')))

运行结果：

另外我们注意到上面的代码中，嵌套了三层for循环，当数据量很大的时候，这样计算很耗时，那么我们来优化一下：

优化版一：

我们注意到，内层for循环主要是为了矩阵计算的，所以我们其实可以直接利用numpy来进行矩阵的运算。

import numpy as npclass FB:def __init__(self, pi, A, B, V):"""初始化模型参数:param pi: 初始状态概率向量:param A: 已学习得到的状态转移概率矩阵，这里直接引用课本中的例子10.2:param B: 已学习得到的概率矩阵，这里直接引用课本中的例子10.2:param V: 已知的观测集合，同样使用例子中的值"""self.pi = np.array(pi)self.A = np.array(A)self.B = np.array(B)self.V = Vdef cal_prob(self, O, opt):if opt == 'f':metrix = self.forward(O)# 计算P(O|λ) 公式10.17return sum(metrix[-1])elif opt == 'b':# 计算P(O|λ) 公式10.21metrix = self.backward(O)return sum(metrix[-1])def forward(self, O):"""前向算法:param O: 已知的观测序列:return: P(O|λ)"""row, col = len(O), self.A.shape[0]alpha_t_plus_1 = np.zeros((row, col), dtype=float)obj_index = self.V.index(O[0])# 初值α 公式10.15alpha_t_plus_1[0][:] = self.pi * self.B[:].T[obj_index]for t, o in enumerate(O[1:]):t += 1# 递推 公式10.16obj_index = self.V.index(o)for i in range(self.A.shape[0]):# 公式10.16alpha_ji = alpha_t_plus_1[t-1][:] @ self.A[:].T[i]alpha_t_plus_1[t][i] = alpha_ji * self.B[i][obj_index]return alpha_t_plus_1def backward(self, O):"""后向算法:param O: 已知的观测序列:return: P(O|λ)"""row, col = len(O), self.A.shape[0]betaT = np.zeros((row+1, col), dtype=float)# 初值β 公式10.19betaT[0][:] = [1] * self.A.shape[0]for t, o in enumerate(O[::-1][1:]):t += 1# 反向递推 公式10.20obj_index = self.V.index(O[t-1])for i in range(self.A.shape[0]):beta_t = self.A[i][:] * self.B[:].T[obj_index] @ betaT[t-1][:].TbetaT[t][i] = beta_tbetaT[-1][:] = [self.pi[i] * self.B[i][self.V.index(O[0])] * betaT[-2][i] for i in range(self.A.shape[0])]return betaT

优化版二：

优化版一中，只省略了最内层循环，即每个状态的每一列。其实近一步观察发现，对于每一行的for循环也可以直接用矩阵计算的方式省略

import numpy as npclass FB:def __init__(self, pi, A, B, V):"""初始化模型参数:param pi: 初始状态概率向量:param A: 已学习得到的状态转移概率矩阵，这里直接引用课本中的例子10.2:param B: 已学习得到的概率矩阵，这里直接引用课本中的例子10.2:param V: 已知的观测集合，同样使用例子中的值"""self.pi = np.array(pi)self.A = np.array(A)self.B = np.array(B)self.V = Vdef cal_prob(self, O, opt):if opt == 'f':metrix = self.forward(O)# 计算P(O|λ) 公式10.17return sum(metrix[-1])elif opt == 'b':# 计算P(O|λ) 公式10.21metrix = self.backward(O)return sum(metrix[-1])def forward(self, O):"""前向算法:param O: 已知的观测序列:return: P(O|λ)"""row, col = len(O), self.A.shape[0]alpha_t_plus_1 = np.zeros((row, col), dtype=float)obj_index = self.V.index(O[0])# 初值α 公式10.15alpha_t_plus_1[0][:] = self.pi * self.B[:].T[obj_index]for t, o in enumerate(O[1:]):t += 1# 递推 公式10.16obj_index = self.V.index(o)alpha_ji = alpha_t_plus_1[t - 1][:].T @ self.Aalpha_t_plus_1[t][:] = alpha_ji * self.B[:].T[obj_index]return alpha_t_plus_1def backward(self, O):"""后向算法:param O: 已知的观测序列:return: P(O|λ)"""row, col = len(O), self.A.shape[0]betaT = np.zeros((row+1, col), dtype=float)# 初值β 公式10.19betaT[0][:] = [1] * self.A.shape[0]for t, o in enumerate(O[::-1][1:]):t += 1# 反向递推 公式10.20obj_index = self.V.index(O[t-1])beta_t = self.A * self.B[:].T[obj_index] @ betaT[t-1][:].TbetaT[t][:] = beta_tbetaT[-1][:] = [self.pi[i] * self.B[i][self.V.index(O[0])] * betaT[-2][i] for i in range(self.A.shape[0])]return betaT

好啦，就是这么一步步递进优化，方法论就是，先用基本的代码将原理实现，帮助理解，然后发现竟然有3个for循环，复杂度且不说，代码首先就不美观，看看有没有可以优化的地方，发现for循环的服务对象是矩阵计算，所以，很自然的想到直接利用矩阵计算的方法一步到位。代码已经放到GitHub上了，我将会持续更新其它算法的实现。

隐马尔科夫模型HMM之前后向算法Python代码实现，包括2个优化版本相关推荐

隐马尔科夫模型(HMM)的无监督学习算法java实现（baum-welch迭代求解），包括串行以及并行实现
HMM的原理就不说了,这里主要说算法的实现. 实际实现起来并不是很困难,前提是你仔细看过hmm的原理,然后很多实现就照着公式写出对应的代码,比如前向算法,后向算法,参数更新都是有明确的公式的,只需要对 ...
隐马尔科夫模型HMM自学（3）
Viterbi Algorithm 本来想明天再把后面的部分写好,可是睡觉今天是节日呢?一时情不自禁就有打开电脑.......... 找到可能性最大的隐含状态序列崔晓源翻译多数情况下,我们都希望 ...
隐马尔科夫模型HMM自学（2）
HMM 定义崔晓源翻译 HMM是一个三元组 (,A,B). the vector of the initial state probabilities; the state transitio ...
隐马尔科夫模型HMM自学（1）
介绍崔晓源翻译我们通常都习惯寻找一个事物在一段时间里的变化规律.在很多领域我们都希望找到这个规律,比如计算机中的指令顺序,句子中的词顺序和语音中的词顺序等等.一个最适用的例子就是天气的预测. 首 ...
【NLP】用于语音识别、分词的隐马尔科夫模型HMM
大家好,今天介绍自然语言处理中经典的隐马尔科夫模型(HMM).HMM早期在语音识别.分词等序列标注问题中有着广泛的应用. 了解HMM的基础原理以及应用,对于了解NLP处理问题的基本思想和技术发展脉络有 ...
第九章隐马尔科夫模型HMM
文章目录 1 隐马尔科夫模型定义 2 概率计算算法 2.1 前向概率 2.2 概率计算 3 学习算法 3.1 EM算法 3.2EM在HMM 4 预测算法 1 隐马尔科夫模型定义隐马尔科夫模型是一个s ...
python地图匹配_基于隐马尔科夫模型(HMM)的地图匹配(Map-Matching)算法
1. 摘要本篇博客简单介绍下用隐马尔科夫模型(Hidden Markov Model, HMM)来解决地图匹配(Map-Matching)问题.转载请注明网址. 2. Map-Matching(MM ...
隐马尔科夫模型 HMM 与语音识别 speech recognition (1)：名词解释
0.引言想在 CSDN 上看一下隐马尔科夫模型,简称HMM(Hidden Markov Model)的例子,找了几篇博文,却发现大部分都是转载的,转载的还没有出处,文中的表述与逻辑也看的人晕头转向, ...
隐马尔科夫模型HMM（三）鲍姆-韦尔奇算法求解HMM参数
在本篇我们会讨论HMM模型参数求解的问题,这个问题在HMM三个问题里算是最复杂的.在研究这个问题之前,建议先阅读这个系列的前两篇以熟悉HMM模型和HMM的前向后向算法,以及EM算法原理总结,这些在本篇 ...

隐马尔科夫模型HMM之前后向算法Python代码实现，包括2个优化版本

隐马尔科夫模型HMM之前后向算法Python代码实现，包括2个优化版本相关推荐

最新文章

热门文章