HMM模型

import numpy as npclass HMM:def __init__(self,A=None,B=None,Pi=None,O = None):if A:self.A = np.array(A)else:self.A = Noneif Pi:self.Pi = np.array(Pi)self.i_num = len(Pi)else:self.Pi = Noneself.i_num = Noneif B:self.B = np.array(B)self.o_num = np.array(B).shape[-1]else:self.B = Noneself.o_num = Noneself.O = Oif O:self.T = len(O)else:self.T = Noneself.ati = Noneself.bti = Nonedef update_ati(self):self.ati = np.zeros((self.T,self.i_num))for i in range(self.i_num):self.ati[0][i] = self.Pi[i]*self.B[i][self.O[0]]for t in range(self.T-1):for i in range(self.i_num):self.ati[t+1][i] = np.sum(self.ati[t]*self.A[:,i])*self.B[i][self.O[t+1]]def update_bti(self):self.bti = np.zeros((self.T,self.i_num))for i in range(self.i_num):self.bti[self.T-1][i] = 1for t in range(self.T-2,-1,-1):for i in range(self.i_num):self.bti[t][i] = np.sum(self.A[i,:]*self.B[:,self.O[t+1]]*self.bti[t+1])def comput_ri(self, t, i):return self.ati[t,i]*self.bti[t,i]/np.sum(self.ati[t]*self.bti[t])def comput_ei(self, t, i, j):return self.ati[t,i]*self.A[i,j]*self.B[j][self.O[t+1]]*self.bti[t+1][j] /\np.sum(self.ati[t]*self.A*self.B[:,self.O[t+1]]*self.bti[t+1])#--------------概率計算問題--------------------def predict_whole_p_o(self):self.update_ati()return np.sum(self.ati[self.T-1])def predict_q_give_i_t(self,t,i):self.update_ati()self.update_bti()return self.comput_ri(t-1,i-1)def predict_p_t_t_1_q(self,t,i,j):self.update_ati()self.update_bti()return self.comput_ei(t-1,i-1,j-1)#-------Baum-Welch學習算法--------------------def Baum_welch(self,O,i_num,max_iter=20):self.T = len(O)self.O = np.array(O)self.i_num = i_numself.A = np.zeros((i_num,i_num))+ 1/i_numself.B = np.zeros((i_num,len(set(O)))) + 1/len(set(O))self.Pi = np.zeros(i_num) + 1/i_numfor n in range(max_iter):self.update_ati()self.update_bti()for i in range(i_num):self.Pi[i] = self.comput_ri(0,i)for j in range(i_num):sum_et = 0sum_rt = 0for t in range(self.T-1):sum_et += self.comput_ei(t,i,j)sum_rt += self.comput_ri(t,i)self.A[i][j] = sum_et/sum_rtfor j in range(i_num):for k in range(len(set(O))):sum_rt1 = 0sum_rt2 = 0for t in range(self.T):rtj = self.comput_ri(t,j)if self.O[t] == k:sum_rt1 += rtjsum_rt2 += rtjself.B[j][k] = sum_rt1/sum_rt2#-------維特比預測算法--------------------def Viterbi(self):p_marix = np.zeros((self.T,self.i_num))rout = [[0] for i in range(self.i_num)]for i in range(self.i_num):p_marix[0][i] = self.Pi[i]*self.B[i][self.O[0]]for t in range(1,self.T):for i in range(self.i_num):max_p = 0max_rout = 0index_rout = 0for j in range(self.i_num):max_p = max(max_p,p_marix[t-1][j]*self.A[j][i]*self.B[i][self.O[t]])if p_marix[t-1][j]*self.A[j][i] > max_rout:max_rout = p_marix[t-1][j]*self.A[j][i]index_rout = jp_marix[t][i] = max_prout[i].append(index_rout)P_max = max(p_marix[self.T-1])max_index = np.argmax(p_marix[self.T-1])rout[max_index].append(max_index)return P_max , [x+1 for x in rout[max_index][1:]]def main():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,1,1,0,1,0,0]HMM_test = HMM(A=A,B=B,Pi=Pi,O=O)print(HMM_test.predict_whole_p_o())print(HMM_test.predict_q_give_i_t(4,3))print(HMM_test.predict_p_t_t_1_q(3, 1, 2))print(HMM_test.ati)print(HMM_test.bti)P_max, rout = HMM_test.Viterbi()print(P_max, rout)HMM_test_predict = HMM()HMM_test_predict.Baum_welch(O=O,i_num=3,max_iter=50)P_max, rout = HMM_test_predict.Viterbi()print(P_max, rout)if __name__ == '__main__':main()###########result################### /usr/bin/python3 /Users/zhengyanzhao/PycharmProjects/tongjixuexi/shixian2/HMM.py0.0019524748085880934 0.21041489951972903 0.09239969242123368[[0.1 0.16 0.28 ][0.077 0.1104 0.0606 ][0.04187 0.035512 0.052836 ][0.0210779 0.02518848 0.01382442][0.01043019 0.01257429 0.00548198][0.00504189 0.00400711 0.00586943][0.00244848 0.00286366 0.00157461][0.00119913 0.00095756 0.00146621][0.00059004 0.00046339 0.00089905]][[0.00364823 0.00395888 0.00340797][0.0079628 0.00754294 0.00835972][0.01479809 0.01691621 0.01385703][0.03194863 0.03446953 0.02971769][0.07018853 0.06518296 0.07310601][0.132638 0.140463 0.122819 ][0.2939 0.2588 0.3123 ][0.54 0.49 0.57 ][1. 1. 1. ]]2.667167999999999e-06 [3, 3, 3, 3, 2, 3, 2, 3, 3]1.0490981222700444e-07 [3, 1, 3, 1, 1, 3, 1, 3, 3]Process finished with exit code 0

總結

以上是生活随笔為你收集整理的统计学习方法第十章作业：HMM模型—概率计算问题、Baum-Welch学习算法、维特比预测算法 代码实现的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。