1. 维特比算法概述
维特比算法是一个通用的解码算法，是基于动态规划的求序列最短路径的方法。
第一个局部状态是在时刻$t$隐藏状态为$i$所有可能的状态转移路径$i_1,i_2,...i_t$中的概率最大值。记为 $\delta_t(i)$:$$\delta_t(i) = \max_{i_1,i_2,...i_{t-1}}\;P(i_t=i, i_1,i_2,...i_{t-1},o_t,o_{t-1},...o_1|\lambda),\; i =1,2,...N$$
由$\delta_t(i)$的定义可以得到$\delta$的递推表达式：$$\begin{align} \delta_{t+1}(i) & = \max_{i_1,i_2,...i_{t}}\;P(i_{t+1}=i, i_1,i_2,...i_{t},o_{t+1},o_{t},...o_1|\lambda) \\ & = \max_{1 \leq j \leq N}\;[\delta_t(j)a_{ji}]b_i(o_{t+1})\end{align}$$
第二个局部状态由第一个局部状态递推得到。我们定义在时刻$t$隐藏状态为$i$的所有单个状态转移路径$(i_1,i_2,...,i_{t-1},i)$中概率最大的转移路径中第$t-1$个节点的隐藏状态为$\Psi_t(i)$,其递推表达式可以表示为：$$\Psi_t(i) = arg \; \max_{1 \leq j \leq N}\;[\delta_{t-1}(j)a_{ji}]$$
有了这两个局部状态，我们就可以从时刻0一直递推到时刻$T$，然后利用$\Psi_t(i)$记录的前一个最可能的状态节点回溯，直到找到最优的隐藏状态序列。
2. 维特比算法流程总结
现在我们来总结下维特比算法的流程：
输入：HMM模型$\lambda = (A, B, \Pi)$，观测序列$O=(o_1,o_2,...o_T)$
输出：最有可能的隐藏状态序列$I^*= \{i_1^*,i_2^*,...i_T^*\}$
1）初始化局部状态：$$\delta_1(i) = \pi_ib_i(o_1),\;i=1,2...N$$$$\Psi_1(i)=0,\;i=1,2...N$$
2) 进行动态规划递推时刻$t=2,3,...T$时刻的局部状态：$$\delta_{t}(i) = \max_{1 \leq j \leq N}\;[\delta_{t-1}(j)a_{ji}]b_i(0_{t}),\;i=1,2...N$$$$\Psi_t(i) = arg \; \max_{1 \leq j \leq N}\;[\delta_{t-1}(j)a_{ji}],\;i=1,2...N$$
3) 计算时刻$T$最大的$\delta_{T}(i)$,即为最可能隐藏状态序列出现的概率。计算时刻$T$最大的$\Psi_t(i)$,即为时刻$T$最可能的隐藏状态。$$P* = \max_{1 \leq j \leq N}\delta_{T}(i)$$$$i_T^* = arg \; \max_{1 \leq j \leq N}\;[\delta_{T}(i)]$$
4) 利用局部状态$\Psi(i)$开始回溯。对于$t=T-1,T-2,...,1$：$$i_t^* = \Psi_{t+1}(i_{t+1}^*)$$
最终得到最有可能的隐藏状态序列$I^*= \{i_1^*,i_2^*,...i_T^*\}$
3. HMM维特比算法求解实例
下面我们仍然用隐马尔科夫模型HMM中盒子与球的例子来看看HMM维特比算法求解。
我们的观察集合是:$$V=\{红，白\}，M=2$$
我们的状态集合是：$$Q =\{盒子1，盒子2，盒子3\}， N=3 $$
而观察序列和状态序列的长度为3.
初始状态分布为：$$\Pi = (0.2,0.4,0.4)^T$$
状态转移概率分布矩阵为：
$$A = \left( \begin{array} {ccc} 0.5 & 0.2 & 0.3 \\ 0.3 & 0.5 & 0.2 \\ 0.2 & 0.3 &0.5 \end{array} \right) $$
观测状态概率矩阵为：
$$B = \left( \begin{array} {ccc} 0.5 & 0.5 \\ 0.4 & 0.6 \\ 0.7 & 0.3 \end{array} \right) $$
球的颜色的观测序列:$$O=\{红，白，红\}$$
首先需要得到三个隐藏状态在时刻1时对应的各自两个局部状态，此时观测状态为1：
$$\delta_1(1) = \pi_1b_1(o_1) = 0.2 \times 0.5 = 0.1$$
$$\delta_1(2) = \pi_2b_2(o_1) = 0.4 \times 0.4 = 0.16$$
$$\delta_1(3) = \pi_3b_3(o_1) = 0.4 \times 0.7 = 0.28$$
$$\Psi_1(1)=\Psi_1(2) =\Psi_1(3) =0$$
现在开始递推三个隐藏状态在时刻2时对应的各自两个局部状态，此时观测状态为2：
$$\delta_2(1) = \max_{1\leq j \leq 3}[\delta_1(j)a_{j1}]b_1(o_2) = \max_{1\leq j \leq 3}[0.1 \times 0.5, 0.16 \times 0.3, 0.28\times 0.2] \times 0.5 = 0.028$$
$$\Psi_2(1)=3$$
$$\delta_2(2) = \max_{1\leq j \leq 3}[\delta_1(j)a_{j2}]b_2(o_2) = \max_{1\leq j \leq 3}[0.1 \times 0.2, 0.16 \times 0.5, 0.28\times 0.3] \times 0.6 = 0.0504$$
$$\Psi_2(2)=3$$
$$\delta_2(3) = \max_{1\leq j \leq 3}[\delta_1(j)a_{j3}]b_3(o_2) = \max_{1\leq j \leq 3}[0.1 \times 0.3, 0.16 \times 0.2, 0.28\times 0.5] \times 0.3 = 0.042$$
$$\Psi_2(3)=3$$
继续递推三个隐藏状态在时刻3时对应的各自两个局部状态，此时观测状态为1：
$$\delta_3(1) = \max_{1\leq j \leq 3}[\delta_2(j)a_{j1}]b_1(o_3) = \max_{1\leq j \leq 3}[0.028 \times 0.5, 0.0504 \times 0.3, 0.042\times 0.2] \times 0.5 = 0.00756$$
依最大概率计算结果： $$\Psi_3(1)=2$$
$$\delta_3(2) = \max_{1\leq j \leq 3}[\delta_2(j)a_{j2}]b_2(o_3) = \max_{1\leq j \leq 3}[0.028 \times 0.2, 0.0504\times 0.5, 0.042\times 0.3] \times 0.4 = 0.01008$$
依最大概率计算结果：$$\Psi_3(2)=2$$
$$\delta_3(3) = \max_{1\leq j \leq 3}[\delta_2(j)a_{j3}]b_3(o_3) = \max_{1\leq j \leq 3}[0.028 \times 0.3, 0.0504 \times 0.2, 0.042\times 0.5] \times 0.7 = 0.0147$$
依最大概率计算结果：$$\Psi_3(3)=3$$
此时已经到最后的时刻，我们开始准备回溯。此时最大概率为$\delta_3(3)$,从而得到$i_3^* =3$
由于$\Psi_3(3)=3$,所以$i_2^* =3$, 而又由于$\Psi_2(3)=3$,所以$i_1^* =3$。从而得到最终的最可能的隐藏状态序列为：$(3,3,3)$
代码实例
1 import numpy as np 2 def viterbi(trainsition_probability,emission_probability,pi,obs_seq): 3 #转换为矩阵进行运算 4 trainsition_probability=np.array(trainsition_probability) 5 emission_probability=np.array(emission_probability) 6 pi=np.array(pi) 7 obs_seq = [0, 2, 3] 8 # 最后返回一个Row*Col的矩阵结果 9 Row = np.array(trainsition_probability).shape[0] 10 Col = len(obs_seq) 11 #定义要返回的矩阵 12 F=np.zeros((Row,Col)) 13 #初始状态 14 F[:,0]=pi*np.transpose(emission_probability[:,obs_seq[0]]) 15 for t in range(1,Col): 16 list_max=[] 17 for n in range(Row): 18 list_x=list(np.array(F[:,t-1])*np.transpose(trainsition_probability[:,n])) 19 #获取最大概率 20 list_p=[] 21 for i in list_x: 22 list_p.append(i*10000) 23 list_max.append(max(list_p)/10000) 24 F[:,t]=np.array(list_max)*np.transpose(emission_probability[:,obs_seq[t]]) 25 return F 26 27 if __name__=='__main__': 28 #隐藏状态 29 invisible=['Sunny','Cloud','Rainy'] 30 #初始状态 31 pi=[0.63,0.17,0.20] 32 #转移矩阵 33 trainsion_probility=[[0.5,0.375,0.125],[0.25,0.125,0.625],[0.25,0.375,0.375]] 34 #发射矩阵 35 emission_probility=[[0.6,0.2,0.15,0.05],[0.25,0.25,0.25,0.25],[0.05,0.10,0.35,0.5]] 36 #最后显示状态 37 obs_seq=[0,2,3] 38 #最后返回一个Row*Col的矩阵结果 39 F=viterbi(trainsion_probility,emission_probility,pi,obs_seq) 40 print(F)
结果：
[[ 0.378 0.02835 0.00070875]
[ 0.0425 0.0354375 0.00265781]
[ 0.01 0.0165375 0.01107422]]
每列代表Dry,Damp,Soggy的概率，每行代表Sunny,Cloud,Rainy，所以可以看出最大概率的天气为｛Sunny，Cloud，Rainy｝