引言
viterbi算法简化最有可能的天气序列的运算过程,forward算法简化该该观察值的概率。
问题描述
你在中国,你朋友F在美国,F的作息有walk, shop, clean,但这选择跟天气有关,我们又知道Rainy的概率比Sunny的概率大
这是初始概率这是天气转移矩阵
这是在相应天气下发生相应事的概率分布
然后,F这三天是walk,walk,shop,问最有可能的天气序列问题分析
同样的,我们先用穷举法来算,即
Sunny Sunny Sunny
Sunny Sunny RainySunny Rainy SunnySunny Rainy RainyRainy Sunny SunnyRainy Sunny RainyRainy Rainy SunnyRainy 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)