# viterbi algorithm(pi,A,B)
import numpy as np
# -*- codeing:utf-8 -*-
__author__ = 'Vera'
Hmm_ob= [0, 1,2]
Hmm_hdn=['runny', 'rainy']
A=np.array([[0.7, 0.3], [0.4, 0.6]])
B=np.array([[0.1, 0.4, 0.5], [0.6, 0.3, 0.1]])
Hmm_pi=[0.6, 0.4]
real_ob= [0, 1, 2,1,0,0,2,1,1,2]
def HMM_compute (obser):
N=len(obser) # N means the time
max_sita=np.zeros((len(Hmm_hdn),len(obser)))#This is the max seta for each state
# Each max state should come from some
max_index = np.zeros((len(Hmm_hdn),len(obser)),dtype=int)
path=np.array(len(obser),dtype=int)
prob=np.zeros(len(Hmm_pi))
for i in range(len(Hmm_hdn)):
max_sita[i][0]=Hmm_pi[i]*B[i][obser[0]]
max_index[i][0]=i
for k in range(1,N):
for i in range(len(Hmm_pi)):
for j in range(len(Hmm_hdn)):
prob[j]=max_sita[j][k-1]*A[j][i]*B[i][obser[k]]
a=np.where(prob == np.max(prob))
max_sita[i][k]=np.max(prob)
max_index[i][k]=a[0]
# get path
path=np.zeros(N)
a = np.where(max_sita[:,N-1] == np.max(max_sita[:,N-1]))
path[N - 1] = a[0]
for i in range(N - 1):
aa = int(path[N - 1 - i])
path[N - 2 - i] = max_index[aa][N - 1 - i]
return path
aaa=HMM_compute(real_ob)
for i in range(len(aaa)):
print(Hmm_hdn[int(aaa[i])])
参考了部分代码,但是找不到了。。。