隐马尔可夫模型二(公式推导)

• • 概率问题
    • 直接计算法
  • 前向后向算法
    • • 前向算法
      • 后向算法
    • 一些期望
  • 学习问题
    • Baum-Welch算法
    • Baum-Welch参数估计公式
  • 预测算法
    • 近似算法
    • 维比特算法
  • 参考文献

前面一篇介绍了隐马尔科夫模型的基本的一些概念，篇主要介绍三个问题的具体解决方法。如果对于概念不太理解的可以参考前一篇博客HMM模型基本概念，本篇博客主要介绍对于三个问题的主要推倒，内容主要基于统计学习方法这本书，但是在上面加上了一些自己的理解。下面一一介绍三个问题以及解决的办法。

概率问题

给定模型λ=(A,B,π)$\lambda =(A,B,\pi )$和观测序列O=(o1,o2,...,oT)$O=(o_1,o_2,...,o_T)$，计算在模型λ$\lambda$下观测序列O$O$出现的概率P(O|λ)$P(O|\lambda )$。

直接计算法

直接计算法说白了就是暴力计算每一种情况的可能。对于所有可能的状态序列I$I$求和，得到观测序列O$O$的概率P(O|λ)$P(O|\lambda )$，即:

P(O|λ)=∑IP(O|I,λ)P(I|λ)$$P(O|\lambda )=\sum _{I}P(O|I,\lambda )P(I|\lambda )$$

=∑i1,i2,...,iTπi1bi1(o1)ai1i2bi2(o2)...aiT−1iTbiT(oT)$$=\sum _{i_1,i_2,...,i_T}{\pi }_{i_1}{b}_{i_1}(o_1)a_{i_1{i}_2}{b}_{i_2}(o_2)...a_{i_{T-1}{i}_T}{b}_{i_T}(o_T)$$

这种计算的缺点在于计算量很大，时间复杂度为O(TNT)$O(T{N}^T)$。

前向后向算法

前向后向算法的核心是利用动态规划的思想减少计算的时间复杂度。

图1

前向算法

前向概率 给定隐马尔可夫模型λ$\lambda$，定义到时刻t$t$部分观测序列为o1,o2,...,ot$o_1,o_2,...,o_t$且状态为qi$q_i$的概率为前向概率，记作

αt(i)=P(o1,o2,...,ot,it=qi|λ)(1)$${\alpha }_t(i)=P(o_1,o_2,...,o_t,i_t=q_i|\lambda )(1)$$

然后可以递推求出前向概率αt(i)${\alpha }_t(i)$以及观测序列P(O|λ)$P(O|\lambda )$

盒子	1	2	3	4
红球数	5	3	6	8
白球数	5	7	4	2

那么这个前向概率到底是什么意思呢？估计好多人还是看的一脸懵逼。还是以之前的盒子与球模型,观测序列为O={红,红,白,白,红}$O=\{红,红,白,白,红\}$，假设t=3,i=1$t=3,i=1$,后面的序列我们不知道,那么αt(i)=P(O={红,红,白},i=1|λ)${\alpha }_t(i)=P(O=\{红,红,白\},i=1|\lambda )$。即前面观测序列为[红,红,白],第三次丑的白色球是从盒子1中抽出的概率。

下面是对前向算法的形式化推导。

输入：隐马尔科夫模型λ$\lambda$,观测序列为O$O$;

输出：观测序列概率P(O|λ)$P(O|\lambda )$;

1.初值

α1(i)=πibi(o1)(2)$${\alpha }_1(i)={\pi }_i{b}_i(o_1)(2)$$

即求第一个观测值对应的状态为i$i$的概率

2.递推，对t=1,2,…,T-1,

αt+1(i)=[∑j=1Nαt(j)aji]bi(ot+1),i=1,2,....,N(3)$${\alpha }_{t+1}(i)=[\sum _{j=1}^N{\alpha }_t(j)a_{ji}]b_i(o_{t+1}),i=1,2,....,N(3)$$

3.终止

P(O|λ)=∑i=1NαT(i)(4)$$P(O|\lambda )=\sum _{i=1}^N{\alpha }_T(i)(4)$$

上面这段公式推导的思想为，先求观测值为o1$o_1$的概率，然后在前面的基础上求观测值为o2$o_2$对应的概率,依次递推，最后求观测值为ot$o_t$的概率。然后每一种概率对应了不同的状态，αt(i)${\alpha }_t(i)$则表示在状态为i$i$时的概率，而P(o1,o2,...,ot,iT=qi|λ)$P(o_1,o_2,...,o_t,i_T=q_i|\lambda )$，然后然后对i$i$求和就能够得到上述(3)式。其递推过程图如下：

图2

前向算法就是利用这种方式将时间复杂度从O(TNT)$O(T{N}^T)$降低到O(TN2)$O(T{N}^2)$,至于降低的原因是减少直接利用了前面的计算结果，避免了每一次都需要重新计算。

还是以盒子与球模型为例，λ=(A,B,π)$\lambda =(A,B,\pi )$,状态集合Q={1,2,3}$Q=\{1,2,3\}$,观测集合V={红,白}$V=\{红,白\}$.

设T=3$T=3$,O={红,白,红}$O=\{红,白,红\}$，求P(O|λ)$P(O|\lambda )$

1.计算初值

α1(1)=π1b1(o1)=0.2×0.5=0.10$${\alpha }_1(1)={\pi }_1{b}_1(o_1)=0.2\times 0.5=0.10$$

α1(2)=π2b2(o1)=0.4×0.4=0.16$${\alpha }_1(2)={\pi }_2{b}_2(o_1)=0.4\times 0.4=0.16$$

α1(3)=π3b3(o1)=0.4×0.7=0.28$${\alpha }_1(3)={\pi }_3{b}_3(o_1)=0.4\times 0.7=0.28$$

2.递推计算

α2(1)=[∑j=13α1(j)aj1]b1(o2)$${\alpha }_2(1)=[\sum _{j=1}^3{\alpha }_1(j)a_{j_1}]b_1(o_2)$$

=(0.10×0.5+0.16×0.3+0.28×0.2)=0.154×0.5=0.77$$=(0.10\times 0.5+0.16\times 0.3+0.28\times 0.2)=0.154\times 0.5=0.77$$

α2(2)=[∑j=13α1(j)aj2]b2(o2)$${\alpha }_2(2)=[\sum _{j=1}^3{\alpha }_1(j)a_{j2}]b_2(o_2)$$

α2(3)=[∑j=13α1(j)aj3]b3(o2)$${\alpha }_2(3)=[\sum _{j=1}^3{\alpha }_1(j)a_{j3}]b_3(o_2)$$

同理可以得到

α3(1)=[∑j=13α2(j)aj1]b1(o3)=0.04187$${\alpha }_3(1)=[\sum _{j=1}^3{\alpha }_2(j)a_{j1}]b_1(o_3)=0.04187$$

α3(2)=[∑j=13α2(j)aj2]b2(o3)=0.03551$${\alpha }_3(2)=[\sum _{j=1}^3{\alpha }_2(j)a_{j2}]b_2(o_3)=0.03551$$

α3(3)=[∑j=13α2(j)aj3]b3(o3)=0.05284$${\alpha }_3(3)=[\sum _{j=1}^3{\alpha }_2(j)a_{j3}]b_3(o_3)=0.05284$$

3.终止

P(O|λ)=∑j=13α3(j)=0.13022$$P(O|\lambda )=\sum _{j=1}^3{\alpha }_3(j)=0.13022$$

前向算法清楚了，其实后向概率也就清楚了。其本质就是和前向算法的思想是一样的，只不过方向相反，从后往前计算。以下图为例,先计算的是qt$q_t$的概率，然后在计算qt−1$q_{t-1}$的概率，以此类推。下面直接给出定义以及推导。

后向算法

后向概率$$给定隐马尔科夫模型λ$\lambda$,定义在t$t$时刻状态为qt$q_t$，从t+1$t+1$到T$T$的部分观测序列为ot+1,ot+2,...,oT$o_{t+1},o_{t+2},...,o_T$的概率为后向概率，记作

βt(i)=P(ot+1,ot+2,...,oT|it=qi,λ)$${\beta }_t(i)=P(o_{t+1},o_{t+2},...,o_T|i_t=q_i,\lambda )$$

然后可以利用递推从后向前求解得到P(O|λ)$P(O|\lambda )$，计算过程如下：

1.首先，对于最终时刻的所有状态qi$q_i$规定βT(i)=1${\beta }_T(i)=1$,即：

βT(i)=1,1,2,...,N(5)$${\beta }_T(i)=1,1,2,...,N(5)$$

本来应该对于这一步应该像前向算法一样存在初始状态概率，但是后向算法将初始概率放到最后一步计算，所有令βT(i)${\beta }_T(i)$规定为=1

2.对于t=T−1,T−2,...,1$t=T-1,T-2,...,1$

βt(i)=[∑j=1Nβt+1(j)aji]bi(ot)(6)$${\beta }_t(i)=[\sum _{j=1}^N{\beta }_{t+1}(j)a_{j_i}]b_i(o_t)(6)$$

这一步和书上不一样的，个人感觉这一步书上的写法是错误的。

3.终止

P(O|λ)=∑i=1Nπibi(o1)β1(i)(7)$$P(O|\lambda )=\sum _{i=1}^N{\pi }_i{b}_i(o_1){\beta }_1(i)(7)$$

并且可以将前向后向概率统一：

P(O|λ)=∑i=1N∑j=1Nαt(i)aijbj(ot+1)βt+1(j)(8)$$P(O|\lambda )=\sum _{i=1}^N\sum _{j=1}^N{\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)(8)$$

最后化简能够得到

P(O|λ)=∑i=1N∑j=1Nαt(i)aijbj(ot+1)βt+1(j)(9)$$P(O|\lambda )=\sum _{i=1}^N\sum _{j=1}^N{\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)(9)$$

使用前向后向算法是为了更好的描述后面两个问题，其实两个算法本质就是一回事o(╯□╰)o

一些期望

这一部分主要是简化一些符号，为后面的计算做准备。

1.给定模型λ$\lambda$和观测O$O$，在时刻t$t$处于状态qi$q_i$的概率，记为：

γi=P(it=qi|O,λ)$${\gamma }_i=P(i_t=q_i|O,\lambda )$$

可以通过前向后向概率计算。

γi(i)=P(it=qi|O,λ)P(O|λ)$${\gamma }_i(i)=\frac{P(i_t=q_i|O,\lambda )}{P(O|\lambda )}$$

而

P(it=qi|O,λ)=αt(i)βt(i)$$P(i_t=q_i|O,\lambda )={\alpha }_t(i){\beta }_t(i)$$

所以有

γt(i)=αt(i)βt(i)∑Ni=1αt(i)βt(i)(9)$${\gamma }_t(i)=\frac{{\alpha }_t(i){\beta }_t(i)}{\sum _{i=1}^N{\alpha }_t(i){\beta }_t(i)}(9)$$

2.给定模型λ$\lambda$和观测O$O$,在时刻t$t$处于状态qi$q_i$且在时刻t+1$t+1$处于状态qj$q_j$的概率，记

ξt(i,j)=P(it=qi,it+1=qj|O,λ)(10)$${\xi }_t(i,j)=P(i_t=q_i,i_{t+1}=q_j|O,\lambda )(10)$$

所以有

ξt(i,j)=αt(i)aijbj(ot+1)βt+1(j)∑Ni=1∑Nj=1αt(i)aijbj(ot+1)βt+1(j)(11)$${\xi }_t(i,j)=\frac{{\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)}{\sum _{i=1}^N\sum _{j=1}^N{\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)}(11)$$

3.一些有用的期望

1)在观测O$O$下状态i$i$出现的期望值

∑t=1Tγt(i)(12)$$\sum _{t=1}^T{\gamma }_t(i)(12)$$

2)在观测O$O$下由状态i$i$转移的期望值

如上图所示，从状态t$t$出发，到其他的状态则称为状态转移，由于最后一个状态是不能转移到下一个状态因此,状态转移的期望值为：

∑t=1T−1γt(i)(13)$$\sum _{t=1}^{T-1}{\gamma }_t(i)(13)$$

3).在观测为O$O$下，由状态i$i$转移到状态j$j$的期望值为：

∑t=1T−1ξt(i,j)(14)$$\sum _{t=1}^{T-1}{\xi }_t(i,j)(14)$$

学习问题

学习问题是为了计算模型参数。在已经给定了观测序列O$O$，根据是否给定状态序列I$I$可以分为监督学习方法和非监督学习方法。监督学习的方法可以利用极大似然，非监督学习主要是利用Baum-Welch算法。

Baum-Welch算法

给定训练集数据为S个长度为T$T$的观测序列O={O1,O2,...,Os}$O=\{O_1,O_2,...,O_s\}$，而没有对应的状态序列，目标是学习隐马尔可夫模型λ=(A,B,π)$\lambda =(A,B,\pi )$,观测序列为O$O$，状态序列为I$I$，那么我们可以将P(O|λ)$P(O|\lambda )$变成包含隐变量的概率模型：

P(O|λ)=∑IP(O|I,λ)P(I|λ)(15)$$P(O|\lambda )=\sum _{I}P(O|I,\lambda )P(I|\lambda )(15)$$

然后可以通过EM算法学习参数，步骤如下：

1.确定完全数据的对数似然函数

观测序列数据为O=(o1,o2,...,oT)$O=(o_1,o_2,...,o_T)$,状态序列数据(隐数据)为I=(i1,i2,...,it)$I=(i_1,i_2,...,i_t)$，完全数据为(O,I)$(O,I)$,所以完全数据的对数似然函数为logP(O,I|λ)$logP(O,I|\lambda )$

2.EM算法的E步，求Q$Q$函数Q(λ,λ¯)$Q(\lambda ,\overline{\lambda })$

Q(λ,λ¯)=∑IP(O,I|λ¯)logP(O,I|λ)(16)$$Q(\lambda ,\overline{\lambda })=\sum _{I}P(O,I|\overline{\lambda })logP(O,I|\lambda )(16)$$

其中

P(O,I|λ)=πi1bi1(o1)ai1i2bi2(o2)....aiT−1iTbiT(oT)$$P(O,I|\lambda )={\pi }_{i_1}{b}_{i_1}(o_1)a_{i_1{i}_2}{b}_{i_2}(o_2)....a_{i_{T-1}{i}_T}{b}_{i_T}(o_T)$$

所以有：

Q(λ,λ¯)=∑IP(O,I|λ¯){logπi1+∑t=1T−1logaitit+1+∑t=1Tlogbit(ot)}(17)$$Q(\lambda ,\overline{\lambda })=\sum _{I}P(O,I|\overline{\lambda })\left\{log{\pi }_{i_1}+\sum _{t=1}^{T-1}log{a}_{i_t{i}_{t+1}}+\sum _{t=1}^{T}log{b}_{i_t}(o_t)\right\}(17)$$

这里不得不吐槽一下李航统计学习方法的符号写法，看了半天才明白到底是怎么写的。

3.EM的M步，极大化Q$Q$函数

可以将上述Q$Q$函数拆分为三项，其中第一项为

∑IP(O,I|λ¯)logπi1$$\sum _{I}P(O,I|\overline{\lambda })log{\pi }_{i_1}$$

并且有∑Ni=1πi=1$\sum _{i=1}^N{\pi }_i=1$,然后利用拉格朗日乘子法，写出拉格朗日函数：

∑i=1NlogπiP(O,i1=i|λ¯)+γ(∑i=1Nπi−1)$$\sum _{i=1}^{N}log{\pi }_{i}P(O,i_1=i|\overline{\lambda })+\gamma \left(\sum _{i=1}^N{\pi }_i-1\right)$$

对齐求偏导并且令结果为0，于是有

∂∂πi[∑i=1NπiP(O,i1=i|λ¯)+γ{∑i=1Nπi−1}]=0(17)$$\frac{\mathrm{\partial }}{\mathrm{\partial }{\pi }_i}\left[\sum _{i=1}^N{\pi }_{i}P(O,i_1=i|\overline{\lambda })+\gamma \left\{\sum _{i=1}^N{\pi }_i-1\right\}\right]=0(17)$$

得

P(O,i1=i|λ¯)+γπi=0(18)$$P(O,i_1=i|\overline{\lambda })+\gamma {\pi }_i=0(18)$$

两边同时对i$i$求和有

γ=−P(O|λ¯)$$\gamma =-P(O|\overline{\lambda })$$

带入到上一步有

πi=P(O,i1=i|λ¯)P(O|λ¯)(19)$${\pi }_i=\frac{P(O,i_1=i|\overline{\lambda })}{P(O|\overline{\lambda })}(19)$$

第二项可以写成

∑IP(O,I|λ¯)(∑t=1T−1logaitit+1)=∑i=1N∑j=1N∑t=1T−1P(O,it=i,it+1=j|λ¯)logaij(20)$$\sum _{I}P(O,I|\overline{\lambda })\left(\sum _{t=1}^{T-1}log{a}_{i_t{i}_{t+1}}\right)=\sum _{i=1}^N\sum _{j=1}^N\sum _{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })log{a}_{ij}(20)$$

和上面的类似，有约束条件∑Nj=1aij=1$\sum _{j=1}^N{a}_{ij}=1$和∑Ni=1aij=1$\sum _{i=1}^N{a}_{ij}=1$的拉格朗日乘子法可求出

∑i=1N∑j=1N∑t=1T−1P(O,it=i,it+1=j|λ¯)logaij+γ(∑j=1Naij−1)$$\sum _{i=1}^N\sum _{j=1}^N\sum _{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })log{a}_{ij}+\gamma (\sum _{j=1}^N{a}_{ij-1})$$

两边同时对aij$a_{ij}$求导得

∑t=1T−1P(O,it=i,it+1=j|λ¯)+γaij=0(21)$$\sum _{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })+\gamma a_{ij}=0(21)$$

两边同时对j$j$求和有

γ=−∑t=1T−1P(O,it=i|λ¯)(22)$$\gamma =-\sum _{t=1}^{T-1}P(O,i_t=i|\overline{\lambda })(22)$$

将式(22)带入(21)有

aij=∑T−1t=1P(O,it=i,it+1=j|λ¯)∑t=1T−1P(O,it=i,|λ¯)(23)$$a_{ij}=\frac{\sum _{t=1}^{T-1}P(O,i_t=i,i_{t+1}=j|\overline{\lambda })}{\sum _{t=1}T-1P(O,i_t=i,|\overline{\lambda })}(23)$$

然后是对bj(k)$b_j(k)$的计算

∑IP(O,I|λ¯)(∑t=1Tlogbit(ot))=∑j=1N∑t=1TP(O,it=j|λ¯)logbj(ot)(24)$$\sum _{I}P(O,I|\overline{\lambda })\left(\sum _{t=1}^{T}log{b}_{i_t}(o_t)\right)=\sum _{j=1}^N\sum _{t=1}^{T}P(O,i_t=j|\overline{\lambda })log{b}_j(o_t)(24)$$

并且约束条件为∑Mk=1bj(k)=1$\sum _{k=1}^M{b}_j(k)=1$

所以构造的拉格朗日函数为

∑j=1N∑t=1TP(O,it=j|λ¯)logbj(ot)+γ(∑k=1Mbj(k)−1)=0$$\sum _{j=1}^N\sum _{t=1}^{T}P(O,i_t=j|\overline{\lambda })log{b}_j(o_t)+\gamma (\sum _{k=1}^M{b}_j(k)-1)=0$$

对bj(k)$b_j(k)$求导有

∂∂bj(k)[∑j=1N∑t=1TP(O,it=j|λ¯)logbj(ot)+γ(∑k=1Mbj(k)−1)]=0$$\frac{\mathrm{\partial }}{\mathrm{\partial }{b}_j(k)}\left[\sum _{j=1}^N\sum _{t=1}^{T}P(O,i_t=j|\overline{\lambda })log{b}_j(o_t)+\gamma (\sum _{k=1}^M{b}_j(k)-1)\right]=0$$

注意，只有当ot=k$o_t=k$时，偏导才不为0所以有

∑t=1TP(O,it=j|λ¯)I(ot=k)+γbj(k)=0(25)$$\sum _{t=1}^{T}P(O,i_t=j|\overline{\lambda })I(o_t=k)+\gamma b_j(k)=0(25)$$

所以同时k求和,并且我们很容易得出∑Mk=1I(ot=k)=1$\sum _{k=1}^{M}I(o_t=k)=1$

γ=−∑t=1TP(O,ot=j|λ¯)(26)$$\gamma =-\sum _{t=1}^{T}P(O,o_t=j|\overline{\lambda })(26)$$

所以将式子(26)带入(25)则有

bj(k)=∑Tt=1P(O,it=j|λ¯)I(ot=k)∑Tt=1P(O,it=j|λ¯)$$b_j(k)=\frac{\sum _{t=1}^{T}P(O,i_t=j|\overline{\lambda })I(o_t=k)}{\sum _{t=1}^{T}P(O,i_t=j|\overline{\lambda })}$$

Baum-Welch参数估计公式

aij=∑T−1t=1ξt(i,j)∑T−1t=1γt(i)$$a_{ij}=\frac{\sum _{t=1}^{T-1}{\xi }_t(i,j)}{\sum _{t=1}^{T-1}{\gamma }_t(i)}$$

bj(k)=∑Tt=1,ot=kγt(j)∑Tt=1γt(j)$$b_j(k)=\frac{\sum _{t=1,o_t=k}^T{\gamma }_t(j)}{\sum _{t=1}^T{\gamma }_t(j)}$$

πi=γ1(i)$${\pi }_i={\gamma }_1(i)$$

预测算法

回到上一篇博客的内容，预测问题是个什么问题呢？预测问题也叫做解码问题，即给定隐马尔科夫模型参数λ=(A,B,π)$\lambda =(A,B,\pi )$,以及观测序列O=(o1,o2,o3,...,ot)$O=(o_1,o_2,o_3,...,o_t)$,求P(I|O)$P(I|O)$,即最有可能出现的状态

解决上面的问题主要有两种方法，一种是近似算法，另外一种是维比特算法。

近似算法

近似算法的思想其实很简单，在每一个t$t$时刻选择最可能出现的状态i∗t$i_t^{\ast }$，从而得到一个近似状态I∗=(i∗1,i∗2,...,i∗T)$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$。

在t时刻处于qi$q_i$的状态的概率为γt(i)${\gamma }_t(i)$

所以在每一个t时刻最有可能的状态为

i∗t=argmax[γt(i)],t=1,2,3..T,1≤t≤N$$i_t^{\ast }=argmax[{\gamma }_t(i)],t=1,2,\mathrm{3..}T,1\le t\le N$$

从而得到状态序列I∗=(i∗1,i∗2,...,i∗T)$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$.但是这种方法只保证了每一个t时刻最有可能的状态，不能保证整体，有点贪心的思想在里面

维比特算法

维比特算法其实就是利用动态规划的思想来求最大路径，类似于以前学的利用动态规划解决最短路径问题。文字描述什么的感觉不容易理解，直接上例子感觉跟容易理解。但是在说明例子之前，我们先定义两个符号δ$\delta$和ψ$\psi$方便后面的计算。

定义在时刻t$t$状态为i$i$的所有单个路径(i1,i2,...,iT)$(i_1,i_2,...,i_T)$中概率最大的值为

δt(i)=maxi1,i2,...,iTP(it=i|λ),i=1,2,...,N,(28)$${\delta }_t(i)=\underset{i_1,i_2,...,i_T}{max}P(i_t=i|\lambda ),i=1,2,...,N,(28)$$

定义在时刻t$t$状态为i$i$的所有单个儿路径(i1,i2,...,it−1,i)$(i_1,i_2,...,i_{t-1},i)$中最大概率的路径的第t−1$t-1$个节点为

ψt(i)=argmax1≤j≤Nψt1(j)aji$${\psi }_t(i)=arg\underset{1\le j\le N}{max}{\psi }_{t_1}(j)a_{ji}$$

下面直接符号化维比特算法具体过程

输入：模型λ=(A,B,π)$\lambda =(A,B,\pi )$和观测O=(o1,o2,...,oT)$O=(o_1,o_2,...,o_T)$

输出：最优路径I∗=(i∗1,i∗2,...,i∗T)$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$

1)初始化

δ1(i)=πibi(o1),i=1,2,...,N$${\delta }_1(i)={\pi }_i{b}_i(o_1),i=1,2,...,N$$

ψ1(i)=0$${\psi }_1(i)=0$$

2)递推：对于t=1,2,3...,T$t=1,2,\mathrm{3...},T$

δ1(i)=max1≤j≤N[δt−1(i)aji]bi(ot),i=1,2,...,N$${\delta }_1(i)=\underset{1\le j\le N}{max}[{\delta }_{t-1}(i)a_{ji}]b_i(o_t),i=1,2,...,N$$

ψi(t)=argmax1≤i≤N[δt−1aji],i=1,2,...,N$${\psi }_i(t)=arg\underset{1\le i\le N}{max}[{\delta }_{t-1}{a}_{ji}],i=1,2,...,N$$

3).终止

P∗=max1≤i≤NδT(i)$$P^{\ast }=\underset{1\le i\le N}{max}{\delta }_T(i)$$

i∗t=argmax1≤i≤N[δT(i)]$$i_t^{\ast }=arg\underset{1\le i\le N}{max}[{\delta }_T(i)]$$

4).最优路径回溯。对于t=T−1,T−2,...,1$t=T-1,T-2,...,1$

I∗=(i∗1,i∗2,i∗3...,i∗T)$$I^{\ast }=(i_1^{\ast },i_2^{\ast },i_3^{\ast }...,i_T^{\ast })$$

那么具体过程是什么样的呢？

还是以上盒子与球的模型为例

已知观测序列为O=(红,白,红)$O=(红,白,红)$，试求最优状态序列I∗=(i∗1,i∗2,...,i∗T)$I^{\ast }=(i_1^{\ast },i_2^{\ast },...,i_T^{\ast })$

(1)初始化，代入公式有

δ1(1)=0.10,δ1(2)=0.16,δ1(3)=0.28$${\delta }_1(1)=0.10,{\delta }_1(2)=0.16,{\delta }_1(3)=0.28$$

记ψ1(i)=0,i=1,2,3${\psi }_1(i)=0,i=1,2,3$

最优路径求求截图如上

(2)t=2的时候如何计算呢，这里以t=2,i=1$t=2,i=1$为例

δ2(1)=max1≤j≤3[δ1(j)aj1b1(o2)]$${\delta }_2(1)=\underset{1\le j\le 3}{max}[{\delta }_1(j)a_{j1}{b}_1(o_2)]$$

=maxj{0.1×0.5,0.16×0.3,0.28×0.2}×0.5$$=\underset{j}{max}\{0.1\times 0.5,0.16\times 0.3,0.28\times 0.2\}\times 0.5$$

=0.028$$=0.028$$

根据计算得到当前路径从3到1的概率最大，所以有ψ2(1)=3${\psi }_2(1)=3$根据上面的公式计算得到：

δ2(2)=0.0504,ψ2(2)=3$${\delta }_2(2)=0.0504,{\psi }_2(2)=3$$

δ2(3)=0.042,ψ2(3)=3$${\delta }_2(3)=0.042,{\psi }_2(3)=3$$

δ3(1)=0.0.0756,ψ3(1)=2$${\delta }_3(1)=\mathrm{0.0.0756},{\psi }_3(1)=2$$

δ3(2)=0.01008,ψ3(2)=2$${\delta }_3(2)=0.01008,{\psi }_3(2)=2$$

δ3(3)=0.0147,ψ3(3)=3$${\delta }_3(3)=0.0147,{\psi }_3(3)=3$$

所以最优路径概率为P∗=0.0147$P^{\ast }=0.0147$

(3)倒推最优路径

最优路径最有一个状态对应的是3，而ψ3(3)=3${\psi }_3(3)=3$，所以第二个状态为3，而ψ2(3)=3${\psi }_2(3)=3$，因此第一个状态也为3

所以最优路径为

I∗=(i∗1,i∗2,i∗3)=(3,3,3)$$I^{\ast }=(i_1^{\ast },i_2^{\ast },i_3^{\ast })=(3,3,3)$$

HMM就这样水完了，感觉这篇博客写得像个草稿，有时间再将思想凝练下吧，暂时就这样窘o(╯□╰)o

参考文献

1.隐马尔可夫(HMM)、前/后向算法、Viterbi算法 再次总结

2.统计学习方法，李航