HMM：前向算法实例

看本篇文章，假设你已经知道HMM中的前向算法相关概念
如果不知道，推荐先学习HMM学习最佳范例中相关文章

这个问题来自于HMM学习最佳范例五：前向算法5
只不过再手动算一下，加深一下自己的理解

已知隐马尔科夫模型如下：

1、隐藏状态 (天气)：Sunny，Cloudy，Rainy；
2、观察状态（海藻湿度）：Dry，Dryish，Damp，Soggy；
3、初始状态概率： Sunny（0.63）， Cloudy（0.17）， Rainy（0.20）；
4、状态转移矩阵：

5、混淆矩阵：

M=4(海藻湿度状态数)
N=3(天气状态数）
T=3(观察天数）

$$状态转移矩阵:A= \begin{bmatrix} 0.500 & 0.375 & 0.125 \ 0.250 & 0.125 & 0.625 \ 0.250 & 0.375 & 0.375 \ \end{bmatrix}$$
$$混淆矩阵:B= \begin{bmatrix} 0.60 & 0.20 & 0.15 & 0.05 \ 0.25 & 0.25 & 0.25 & 0.25 \ 0.05 & 0.10 & 0.35 & 0.50 \ \end{bmatrix}$$
$$初始向量:\Pi = \begin{bmatrix} 0.63 &0.17& 0.20\ \end{bmatrix}$$
$$观察序列1,3,4 (“Dry,Damp,Soggy”)$$

计算

用一个$\alpha =N\times T$的矩阵来保存局部概率
$\alpha_{ij}表示t=i的时候，a=j（例如j=1表示天气为“sunny”）$的概率，

1.先计算t=1的局部概率，我们需要用到的初始状态向量$\Pi$，
由于第一天的海藻湿度为1（“Dry”） ，的第一列:
$\alpha_{11}=\Pi_1*B_{11}=0.60*0.63=0.378$
$\alpha_{21}=\Pi_2*B_{21}=0.25*0.17=0.0425$
$\alpha_{31}=\Pi_3*B_{31}=0.05*0.20=0.01$

2.计算t=2的局部概率，当计算t>1的时候，会用到t-1的局部概率
以计算$\alpha_{21}$为例子，首先我们要计算，从t=1的天气，通过状态转移矩阵，得出的通过第一天的天气推测出第二天是“sunny”概率：
$P(第二天sunny)=P(第一天sunny)\times P(第二天sunny|第一天sunny)+P(第一天cloudy) \times P(第二天sunny|第一天cloudy)+P(第一天rainy) \times P(第二天sunny|第一天rainy)$
$=\alpha_{11} \times A_{11}+\alpha_{21} \times A_{21}+\alpha_{31} \times A_{31}=0.202125$
然后，再根据第二天的海藻状态为3（“Damp”），得出

$$\alpha_{21}=P(第二天sunny) \times P(天气sunny|海藻Damp)= 0.202125 \times 0.150000 = 0.030319$$
$\alpha_{22}=(\alpha_{11} \times A_{12}+\alpha_{22} \times

python实现如下：

t=3
m=4
n=3
pi=[0.63 ,0.17 ,0.2]
o=[1,3,4]
a=[
[0.5,0.375 ,0.125],
[0.25 ,0.125 ,0.625],
[0.25 ,0.375 ,0.375]]

b=[
[0.6 ,0.2 ,0.15 ,0.05],
[0.25 ,0.25 ,0.25 ,0.25],
[0.05 ,0.1 ,0.35 ,0.5]]


alpha=[[0]*t for x in range(t)]

for i in range(t):
    alpha[0][i]=round(pi[i]*b[i][o[0]-1],10)
    print "alpha[0][%d]=%f*%f=%f" %(i,pi[i],b[i][o[0]-1],alpha[0][i])

for i in range(1,t):
    for j in range(0,n):
        sum=0;
        symStr=str();
        for k in range(0,n):
            sum+=alpha[i-1][k]*a[k][j]
            symStr+="%f*%f+" % (alpha[i-1][k],a[k][j]) 
        print "sum=%s=%f" % (symStr[:-1],sum)
        alpha[i][j]=sum*b[j][o[i]-1]
        print "alpha[%d][%d]=%f*%f=%f" % (i,j,sum,b[j][o[i]-1],alpha[i][j])
prob=0;
for i in range(n):
    prob+=alpha[t-1][i]
print "prib=%f" % prob

输出结果：

alpha[0][0]=0.630000*0.600000=0.378000
alpha[0][1]=0.170000*0.250000=0.042500
alpha[0][2]=0.200000*0.050000=0.010000
sum=0.378000*0.500000+0.042500*0.250000+0.010000*0.250000=0.202125
alpha[1][0]=0.202125*0.150000=0.030319
sum=0.378000*0.375000+0.042500*0.125000+0.010000*0.375000=0.150812
alpha[1][1]=0.150812*0.250000=0.037703
sum=0.378000*0.125000+0.042500*0.625000+0.010000*0.375000=0.077563
alpha[1][2]=0.077563*0.350000=0.027147
sum=0.030319*0.500000+0.037703*0.250000+0.027147*0.250000=0.031372
alpha[2][0]=0.031372*0.050000=0.001569
sum=0.030319*0.375000+0.037703*0.125000+0.027147*0.375000=0.026263
alpha[2][1]=0.026263*0.250000=0.006566
sum=0.030319*0.125000+0.037703*0.625000+0.027147*0.375000=0.037534
alpha[2][2]=0.037534*0.500000=0.018767
prib=0.026901