1、隐马尔可夫模型中的 Viterbi 算法先用一句话来简单描述一下：给出一个观测序列 o1,o2,o3 .，我们希望找到观测序列背后的隐藏状态序列 s1, s2, s3, .；Viterbi 以它的发明者名字命名，正是这样一种由动态规划的方法来寻找出现概率最大的隐藏状态序列（被称为 Viterbi 路径）的算法。这里需要抄一点有关隐马可夫序列（HMM，Hidden Markov Model）的书页来解释一下观测序列和隐藏状态序列。首先从最简单的离散 Markov 过程入手，我们知道，Markov 随机过程具有如下的性质：在任意时刻，从当前状态转移到下一个状态的概率与当前状态之前的那些状态没有关系
2、。所以，我们可以用一个状态转移概率矩阵来描述它。假设我们有 n 个离散状态 S1, S2,Sn，我们可以构造一个矩阵 A，矩阵中的元素 aij 表示从当前状态 Si 下一时刻迁移到 Sj 状态的概率。但是在很多情况下，Markov 模型中的状态是我们观察不到的。例如，容器与彩球的模型：有若干个容器，每个容器中按已知比例放入各色的彩球（这样，选择了容器后，我们可以用概率来预测取出各种彩球的可能性） ；我们做这样的实验，实验者从容器中取彩球先选择一个容器，再从中抓出某一个球，只给观察者看球的颜色；这样，每次取取出的球的颜色是可以观测到的，即 o1, o2,，但是每次选择哪个容器是不暴露给观察者的，
3、容器的序列就组成了隐藏状态序列 S1, S2,Sn。这是一个典型的可以用HMM 描述的实验。HMM 有几个重要的任务，其中之一就是期望通过观察序列来猜测背后最有可能的隐藏序列。在上面的例子中，就是找到我们在实验中最有可能选择到的容器序列。Viterbi 正是用来解决这个问题的算法。HMM 另外两个任务是：a) 给定一个 HMM，计算一个观测序列出现的可能性；b)已知一个观测序列，HMM 参数不定，如何优化这些参数使得观测序列的出现概率最大。解决前一个问题可以用与 Viberbi 结构非常类似的 Forward 算法来解决（实际上在下面合二为一） ，而后者可以用 Baum-Welch/EM 算法
4、来迭代逼近。从 Wiki 上抄一个例子来说明 Viterbi 算法。假设你有一个朋友在外地，每天你可以通过电话来了解他每天的活动。他每天只会做三种活动之一Walk, Shop, Clean。你的朋友从事哪一种活动的概率与当地的气候有关，这里，我们只考虑两种天气Rainy, Sunny 。我们知道，天气与运动的关系如下：Rainy SunnyWalk 0.1 0.6Shop 0.4 0.3Clean 0.5 0.1例如，在下雨天出去散步的可能性是 0.1。而天气之前互相转换的关系如下， （从行到列）Rainy SunnyRainy 0.7 0.3Sunny 0.4 0.6例如，从今天是晴天而明天
5、就开始下雨的可能性是 0.4 。同时为了求解问题我们假设初始情况：通话开始的第一天的天气有 0.6 的概率是 Rainy，有 0.4 概率是 Sunny。OK ，现在的问题是，如果连续三天，你发现你的朋友的活动是： Walk-Shop-Clean；那么，如何判断你朋友那里这几天的天气是怎样的？解决这个问题的 python 伪代码如下：def forward_viterbi(y, X, sp, tp, ep):T = for state in X:# prob. V. path V. prob.Tstate = (spstate, state, spstate)for output in y:U
6、 = for next_state in X:total = 0argmax = Nonevalmax = 0for source_state in X:(prob, v_path, v_prob) = Tsource_statep = epsource_stateoutput * tpsource_statenext_stateprob *= pv_prob *= ptotal += probif v_prob valmax:argmax = v_path + next_statevalmax = v_probUnext_state = (total, argmax, valmax)T =
7、U# apply sum/max to the final states:total = 0argmax = Nonevalmax = 0for state in X:(prob, v_path, v_prob) = Tstatetotal += probif v_prob valmax:argmax = v_pathvalmax = v_probreturn (total, argmax, valmax)几点说明：1、算法对于每一个状态要记录一个三元组：(prob, v_path, v_prob)，其中，prob 是从开始状态到当前状态所有路径（不仅仅 是最有可能的 viterbi 路径）的
8、概率加在一起的结果（作为算法附产品，它可以输出一个观察序列在给定 HMM 下总的出现概率，即 forward 算法的输出） ，v_path 是从开始状态一直到当前状态的 viterbi 路径，v_prob 则是该路径的概率。 2、算法开始，初始化 T （T 是一个 Map，将每一种可能状态映射到上面所说的三元组上） 3、三重循环，对每个一活动 y，考虑下一步每一个可能的状态 next_state，并重新计算若从 T 中的当前状态 state 跃迁到 next_state 概率会有怎样的变化。跃迁主要考虑天气转移（tpsource_statenext_state ）与该天气下从事某种活动（eps
9、ource_stateoutput）的联合概率。所有下一步状态考虑完后，要从 T 中找出最优的选择 viterbi 路径即概率最大的 viterbi 路径，即上面更新 Map U 的代码 Unext_state = (total, argmax, valmax)。 4、算法最后还要对 T 中的各种情况总结，对 total 求和，选择其中一条作为最优的 viterbi 路径。 5、算法输出四个天气状态，这是因为，计算第三天的概率时，要考虑天气转变到下一天的情况。 6、通过程序的输出可以帮助理解这一过程： observation=Walknext_state=Sunnystate=Sunnyp=0
10、.36triple=(0.144,Sunny-,0.144)state=Rainyp=0.03triple=(0.018,Rainy-,0.018)Update USunny=(0.162,Sunny-Sunny-,0.144)next_state=Rainystate=Sunnyp=0.24triple=(0.096,Sunny-,0.096)state=Rainyp=0.07triple=(0.042,Rainy-,0.042)Update URainy=(0.138,Sunny-Rainy-,0.096)observation=Shopnext_state=Sunnystate=Sunn
11、yp=0.18triple=(0.02916,Sunny-Sunny-,0.02592)state=Rainyp=0.12triple=(0.01656,Sunny-Rainy-,0.01152)Update USunny=(0.04572,Sunny-Sunny-Sunny-,0.02592)next_state=Rainystate=Sunnyp=0.12triple=(0.01944,Sunny-Sunny-,0.01728)state=Rainyp=0.28triple=(0.03864,Sunny-Rainy-,0.02688)Update URainy=(0.05808,Sunny
12、-Rainy-Rainy-,0.02688)observation=Cleannext_state=Sunnystate=Sunnyp=0.06triple=(0.0027432,Sunny-Sunny-Sunny-,0.0015552)state=Rainyp=0.15triple=(0.008712,Sunny-Rainy-Rainy-,0.004032)Update USunny=(0.0114552,Sunny-Rainy-Rainy-Sunny-,0.004032)next_state=Rainystate=Sunnyp=0.04triple=(0.0018288,Sunny-Sun
13、ny-Sunny-,0.0010368)state=Rainyp=0.35triple=(0.020328,Sunny-Rainy-Rainy-,0.009408)Update URainy=(0.0221568,Sunny-Rainy-Rainy-Rainy-,0.009408)final triple=(0.033612,Sunny-Rainy-Rainy-Rainy-,0.009408)所以，最终的结果是，朋友那边这几天最可能的天气情况是 Sunny-Rainy-Rainy-Rainy，它有0.009408 的概率出现。而我们算法的另一个附带的结论是，我们所观察到的朋友这几天的活动序列：
14、Walk-Shop-Clean 在我们的隐马可夫模型之下出现的总概率是 0.033612。附：c+ 主要代码片断void forward_viterbi(const vector map map / initializationInitParameters();map T;for (vector:iterator it = states.begin();it != states.end();+it).viterbi_triple_t foo;foo.prob = sp*it;foo.vpath.push_back(*it);foo.vprob = sp*it;T*it = foo;map U;
15、double total = 0;vector argmax;double valmax = 0;double p = 0;for (vector:const_iterator itob = ob.begin(); itob != ob.end(); +itob).cout :iterator itNextState = states.begin();itNextState != states.end();+itNextState).cout :iterator itSrcState = states.begin();itSrcState != states.end();+itSrcState
16、).cout valmax).foo.vpath.push_back(*itNextState);argmax = foo.vpath;valmax = foo.vprob;U*itNextState = viterbi_triple_t(total, argmax, valmax);cout :iterator itState = states.begin();itState != states.end();+itState).viterbi_triple_t foo = T*itState;total += foo.prob;if (foo.vprob valmax).argmax.swap(foo.vpath);valmax = foo.vprob;vtriple.prob = total;vtriple.vpath = argmax;vtriple.vprob = valmax;cout “final triple=“ vtriple endl;
