隐马尔可夫HMM中viterbi算法

来源：互联网发布：淘宝图文详情怎么制作编辑：程序博客网时间：2024/06/05 09:34

引言

viterbi算法简化最有可能的天气序列的运算过程，forward算法简化该该观察值的概率。

问题描述

你在中国，你朋友F在美国，F的作息有walk, shop, clean，但这选择跟天气有关，我们又知道Rainy的概率比Sunny的概率大
这是初始概率

这是天气转移矩阵

这是在相应天气下发生相应事的概率分布
然后，F这三天是walk，walk，shop，问最有可能的天气序列

问题分析

同样的，我们先用穷举法来算，即

Sunny Sunny Sunny
Sunny Sunny Rainy
Sunny Rainy Sunny
Sunny Rainy Rainy
Rainy Sunny Sunny
Rainy Sunny Rainy
Rainy Rainy Sunny
Rainy Rainy Rainy

分别求出在各个情况下，{walk，walk，shop}的概率，然后比较概率，最大的就是我们想要的结果，我们发现这跟我们forward的应用场景相似，只不过，forward算法是求和，而viterbi算法是求最大，forward算法是所有路径之和，而viterbi算法只是求一条路径。这也就意味着这两算法可以写在一起。

那我们重点讲讲算法吧，这次我是看wiki上的代码，再写，再比较，发现wiki上的代码有一处甚是费解，因此，加以修改，并一目了然了，我把我的思考罗列如下

1）我的思路是把第一天的情况特殊考虑，可wiki上的代码直接整合在一起了，其意义甚是费解，所以，我把第一天特殊考虑了

2）再者，我发现wiki上的代码，明明要求只需要确定三天的天气却出现四天的天气，估计是第一天没用的

3）我修改后代码跟wiki上代码，forward算法最终结果是一样的，但过程中求得值是不一样的

4）用viterbi求得概率是不同的，但wiki上的代码，把第一天去掉之后，剩下的天气序列跟我的是一样

不管怎样，wiki上的代码整体框架还是非常好的（或许还是我改错了），思路如下：

1）我上一篇自己写forward算法，差不多用二维数组来记录整个过程，其实没必要，只需一维足以，用完一行，替代一行，所以才会有字典T和字典U

2）返回三样，分别是该观察序列的最大概率（即之前用forward求出来的），最可能出现的天气序列，该天气序列的概率，同时每次都更新这三个值

3）三重for循环，第一重是第几天，第二重是第几天晴天还是雨天，第三重是前一天是晴天还是雨天下求得概率，再具体求和，求最大值

我稍加修改的Python代码及结果

 1 def forward_viterbi(obs, states, start_p, trans_p, emit_p): 2     T = {} 3     #for the first one 4     for state in states: 5         T[state] = (start_p[state] * emit_p[state][obs[0]], [state], start_p[state]) 6     for output in obs: 7         U = {} 8         #pass the first one 9         if output == obs[0]:10             continue11         for next_state in states:12             total = 013             argmax = None14             valmax = 015             for source_state in states:16                 #get the previous one17                 (prob, v_path, v_prob) = T[source_state]18                 #caculate the commone one19                 p = trans_p[source_state][next_state] * emit_p[next_state][output]20                 prob *= p21                 v_prob *= p22                 #sum up23                 total += prob24                 #find the max25                 if v_prob > valmax:26                     argmax = v_path + [next_state]27                     valmax = v_prob28             U[next_state] = (total, argmax, valmax)29         T = U30     #sum up & find the max one31     total = 032     argmax = None33     valmax = 034     for state in states:35         (prob, v_path, v_prob) = T[state]36         total += prob37         if v_prob > valmax:38             valmax = v_prob39             argmax = v_path40     return (total, argmax, valmax)41 states = ('Rainy', 'Sunny')42 observations = ('walk', 'shop', 'clean')43 start_probability = {'Rainy': 0.6, 'Sunny': 0.4}44 transition_probability = {45    'Rainy' : {'Rainy': 0.7, 'Sunny': 0.3},46    'Sunny' : {'Rainy': 0.4, 'Sunny': 0.6},47    }48 emission_probability = {49    'Rainy' : {'walk': 0.1, 'shop': 0.4, 'clean': 0.5},50    'Sunny' : {'walk': 0.6, 'shop': 0.3, 'clean': 0.1},51    }52 #A simple example of the using algorithm53 def example():54     return forward_viterbi(observations,55                            states,56                            start_probability,57                            transition_probability,58                            emission_probability)59 print example()60 61 #(0.033611999999999996, ['Rainy', 'Rainy', 'Rainy'], 0.05879999999999999)

wiki代码及结果

 1 def forward_viterbi(obs, states, start_p, trans_p, emit_p): 2     T = {} 3     for state in states: 4         ##          prob.           V. path  V. prob. 5         T[state] = (start_p[state], [state], start_p[state]) 6     for output in obs: 7         print T 8         U = {} 9         for next_state in states:10             total = 011             argmax = None12             valmax = 013             for source_state in states:14                 (prob, v_path, v_prob) = T[source_state]15                 p = emit_p[source_state][output] * trans_p[source_state][next_state]16                 #caculate17                 prob *= p18                 v_prob *= p19                 #the different way to gain20                 total += prob21                 if v_prob > valmax:22                     argmax = v_path + [next_state]23                     valmax = v_prob24             U[next_state] = (total, argmax, valmax)25         T = U26     print T27     ## apply sum/max to the final states:28     total = 029     argmax = None30     valmax = 031     for state in states:32         (prob, v_path, v_prob) = T[state]33         total += prob34         if v_prob > valmax:35             argmax = v_path36             valmax = v_prob37     return (total, argmax, valmax)38    39 states = ('Rainy', 'Sunny')40 observations = ('walk', 'shop', 'clean')41 start_probability = {'Rainy': 0.6, 'Sunny': 0.4}42 transition_probability = {43    'Rainy' : {'Rainy': 0.7, 'Sunny': 0.3},44    'Sunny' : {'Rainy': 0.4, 'Sunny': 0.6},45    }46 emission_probability = {47    'Rainy' : {'walk': 0.1, 'shop': 0.4, 'clean': 0.5},48    'Sunny' : {'walk': 0.6, 'shop': 0.3, 'clean': 0.1},49    }50 #A simple example of the using algorithm51 def example():52     return forward_viterbi(observations,53                            states,54                            start_probability,55                            transition_probability,56                            emission_probability)57 print example()58 59 #(0.033611999999999996, ['Sunny', 'Rainy', 'Rainy', 'Rainy'], 0.009407999999999998)

0 0