在HMM这个翻译系列的原文中,作者举了一个前向算法的交互例子,这也是这个系列中比较出彩的地方,但是,在具体运行这个例子的时候,却发现其似乎有点问题。
先说一下如何使用这个交互例子,运行时需要浏览器支持java,我用的是firefox。首先在Set按钮前面的对话框里上观察序列,如“Dry,Damp, Soggy” 或“Dry Damp Soggy”,观察符号间用逗号或空格隔开;然后再点击Set按钮,这样就初始化了观察矩阵;如果想得到一个总的结果,即Pr(观察序列|隐马尔科夫模型),就点旁边的Run按钮;如果想一步一步观察计算过程,即每个节点的局部概率,就单击旁边的Step按钮。
原文交互例子(即天气这个例子)中所定义的已知隐马尔科夫模型如下:
1、隐藏状态 (天气):Sunny,Cloudy,Rainy;
2、观察状态(海藻湿度):Dry,Dryish,Damp,Soggy;
3、初始状态概率: Sunny(0.63), Cloudy(0.17), Rainy(0.20);
4、状态转移矩阵:
weather today
Sunny Cloudy Rainy
weather Sunny 0.500 0.375 0.125
yesterday Cloudy 0.250 0.125 0.625
Rainy 0.250 0.375 0.375
5、混淆矩阵:
observed states
Dry Dryish Damp Soggy
Sunny 0.60 0.20 0.15 0.05
hidden Cloudy 0.25 0.25 0.25 0.25
states Rainy 0.05 0.10 0.35 0.50
为了UMDHMM也能运行这个例子,我们将上述天气例子中的隐马尔科夫模型转化为如下的UMDHMM可读的HMM文件weather.hmm:
--------------------------------------------------------------------
M= 4
N= 3
A:
0.500 0.375 0.125
0.250 0.125 0.625
0.250 0.375 0.375
B:
0.60 0.20 0.15 0.05
0.25 0.25 0.25 0.25
0.05 0.10 0.35 0.50
pi:
0.63 0.17 0.20
--------------------------------------------------------------------
在运行例子之前,如果读者也想观察每一步的运算结果,可以将umdhmm-v1.02目录下forward.c中的void Forward(…)函数替换如下:
--------------------------------------------------------------------
void Forward(HMM *phmm, int T, int *O, double **alpha, double *pprob)
{
int i, j; /* state indices */
int t; /* time index */
double sum; /* partial sum */
/* 1. Initialization */
for (i = 1; i <= phmm->N; i++)
{
alpha[1][i] = phmm->pi[i]* phmm->B[i][O[1]];
printf( "a[1][%d] = pi[%d] * b[%d][%d] = %f * %f = %f\\n",i, i, i, O[i], phmm->pi[i], phmm->B[i][O[1]], alpha[1][i] );
}
/* 2. Induction */
for (t = 1; t < T; t++)
{
for (j = 1; j <= phmm->N; j++)
{
sum = 0.0;
for (i = 1; i <= phmm->N; i++)
{
sum += alpha[t][i]* (phmm->A[i][j]);
printf( "a[%d][%d] * A[%d][%d] = %f * %f = %f\\n", t, i, i, j, alpha[t][i], phmm->A[i][j], alpha[t][i]* (phmm->A[i][j]));
printf( "sum = %f\\n", sum );
}
alpha[t+1][j] = sum*(phmm->B[j][O[t+1]]);
printf( "a[%d][%d] = sum * b[%d][%d]] = %f * %f = %f\\n",t+1, j, j, O[t+1], sum, phmm->B[j][O[t+1]], alpha[t+1][j] );
}
}
/* 3. Termination */
*pprob = 0.0;
for (i = 1; i <= phmm->N; i++)
{
*pprob += alpha[T][i];
printf( "alpha[%d][%d] = %f\\n", T, i, alpha[T][i] );
printf( "pprob = %f\\n", *pprob );
}
}
--------------------------------------------------------------------
替换完毕之后,重新“make clean”,“make all”,这样新的testfor可执行程序就可以输出前向算法每一步的计算结果。
现在我们就用testfor来运行原文中默认给出的观察序列“Dry,Damp,Soggy”,其所对应的UMDHMM可读的观察序列文件test1.seq:
--------------------------------------------------------------------
T=3
1 3 4
--------------------------------------------------------------------
好了,一切准备工作就绪,现在就输入如下命令:
testfor weather.hmm test1.seq > result1
result1就包含了所有的结果细节:
--------------------------------------------------------------------
Forward without scaling
a[1][1] = pi[1] * b[1][1] = 0.630000 * 0.600000 = 0.378000
a[1][2] = pi[2] * b[2][3] = 0.170000 * 0.250000 = 0.042500
a[1][3] = pi[3] * b[3][4] = 0.200000 * 0.050000 = 0.010000
…
pprob = 0.026901
log prob(O| model) = -3.615577E+00
prob(O| model) = 0.026901
…
--------------------------------------------------------------------
黑体部分是最终的观察序列的概率结果,即本例中的Pr(观察序列|HMM) = 0.026901。
但是,在原文中点Run按钮后,结果却是:Probability of this model = 0.027386915。
这其中的差别到底在哪里?我们来仔细观察一下中间运行过程:
在初始化亦t=1时刻的局部概率计算两个是一致的,没有问题。但是,t=2时,在隐藏状态“Sunny”的局部概率是不一致的。英文原文给出的例子的运行结果是:
Alpha = (((0.37800002*0.5) + (0.0425*0.375) + (0.010000001*0.125)) * 0.15) = 0.03092813
而UMDHMM给出的结果是:
--------------------------------------------------------------------
a[1][1] * A[1][1] = 0.378000 * 0.500000 = 0.189000
sum = 0.189000
a[1][2] * A[2][1] = 0.042500 * 0.250000 = 0.010625
sum = 0.199625
a[1][3] * A[3][1] = 0.010000 * 0.250000 = 0.002500
sum = 0.202125
a[2][1] = sum * b[1][3]] = 0.202125 * 0.150000 = 0.030319
--------------------------------------------------------------------
区别就在于状态转移概率的选择上,原文选择的是状态转移矩阵中的第一行,而UMDHMM选择的则是状态转移矩阵中的第一列。如果从原文给出的状态转移矩阵来看,第一行代表的是从前一时刻的状态“Sunny”分别到当前时刻的状态“Sunny”,“Cloudy”,“Rainy”的概率;而第一列代表的是从前一时刻的状态“Sunny”,“Cloudy”,“Rainy”分别到当前时刻状态“Sunny”的概率。这样看来似乎原文的计算过程有误,读者不妨多试几个例子看看,前向算法这一章就到此为止了。
未完待续:维特比算法1
本翻译系列原文请参考:http://www.comp.leeds.ac.uk/roger/HiddenMarkovModels/html_dev/main.html
注:原创文章,转载请注明出处“我爱自然语言处理”:www.52nlp.cn
本文链接地址:https://www.52nlp.cn/hmm-learn-best-practices-five-forward-algorithm-5
今天尝试按照关键字搜索了一下你的网站,确实很有效果,做自己想做并且喜欢做的事情,你的坚持让我很感动。
[回复]
admin 回复:
3 8 月, 2009 at 18:59
恩,你也一样!
[回复]
真的很少能关于hmm这么全面的网站了,谢谢了
[回复]
admin 回复:
14 8 月, 2009 at 21:28
谢谢,欢迎常来看看!
[回复]
确实很全面啊
知道关于HSMM的这方面文章不?
[回复]
admin 回复:
20 10 月, 2009 at 19:58
呵呵,这个不清楚!
[回复]
前辈,多谢你的辛苦翻译,让我对隐马有了更深层次的了解,最近我用HMM来对时间序列建模,有两个问题想向你请教一下,不知道能否留个联系方式?
我的QQ:627431239,多谢~
[回复]
可以发邮件给我,52nlpcn#gmail.com
[回复]
我先看到英文原文,实现了算法,对比了结果,也觉得那个转移矩阵不对,正在怀疑之时,看到了你这篇文章,谢谢!
[回复]
52nlp 回复:
30 12 月, 2009 at 19:58
不客气!尽信书不如无书吧!
[回复]
请问大侠:HMM用于时间序列建立模型,特别是考虑到输入对输出结果的影响,如何解决连续性\小样本和控制输入问题?谢谢
[回复]
52nlp 回复:
18 1 月, 2010 at 19:19
不好意思,这个问题我不太清楚,我对HMM的理解还仅限于离散这一块儿。
[回复]
请问
Forward without scaling
a[1][1] = pi[1] * b[1][1] = 0.630000 * 0.600000 = 0.378000
a[1][2] = pi[2] * b[2][3] = 0.170000 * 0.250000 = 0.042500
a[1][3] = pi[3] * b[3][4] = 0.200000 * 0.050000 = 0.010000
是不是错了?是不是应:
a[1][1] = pi[1] * b[1][1] = 0.630000 * 0.600000 = 0.378000
a[1][2] = pi[2] * b[2][1] = 0.170000 * 0.250000 = 0.042500
a[1][3] = pi[3] * b[3][1] = 0.200000 * 0.050000 = 0.010000
[回复]
52nlp 回复:
14 1 月, 2011 at 13:33
好像没错吧,1,3,4代表的是观察序列。不过好久了,我都有点忘了,呵呵。
[回复]
xiao 回复:
22 12 月, 2011 at 14:07
这个问题我认为博主有点小笔误, yxs 更正的是正确的,即,正确结果如下:
a[1][1] = pi[1] * b[1][1] = 0.630000 * 0.600000 = 0.378000
a[1][2] = pi[2] * b[2][1] = 0.170000 * 0.250000 = 0.042500
a[1][3] = pi[3] * b[3][1] = 0.200000 * 0.050000 = 0.010000
b的第二个角标都是1。因为观察序列1,3,4表示在t=1时的观察值为1,在t=2时的观察值为3,在t=3时的观察值为4,所以在求解t=1时的局部概率时,应该分别取混淆矩阵中的b[1][1], b[2][1],b[3][1]。另外,在翻译混淆矩阵的时候,个人认为混合矩阵稍显恰当。以上仅是个人拙见,请博主指点。再次感谢博主辛勤的努力,您的博客对我很有帮助。
[回复]
52nlp 回复:
27 12 月, 2011 at 23:39
谢谢质疑,主要这个结果是程序给出的,我没有手写,所以不太确定,都有点糊涂了。
[回复]
xiao 回复:
28 12 月, 2011 at 09:56
博主好,谢谢您的认真回复,我的理解是这样的:您程序中的initialization部分,printf( “a[1][%d] = pi[%d] * b[%d][%d] = %f * %f = %f\n”,i, i, i, O[i], phmm->pi[i], phmm->B[i][O[1]], alpha[1][i] );
中的O[i]写的有误,应该是O[1],因为现在是计算t=1时的局部概率,您在这句中后面用的就是phmm->B[i][O[1]],也就是表明现在要乘的是O[1]对应的元素。
同时,您这里的小笔误不影响计算的正确性,因为在输出语句的上一句,也就是真正局部概率的那一句alpha[1][i] = phmm->pi[i]* phmm->B[i][O[1]];您正确使用了B[i][O[1]]。以上仅是个人拙见,请大家讨论。
再次感谢楼主认真细致的学术品格和细心翻译,我学到了很多。
52nlp 回复:
29 12 月, 2011 at 13:53
最近实在太忙了,我元旦的时候再看一下这个问题,谢谢你这么细致的检查!
前辈,我想问一下,UMDHMM前向算法中 初始向量都是1 怎么解释呢?代码如下:
/* 1. Initialization */
for (i = 1; i N; i++)
beta[T][i] = 1.0;
还有关于带有加权参量的向前算法,这个加权变量scale[t] 是一个时间变量,是不是可以解决部分状态转移和时间无关的缺点?
[回复]
52nlp 回复:
29 4 月, 2011 at 22:44
beta[T][i] = 1.0 这个应该是后向算法里的吧?
关于scale[t],抱歉我不太清楚是否能解决这个问题!
[回复]
simon 回复:
3 5 月, 2011 at 15:52
哦 是的 猪头了一回
[回复]
真不错 连程序都出来了
[回复]
英文原文确实有误,状态转移矩阵的选择不对
我用穷举法验证了一下
P(obervation|HMM)= 0.026901406
[回复]
我想用来做手势识别,但是训练以后在计算前向概率的时候概率值很低很低,一般一个序列15的长度计算得到的概率值都是-30多,不知道这样是否正常?谢谢了
[回复]
52nlp 回复:
26 11 月, 2011 at 10:15
不清楚,只要有可比性,应该是正常的吧。
[回复]
袁军 回复:
6 7 月, 2012 at 11:14
我也碰到了相同的问题,不知道您已经解决了吗?
[回复]
看到umdhmm的hmm.c里面对于M的定义是 离散的观测符号的个数,请问:
如果观测的概率分布不是离散的,而是高斯分布,还能不能使用umdhmm?
Thanks.
[回复]
52nlp 回复:
5 12 月, 2011 at 23:26
不清楚,不过即使是连续的分布,对于M也只是其中的一些观察值,代表的是一些抽样的点吧。
[回复]
博主您好,请教您两个问题,1.隐马尔科夫模型多观察序列的概率应该怎么求?2.给定一个隐马尔科夫,用前向算法求一个观察序列的概率,但这个观察序列的T=100,算到后面alpha的值都为零了,是怎么回事?望博主赐教!
[回复]
52nlp 回复:
9 7 月, 2012 at 11:32
可以考虑用scaling 的方法,建议看这个文章的第6节:HMM Scaling
http://www.cs.sjsu.edu/~stamp/RUA/HMM.pdf
[回复]
袁军 回复:
14 7 月, 2012 at 20:31
谢谢博主,这是一种数据下溢问题,需增加比例因子,程序中我调了好久也没弄好,请问博主有程序这方面资料吗?
[回复]
52nlp 回复:
17 7 月, 2012 at 15:59
抱歉,不太清楚。
老师,我按照你说的执行了下,为什么在linux终端会出现:
yuemeng@yuemeng-desktop:~/下载/umdhmm/umdhmm-v1.02$ ./testfor weather.hmm test1.seq
bash: ./testfor: cannot execute binary file
yuemeng@yuemeng-desktop:~/下载/umdhmm/umdhmm-v1.02$
这样的提示呢,请问那地方出错了呢,上面有我程序的路径,请指导下啊
非常感谢
[回复]
sunmoon 回复:
9 5 月, 2016 at 16:30
我遇到和你同样的问题,请问你是怎么解决的
[回复]
嗯,这个问题解决了啊,还是非常感谢你的文章,关于HMM最好的文章
[回复]
感谢博主的文章. 我刚看到这里,感觉文章中的:
Alpha = (((0.37800002*0.5) + (0.0425*0.375) + (0.010000001*0.125)) * 0.15) = 0.03092813 应该是对的,而程序里的(0.5,0.25,0.25)是错误的. 因为是从其他隐藏状态向第一个隐藏状态转移的概率,所以使用第一列来计算,而不是第一行. 不知道我理解的对不对.
[回复]
mike98 回复:
7 3 月, 2013 at 07:27
Alpha = (((0.37800002*0.5) + (0.0425*0.375) + (0.010000001*0.125)) * 0.15) = 0.03092813
肯定是错的。0.375在后面给出的数据中应该为0.325。另外,从转移矩阵上看,行为昨天,列为今天,昨天到今天显然应该取列向量。
[回复]
我用楼主的程序和数据实验,结果居然全是0。不知道怎么回事?
[回复]
William 回复:
16 10 月, 2012 at 01:21
在Mac OS上似乎ReadHMM()无法正确的得到hmm中的数据。验证了一个hmm.pi[1] = 0.000000
[回复]
52nlp 回复:
16 10 月, 2012 at 13:18
是在mac上编译的吗?
[回复]
William 回复:
16 10 月, 2012 at 15:58
是的,我简单的跟了一下,应该是读取hmm数据时,数据类型(double)的问题,可能原来的实现只是在Linux平台。
William 回复:
16 10 月, 2012 at 22:32
好奇怪。。ReadHMM() 中的fscanf(fp, "%lf", &(phmm->A[i][j])); 居然不work。
搞定了,原来是我的weather.hmm格式有问题。多了个tab。。。
[回复]
感谢对原文中错误的指正,当时看到原文就觉得很疑惑感觉他是不是写错了。但是限于自己水平不能确定。
[回复]
原文确实是有些问题,有可能是老外看transition matrix的顺序和我们不一样
[回复]
我在linue下进入目录,make all ,显示结果‘all’ is up to data
[回复]
豹豹 回复:
8 4 月, 2016 at 16:24
make clean
make & make all
[回复]
你好,楼主,,printf( “a[1][%d] = pi[%d] * b[%d][%d] = %f * %f = %f\n”,i, i, i, O[i], phmm->pi[i], phmm->B[i][O[1]], alpha[1][i] );,,我在linux下运行这个时,为什么报错说的a[][]和b[][]没有声明
[回复]
老师,你好,
forward.c: In function ‘Forward’:
forward.c:28: error: stray ‘\342’ in program
forward.c:28: error: stray ‘\200’ in program
forward.c:28: error: stray ‘\234’ in program
forward.c:28: error: ‘a’ undeclared (first use in this function)
forward.c:28: error: (Each undeclared identifier is reported only once
forward.c:28: error: for each function it appears in.)
forward.c:28: error: expected expression before ‘%’ token
forward.c:28: error: ‘pi’ undeclared (first use in this function)
forward.c:28: error: expected expression before ‘%’ token
forward.c:28: error: ‘b’ undeclared (first use in this function)
forward.c:28: error: expected expression before ‘%’ token
forward.c:28: error: expected expression before ‘%’ token
forward.c:28: error: expected expression before ‘%’ token
forward.c:28: error: stray ‘\’ in program
forward.c:28: error: stray ‘\’ in program
forward.c:28: error: stray ‘\342’ in program
forward.c:28: error: stray ‘\200’ in program
forward.c:28: error: stray ‘\235’ in program
make: *** [forward.o] Error 1
这些报错,怎么办
[回复]
52nlp 回复:
28 4 月, 2016 at 20:01
你用的什么编译器?
[回复]