这篇文章主要讲解了“python模拟隐马尔可夫模型的方法是什么”,文中的讲解内容简单清晰,易于学习与理解,下面请大家跟着小编的思路慢慢深入,一起来研究和学习“python模拟隐马尔可夫模型的方法是什么”吧!
import numpy as np class HiddenMarkov: def __init__(self): self.alphas = None self.forward_P = None self.betas = None self.backward_P = None # 前向算法 # Q 是状态集合,里面包含了所有可能的状态 # V 是我们的观测的集合,里面包含了所有可能的观测结果 # A 状态转移概率分布 # B 观测概率分布 # O 观测序列,依次为观测值 # PI 初始概率分布。根据这个先生成初始状态。 def forward(self, Q, V, A, B, O, PI): # 状态序列的大小 N = len(Q) # 观测序列的大小 M = len(O) # 初始化前向概率alpha值 alphas = np.zeros((N, M)) # 时刻数=观测序列数 T = M # 遍历每一个时刻,计算前向概率alpha值 for t in range(T): # 得到序列对应的索引 indexOfO = V.index(O[t]) # 遍历状态序列 for i in range(N): # 初始化alpha初值 if t == 0: # P176 公式(10.15) alphas[i][t] = PI[t][i] * B[i][indexOfO] print('alpha1(%d) = p%db%db(o1) = %f' % (i + 1, i, i, alphas[i][t])) else: # P176 公式(10.16) alphas[i][t] = np.dot([alpha[t - 1] for alpha in alphas], [a[i] for a in A]) * B[i][indexOfO] print('alpha%d(%d) = [sigma alpha%d(i)ai%d]b%d(o%d) = %f' % (t + 1, i + 1, t - 1, i, i, t, alphas[i][t])) # P176 公式(10.17) self.forward_P = np.sum([alpha[M - 1] for alpha in alphas]) self.alphas = alphas # 后向算法 # Q 是状态集合,里面包含了所有可能的状态 # V 是我们的观测的集合,里面包含了所有可能的观测结果 # A 状态转移概率分布 # B 观测概率分布 # O 观测序列,依次为观测值 # PI 初始概率分布。根据这个先生成初始状态。 def backward(self, Q, V, A, B, O, PI): # 状态序列的大小 N = len(Q) # 观测序列的大小 M = len(O) # 初始化后向概率beta值,P178 公式(10.19) betas = np.ones((N, M)) for i in range(N): print('beta%d(%d) = 1' % (M, i + 1)) # 对观测序列逆向遍历 for t in range(M - 2, -1, -1): # 得到序列对应的索引 indexOfO = V.index(O[t + 1]) # 遍历状态序列 for i in range(N): # P178 公式(10.20) betas[i][t] = np.dot( np.multiply(A[i], [b[indexOfO] for b in B]), [beta[t + 1] for beta in betas]) realT = t + 1 realI = i + 1 print('beta%d(%d) = sigma[a%djbj(o%d)beta%d(j)] = (' % (realT, realI, realI, realT + 1, realT + 1), end='') for j in range(N): print("%.2f * %.2f * %.2f + " % (A[i][j], B[j][indexOfO], betas[j][t + 1]), end='') print("0) = %.3f" % betas[i][t]) # 取出第一个值 indexOfO = V.index(O[0]) self.betas = betas # P178 公式(10.21) P = np.dot(np.multiply(PI, [b[indexOfO] for b in B]), [beta[0] for beta in betas]) self.backward_P = P print("P(O|lambda) = ", end="") for i in range(N): print("%.1f * %.1f * %.5f + " % (PI[0][i], B[i][indexOfO], betas[i][0]), end="") print("0 = %f" % P) # 维特比算法:动态规划解隐马尔代夫模型预测问题 # Q 是状态集合,里面包含了所有可能的状态 # V 是我们的观测的集合,里面包含了所有可能的观测结果 # A 状态转移概率分布 # B 观测概率分布 # O 观测序列,依次为观测值 # PI 初始概率分布。根据这个先生成初始状态。 def viterbi(self, Q, V, A, B, O, PI): # 状态序列的大小 N = len(Q) # 观测序列的大小 M = len(O) # 初始化daltas:存当前时刻当前状态的所有单个路径的概率最大值 deltas = np.zeros((N, M)) # 初始化psis:存当前时刻当前状态所有单个路径中概率最大路径的前一时刻结点 psis = np.zeros((N, M)) # 初始化最优路径矩阵,该矩阵维度与观测序列维度相同。这是我们最后的输出。 I = np.zeros((1, M)) # 遍历观测序列 for t in range(M): # 递推从t=2开始 realT = t + 1 # 得到序列对应的索引 indexOfO = V.index(O[t]) for i in range(N): realI = i + 1 if t == 0: # P185 算法10.5 步骤(1) deltas[i][t] = PI[0][i] * B[i][indexOfO] psis[i][t] = 0 print('delta1(%d) = pi%d * b%d(o1) = %.2f * %.2f = %.2f' % (realI, realI, realI, PI[0][i], B[i][indexOfO], deltas[i][t])) print('psis1(%d) = 0' % (realI)) else: # # P185 算法10.5 步骤(2) deltas[i][t] = np.max( np.multiply([delta[t - 1] for delta in deltas], [a[i] for a in A])) * B[i][indexOfO] print( 'delta%d(%d) = max[delta%d(j)aj%d]b%d(o%d) = %.2f * %.2f = %.5f' % (realT, realI, realT - 1, realI, realI, realT, np.max( np.multiply([delta[t - 1] for delta in deltas], [a[i] for a in A])), B[i][indexOfO], deltas[i][t])) # 对于y=f(x),argmax返回取得最大值y时的x psis[i][t] = np.argmax( np.multiply([delta[t - 1] for delta in deltas], [a[i] for a in A])) print('psis%d(%d) = argmax[delta%d(j)aj%d] = %d' % (realT, realI, realT - 1, realI, psis[i][t])) # 得到最优路径的终结点 I[0][M - 1] = np.argmax([delta[M - 1] for delta in deltas]) print('i%d = argmax[deltaT(i)] = %d' % (M, I[0][M - 1] + 1)) # 递归由后向前得到其他结点 for t in range(M - 2, -1, -1): I[0][t] = psis[int(I[0][t + 1])][t + 1] print('i%d = psis%d(i%d) = %d' % (t + 1, t + 2, t + 2, I[0][t] + 1)) # 输出最优路径 print('最优路径是:', "->".join([str(int(i + 1)) for i in I[0]])) # 习题10.1 Q = [1, 2, 3] V = ['红', '白'] 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]] # O = ['红', '白', '红', '红', '白', '红', '白', '白'] O = ['红', '白', '红', '白'] # 习题10.1的例子 PI = [[0.2, 0.4, 0.4]] HMM = HiddenMarkov() HMM.forward(Q, V, A, B, O, PI) print("P(O|λ)={}".format(HMM.forward_P)) # HMM.backward(Q, V, A, B, O, PI) # print("P(O|λ)={}".format(HMM.backward_P)) # HMM.viterbi(Q, V, A, B, O, PI)
结果
前向算法
alpha1(1) = p0b0b(o1) = 0.100000 alpha1(2) = p1b1b(o1) = 0.160000 alpha1(3) = p2b2b(o1) = 0.280000 alpha2(1) = [sigma alpha0(i)ai0]b0(o1) = 0.077000 alpha2(2) = [sigma alpha0(i)ai1]b1(o1) = 0.110400 alpha2(3) = [sigma alpha0(i)ai2]b2(o1) = 0.060600 alpha3(1) = [sigma alpha1(i)ai0]b0(o2) = 0.041870 alpha3(2) = [sigma alpha1(i)ai1]b1(o2) = 0.035512 alpha3(3) = [sigma alpha1(i)ai2]b2(o2) = 0.052836 alpha4(1) = [sigma alpha2(i)ai0]b0(o3) = 0.021078 alpha4(2) = [sigma alpha2(i)ai1]b1(o3) = 0.025188 alpha4(3) = [sigma alpha2(i)ai2]b2(o3) = 0.013824 P(O|λ)=0.06009079999999999
后向算法
beta4(1) = 1 beta4(2) = 1 beta4(3) = 1 beta3(1) = sigma[a1jbj(o4)beta4(j)] = (0.50 * 0.50 * 1.00 + 0.20 * 0.60 * 1.00 + 0.30 * 0.30 * 1.00 + 0) = 0.460 beta3(2) = sigma[a2jbj(o4)beta4(j)] = (0.30 * 0.50 * 1.00 + 0.50 * 0.60 * 1.00 + 0.20 * 0.30 * 1.00 + 0) = 0.510 beta3(3) = sigma[a3jbj(o4)beta4(j)] = (0.20 * 0.50 * 1.00 + 0.30 * 0.60 * 1.00 + 0.50 * 0.30 * 1.00 + 0) = 0.430 beta2(1) = sigma[a1jbj(o3)beta3(j)] = (0.50 * 0.50 * 0.46 + 0.20 * 0.40 * 0.51 + 0.30 * 0.70 * 0.43 + 0) = 0.246 beta2(2) = sigma[a2jbj(o3)beta3(j)] = (0.30 * 0.50 * 0.46 + 0.50 * 0.40 * 0.51 + 0.20 * 0.70 * 0.43 + 0) = 0.231 beta2(3) = sigma[a3jbj(o3)beta3(j)] = (0.20 * 0.50 * 0.46 + 0.30 * 0.40 * 0.51 + 0.50 * 0.70 * 0.43 + 0) = 0.258 beta1(1) = sigma[a1jbj(o2)beta2(j)] = (0.50 * 0.50 * 0.25 + 0.20 * 0.60 * 0.23 + 0.30 * 0.30 * 0.26 + 0) = 0.112 beta1(2) = sigma[a2jbj(o2)beta2(j)] = (0.30 * 0.50 * 0.25 + 0.50 * 0.60 * 0.23 + 0.20 * 0.30 * 0.26 + 0) = 0.122 beta1(3) = sigma[a3jbj(o2)beta2(j)] = (0.20 * 0.50 * 0.25 + 0.30 * 0.60 * 0.23 + 0.50 * 0.30 * 0.26 + 0) = 0.105 P(O|lambda) = 0.2 * 0.5 * 0.11246 + 0.4 * 0.4 * 0.12174 + 0.4 * 0.7 * 0.10488 + 0 = 0.060091 P(O|λ)=[0.0600908]
维特比算法
delta1(1) = pi1 * b1(o1) = 0.20 * 0.50 = 0.10 psis1(1) = 0 delta1(2) = pi2 * b2(o1) = 0.40 * 0.40 = 0.16 psis1(2) = 0 delta1(3) = pi3 * b3(o1) = 0.40 * 0.70 = 0.28 psis1(3) = 0 delta2(1) = max[delta1(j)aj1]b1(o2) = 0.06 * 0.50 = 0.02800 psis2(1) = argmax[delta1(j)aj1] = 2 delta2(2) = max[delta1(j)aj2]b2(o2) = 0.08 * 0.60 = 0.05040 psis2(2) = argmax[delta1(j)aj2] = 2 delta2(3) = max[delta1(j)aj3]b3(o2) = 0.14 * 0.30 = 0.04200 psis2(3) = argmax[delta1(j)aj3] = 2 delta3(1) = max[delta2(j)aj1]b1(o3) = 0.02 * 0.50 = 0.00756 psis3(1) = argmax[delta2(j)aj1] = 1 delta3(2) = max[delta2(j)aj2]b2(o3) = 0.03 * 0.40 = 0.01008 psis3(2) = argmax[delta2(j)aj2] = 1 delta3(3) = max[delta2(j)aj3]b3(o3) = 0.02 * 0.70 = 0.01470 psis3(3) = argmax[delta2(j)aj3] = 2 delta4(1) = max[delta3(j)aj1]b1(o4) = 0.00 * 0.50 = 0.00189 psis4(1) = argmax[delta3(j)aj1] = 0 delta4(2) = max[delta3(j)aj2]b2(o4) = 0.01 * 0.60 = 0.00302 psis4(2) = argmax[delta3(j)aj2] = 1 delta4(3) = max[delta3(j)aj3]b3(o4) = 0.01 * 0.30 = 0.00220 psis4(3) = argmax[delta3(j)aj3] = 2 i4 = argmax[deltaT(i)] = 2 i3 = psis4(i4) = 2 i2 = psis3(i3) = 2 i1 = psis2(i2) = 3 最优路径是: 3->2->2->2
感谢各位的阅读,以上就是“python模拟隐马尔可夫模型的方法是什么”的内容了,经过本文的学习后,相信大家对python模拟隐马尔可夫模型的方法是什么这一问题有了更深刻的体会,具体使用情况还需要大家实践验证。这里是亿速云,小编将为大家推送更多相关知识点的文章,欢迎关注!
免责声明:本站发布的内容(图片、视频和文字)以原创、转载和分享为主,文章观点不代表本网站立场,如果涉及侵权请联系站长邮箱:is@yisu.com进行举报,并提供相关证据,一经查实,将立刻删除涉嫌侵权内容。