首先需要说明的是，本节不是这个系列的翻译，而是作为前向算法这一章的补充，希望能从实践的角度来说明前向算法。除了用程序来解读hmm的前向算法外，还希望将原文所举例子的问题拿出来和大家探讨。
文中所举的程序来自于UMDHMM这个C语言版本的HMM工具包，具体见《几种不同程序语言的HMM版本》(http://www.52nlp.cn/several-different-programming-language-version-of-hmm)。先说明一下UMDHMM这个包的基本情况，在linux环境下，进入umdhmm-v1.02目录，“make all”之后会产生4个可执行文件，分别是：
genseq: 利用一个给定的隐马尔科夫模型产生一个符号序列（Generates a symbol sequence using the specified model sequence using the specified model）
testfor: 利用前向算法计算log Prob(观察序列| HMM模型)（Computes log Prob(observation|model) using the Forward algorithm.）
testvit: 对于给定的观察符号序列及HMM，利用Viterbi 算法生成最可能的隐藏状态序列（Generates the most like state sequence for a given symbol sequence, given the HMM, using Viterbi）
esthmm: 对于给定的观察符号序列，利用BaumWelch算法学习隐马尔科夫模型HMM（Estimates the HMM from a given symbol sequence using BaumWelch）。
这些可执行文件需要读入有固定格式的HMM文件及观察符号序列文件，格式要求及举例如下：
HMM 文件格式：
M= number of symbols
N= number of states
A:
a11 a12 ... a1N
a21 a22 ... a2N
. . . .
. . . .
. . . .
aN1 aN2 ... aNN
B:
b11 b12 ... b1M
b21 b22 ... b2M
. . . .
. . . .
. . . .
bN1 bN2 ... bNM
pi:
pi1 pi2 ... piN
HMM文件举例：
M= 2
N= 3
A:
0.333 0.333 0.333
0.333 0.333 0.333
0.333 0.333 0.333
B:
0.5 0.5
0.75 0.25
0.25 0.75
pi:
0.333 0.333 0.333
观察序列文件格式：
T=seqence length
o1 o2 o3 . . . oT
观察序列文件举例：
T= 10
1 1 1 1 2 1 2 2 2 2
对于前向算法的测试程序testfor来说，运行：
testfor model.hmm（HMM文件） obs.seq（观察序列文件）
就可以得到观察序列的概率结果的对数值，这里我们在testfor.c的第58行对数结果的输出下再加一行输出：
fprintf(stdout, "prob(O| model) = %f\n", proba);
就可以输出运用前向算法计算观察序列所得到的概率值。至此，所有的准备工作已结束，接下来，我们将进入具体的程序解读。
首先，需要定义HMM的数据结构，也就是HMM的五个基本要素，在UMDHMM中是如下定义的（在hmm.h中）：
typedef struct
{
int N; /* 隐藏状态数目;Q={1,2,...,N} /
int M; / 观察符号数目; V={1,2,...,M}/
double A; / 状态转移矩阵A[1..N][1..N]. a[i][j] 是从t时刻状态i到t+1时刻状态j的转移概率 /
double B; / 混淆矩阵B[1..N][1..M]. b[j][k]在状态j时观察到符合k的概率。/
double pi; / 初始向量pi[1..N]，pi[i] 是初始状态概率分布 /
} HMM;
前向算法程序示例如下（在forward.c中）：
/
函数参数说明：
phmm：已知的HMM模型；T：观察符号序列长度；
O：观察序列；alpha：局部概率；pprob：最终的观察概率
*/
void Forward(HMM *phmm, int T, int *O, double **alpha, double pprob)
{
int i, j; / 状态索引 /
int t; / 时间索引 /
double sum; /求局部概率时的中间值 */
/* 1. 初始化：计算t=1时刻所有状态的局部概率alpha： /
for (i = 1; i <= phmm->N; i++)
alpha[1][i] = phmm->pi[i] phmm->B[i][O[1]];
/* 2. 归纳：递归计算每个时间点，t=2，... ，T时的局部概率 /
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]);
alpha[t+1][j] = sum*(phmm->B[j][O[t+1]]);
}
}
/* 3. 终止：观察序列的概率等于T时刻所有局部概率之和/
pprob = 0.0;
for (i = 1; i <= phmm->N; i++)
*pprob += alpha[T][i];
}
下一节我将用这个程序来验证英文原文中所举前向算法演示例子的问题。
网友评论