概率计算算法

直接计算法

给定模型$\lambda =(A,B,\pi )$和观测序列$O=(o_1,o_2,\cdots ,o_T)$，计算观测序列$O$出现的概率$P(O|\lambda )$。最直接的方法就是按概率公式直接计算。通过列举所有可能的长度为$T$的状态序列$I=(i_1,i_2,\cdots ,i_T)$，求各个状态序列$I$与观测序列$O=(o_1,o_2,\cdots ,o_T)$的联合概率$P(O,I|\lambda )$，然后对所有可能的状态序列求和，就得到$P(O|\lambda )$。

状态序列$I=(i_1,i_2,\cdots ,i_T)$的概率为
$$\begin{array}{rlr}P(I|\lambda )& =P(i_1,i_2,\cdots ,i_T|\lambda )& \\ & =P(i_1)\prod _{t=2}^{T}P(i_t|i_1,\cdots ,i_{t-1}，\lambda )& 联合概率分布\\ & =P(i_1)\prod _{t=2}^{T}P(i_t|i_{t-1},\lambda )& 齐次马尔可夫假设\\ & ={\pi }_{i_1}\prod _{t=2}^T{a}_{i_{t-1},i_t}& \\ & ={\pi }_{i_1}{a}_{i_1,i_2}{a}_{i_2,i_3}\cdots a_{i_{T-1},i_T}& \end{array}\text{\hspace{1em}(10.10)}$$

就是由初始概率生成第一个状态$i_1$，然后转移到第二个状态$i_2$，最后到第$T$个状态$i_T$的概率之积。

对于固定的状态序列$I=(i_1,i_2,\cdots ,i_T)$，观测序列$O=(o_1,o_2,\cdots ,o_T)$的概率是：
$$\begin{array}{rlr}P(O|I,\lambda )& =P(o_1,o_2,\cdots ,o_T|i_1,i_2,\cdots ,i_T,\lambda )& \\ & =P(o_T|o_1,o_2,\cdots ,o_{T-1},i_1,i_2,\cdots ,i_T,\lambda )P(o_1,o_2,\cdots ,o_{T-1}|i_1,i_2,\cdots ,i_T,\lambda )& \\ & =P(o_T|i_T,\lambda )P(o_1,o_2,\cdots ,o_{T-1}|i_1,i_2,\cdots ,i_T,\lambda )& 观测独立假设\\ & =P(o_T|i_T,\lambda )P(o_{T-1}|o_1,o_2,\cdots ,o_{T-2},i_1,i_2,\cdots ,i_T,\lambda )P(o_1,o_2,\cdots ,o_{T-2}|i_1,i_2,\cdots ,i_T,\lambda )& \\ & =P(o_T|i_T,\lambda )P(o_{T-1}|i_{T-1},\lambda )P(o_1,o_2,\cdots ,o_{T-2}|i_1,i_2,\cdots ,i_T,\lambda )& \\ & =\prod _{t=1}^{T}P(o_t|i_t)& \\ & =\prod _{t=1}^T{b}_{i_t}(o_t)& \\ & =b_{i_1}(o_1)b_{i_2}(o_2)\cdots b_{i_T}(o_T)& \end{array}\text{\hspace{1em}(10.11)}$$

由状态$i_1$产生观测$o_1$，状态$i_2$产生观测$o_2$，到状态$i_T$产生观测$o_T$的概率之积。

那么$O$和$I$同时出现的概率为：
$$\begin{array}{rl}P(O,I|\lambda )& =P(O|I,\lambda )P(I|\lambda )\\ & =b_{i_1}(o_1)b_{i_2}(o_2)\cdots b_{i_T}(o_T){\pi }_{i_1}{a}_{i_1,i_2}{a}_{i_2,i_3}\cdots a_{i_{T-1},i_T}\\ & ={\pi }_{i_1}{b}_{i_1}(o_1)a_{i_1,i_2}{b}_{i_2}(o_2)\cdots a_{i_{T-1},i_T}{b}_{i_T}(o_T)\end{array}\text{\hspace{1em}(10.12)}$$
上式只是针对一种状态序列$I$，对所有可能的状态序列$I$求和，就可以得到观测序列$O$的概率$P(O|\lambda )$，即
$$\begin{array}{rl}P(O|\lambda )& =\sum _{I}P(O,I|\lambda )\\ & =\sum _{I}P(O|I,\lambda )P(I|\lambda )\\ & =\sum _{i_1,i_2,\cdots ,i_T}{\pi }_{i_1}{b}_{i_1}(o_1)a_{i_1,i_2}{b}_{i_2}(o_2)\cdots a_{i_{T-1},i_T}{b}_{i_T}(o_T)\end{array}\text{\hspace{1em}(10.13)}$$
但是，直接利用上式计算量非常大。其中${\sum }_{i_1,i_2,\cdots ,i_T}$共有$N^T$种可能(长度为$T$，每个位置都有$N$种可能)，而计算${\pi }_{i_1}{b}_{i_1}(o_1)a_{i_1,i_2}{b}_{i_2}(o_2)\cdots a_{i_{T-1},i_T}{b}_{i_T}(o_T)$的时间复杂度为$O(T)$，所以整体的时间复杂度是$O(T{N}^T)$阶的，这种算法实际上不可行。

下面介绍计算观测序列概率$P(O|\lambda )$的有效算法：前向-后向算法(forward-backward algorithm)。

前向算法

定义10.2(前向概率) 给定隐马尔可夫模型$\lambda$，定义到时刻$t$部分观测序列为$o_1,o_2,\cdots ,o_t$且状态为$q_i$的概率为前向概率，记作
$${\alpha }_t(i)=P(o_1,o_2,\cdots ,o_t,i_t=q_i|\lambda )\text{\hspace{1em}(10.14)}$$
这个公式中$\alpha$的下标变成了观测序列，参数变成了状态。

从这个定义可以看到，显然我们可以基于时刻$a_{t-1}$来表示${\alpha }_t$，即可以递推地求得前向概率${\alpha }_t(i)$及观测序列概率$P(O|\lambda )$。

算法10.2(观测序列概率的前向算法)

输入：隐马尔可夫模型$\lambda$，观测序列$O$；

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

(1) 初值
$${\alpha }_1(i)={\pi }_i{b}_i(o_1),i=1,2,\cdots ,N\text{\hspace{1em}(10.15)}$$
在$t=1$时，只看到观测$o_1$，是初始时刻的状态$i_1=q_i$和观测$o_1$的联合概率，即初始状态为$i$的概率${\pi }_i$乘以 状态$i$产生观测$o_1$的概率$b_i(o_1)$。

(2) 递推 对$t=1,2,\cdots ,T-1$，

计算到时刻$t+1$部分观测序列为$o_1,o_2,\cdots ,o_t,o_{t+1}$且在时刻$t+1$处处于状态$q_i$的前向概率。

公式推导如下：
$$\begin{array}{rlr}{\alpha }_{t+1}(i)& =P(o_1,o_2,\cdots ,o_t,o_{t+1},i_{t+1}=q_i|\lambda )& \\ & =\sum _{j=1}^{N}P(o_1,o_2,\cdots ,o_t,o_{t+1},i_t=q_j,i_{t+1}=q_i|\lambda )& 引入i_t=q_j\\ & =\sum _{j=1}^{N}P(o_{t+1}|o_1,o_2,\cdots ,o_t,i_t=q_j,i_{t+1}=q_i,\lambda )P(o_1,o_2,\cdots ,o_t,i_t=q_j,i_{t+1}=q_i|\lambda )& \\ & =\sum _{j=1}^{N}P(o_{t+1}|i_{t+1}=q_i,\lambda )P(o_1,o_2,\cdots ,o_t,i_t=q_j,i_{t+1}=q_i|\lambda )& 观测独立假设\\ & =\sum _{j=1}^{N}P(o_{t+1}|i_{t+1}=q_i,\lambda )P(i_{t+1}=q_i|o_1,o_2,\cdots ,o_t,i_t=q_j,i_{t+1}=q_i,\lambda )P(o_1,o_2,\cdots ,o_t,i_t=q_j|\lambda )& \\ & =\sum _{j=1}^{N}P(o_{t+1}|i_{t+1}=q_i,\lambda )P(i_{t+1}=q_i|i_t=q_j,\lambda )P(o_1,o_2,\cdots ,o_t,i_t=q_j|\lambda )& 齐次马尔可夫假设\\ & =\left[\sum _{j=1}^{N}P(i_{t+1}=q_i|i_t=q_j,\lambda )P(o_1,o_2,\cdots ,o_t,i_t=q_j|\lambda )\right]P(o_{t+1}|i_{t+1}=q_i)& \\ & =\left[\sum _{j=1}^N{a}_{ji}{\alpha }_t(j)\right]b_i(o_{t+1})& \end{array}$$
把上式整理一下，就得到了书中的形式：

$${\alpha }_{t+1}(i)=\left[\sum _{j=1}^N{\alpha }_t(j)a_{ji}\right]b_i(o_{t+1}),i=1,2,\cdots ,N\text{\hspace{1em}(10.16)}$$

(3) 终止

因为
$${\alpha }_T(i)=P(o_1,o_2,\cdots ,o_T,i_T=q_i|\lambda )$$
所以
$$P(O|\lambda )=\sum _{i=1}^{N}P(O,i_t=q_i|\lambda )=\sum _{i=1}^N{\alpha }_T(i)\text{\hspace{1em}(10.17)}$$
如图10.2所示，前向算法实际是基于“状态序列的路径结构”递推计算$P(O|\lambda )$的算法。其高效的关键在于记住并利用前一步计算的前向概率，避免了重复计算，然后利用路径结构将前向概率递推到全局，得到$P(O|\lambda )$。这样，利用前向概率计算$P(O|\lambda )$的计算量是$O(N^{2}T)$阶的，而不是直接计算的$O(T{N}^T)$。

在每个时间步$t$，每个状态都需要与前一个时间步$t-1$的$N$个状态的结果相乘。而每个时间步$t$，都有$N$个状态，所以是$N^2$，总共有$T$个时间步，所以总量是$O(N^{2}T)$。

后向算法

定义 10.3(后向概率) 给定隐马尔可夫模型$\lambda$，定义在时刻$t$状态为$q_i$的条件下，从$t+1$到$T$的部分观测序列为$o_{t+1},o_{t+2},\cdots ,o_T$的概率为后向概率，记作
$${\beta }_t(i)=P(o_{t+1},o_{t+2},\cdots ,o_T|i_t=q_i,\lambda )\text{\hspace{1em}(10.18)}$$
这里和前向概率有点不同，隐藏状态$i_t=q_i$作为条件观测$o_{t+1:T}$，下面也推导一下。

同样可以用递推的方法求得后向概率${\beta }_t(i)$及观测序列$P(O|\lambda )$。

算法 10.3(观测序列概率的后向算法)

输入： 隐马尔可夫模型$\lambda$，观测序列$O$；

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

(1)
$${\beta }_T(i)=1,i=1,2,\cdots ,N\text{\hspace{1em}(10.19)}$$
初始概率为$1$。

(2) 对$t=T-1,T-2,\cdots 1$

我们希望通过${\beta }_{t+1}(j)=P(o_{t+2},\cdots ,o_T|i_{t+1}=q_j,\lambda )$来表示${\beta }_t(i)$。
$$\begin{array}{rlr}{\beta }_t(i)& =P(o_{t+1},o_{t+2},\cdots ,o_T|i_t=q_i,\lambda )& \\ & =\sum _{j=1}^{N}P(o_{t+1},o_{t+2},\cdots ,o_T,i_{t+1}=q_j|i_t=q_i,\lambda )& \\ & =\sum _{j=1}^{N}P(o_{t+1},o_{t+2},\cdots ,o_T,|i_{t+1}=q_j,i_t=q_i,\lambda )P(i_{t+1}=q_j|i_t=q_i,\lambda )& \\ & =\sum _{j=1}^{N}P(o_{t+1},o_{t+2},\cdots ,o_T,|i_{t+1}=q_j,\lambda )P(i_{t+1}=q_j|i_t=q_i,\lambda )& 阻隔\\ & =\sum _{j=1}^{N}P(o_{t+1}|o_{t+2},\cdots ,o_T,i_{t+1}=q_j,\lambda )P(o_{t+2},\cdots ,o_T|i_{t+1}=q_j,\lambda )P(i_{t+1}=q_j|i_t=q_i,\lambda )& \\ & =\sum _{j=1}^{N}P(o_{t+1}|i_{t+1}=q_j,\lambda )P(o_{t+2},\cdots ,o_T|i_{t+1}=q_j,\lambda )P(i_{t+1}=q_j|i_t=q_i,\lambda )& 观测独立假设\\ & =\sum _{j=1}^N{b}_j(o_{t+1}){\beta }_{t+1}(j)a_{ij}& \\ & =\sum _{j=1}^N{a}_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)i=1,2,\cdots ,N& \end{array}\text{\hspace{1em}(10.20)}$$

(3)
$$\begin{array}{rl}P(O|\lambda )& =P(o_1,\cdots ,o_T|\lambda )\\ & =\sum _{i=1}^{N}P(o_1,\cdots ,o_T,i_1=q_i|\lambda )\\ & =\sum _{i=1}^{N}P(o_1,\cdots ,o_T|i_1=q_i,\lambda )P(i_1=q_i|\lambda )\\ & =\sum _{i=1}^{N}P(o_1|o_2,\cdots ,o_T,i_1=q_i,\lambda )P(o_2,\cdots ,o_T|i_1=q_i,\lambda )P(i_1=q_i|\lambda )\\ & =\sum _{i=1}^{N}P(o_1|i_1=q_i,\lambda )P(o_2,\cdots ,o_T|i_1=q_i,\lambda )P(i_1=q_i|\lambda )\\ & =\sum _{i=1}^N{b}_i(o_1){\beta }_1(i){\pi }_i\\ & =\sum _{i=1}^N{\pi }_i{b}_i(o_1){\beta }_1(i)\end{array}\text{\hspace{1em}(10.21)}$$
利用前向概率和后向概率的定义可以将观测序列概率$P(O|\lambda )$统一写成
$$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),t=1,2,\cdots ,T-1\text{\hspace{1em}(10.22)}$$
这是公式怎么来的呢，我们也推导一下。为了简便，省去固定参数$\lambda$。

$$\begin{array}{rlr}P(O)& =P(o_1,\cdots ,o_T)& \\ & =\sum _{i=1}^N\sum _{j=1}^{N}P(o_1,\cdots ,o_T,i_t=q_i,i_{t+1}=q_j)& \\ & =\sum _{i=1}^N\sum _{j=1}^{N}P(o_1,\cdots ,o_t,o_{t+2},\cdots ,o_T,i_t=q_i,i_{t+1}=q_j)P(o_{t+1}|o_1,\cdots ,o_t,o_{t+2},\cdots ,o_T,i_t=q_i,i_{t+1}=q_j)& \\ & =\sum _{i=1}^N\sum _{j=1}^{N}P(o_1,\cdots ,o_t,o_{t+2},\cdots ,o_T,i_t=q_i,i_{t+1}=q_j)P(o_{t+1}|i_{t+1}=q_j)& 观测独立假设\\ & =\sum _{i=1}^N\sum _{j=1}^{N}P(o_1,\cdots ,o_t,o_{t+2},\cdots ,o_T,i_t=q_i,i_{t+1}=q_j)b_j(o_{t+1})& \\ & =\sum _{i=1}^N\sum _{j=1}^{N}P(o_1,\cdots ,o_t,i_t=q_i,i_{t+1}=q_j)P(o_{t+2},\cdots ,o_T|o_1,\cdots ,o_t,i_t=q_i,i_{t+1}=q_j)b_j(o_{t+1})& \\ & =\sum _{i=1}^N\sum _{j=1}^{N}P(o_1,\cdots ,o_t,i_t=q_i,i_{t+1}=q_j)P(o_{t+2},\cdots ,o_T|i_{t+1}=q_j)b_j(o_{t+1})& D-划分\\ & =\sum _{i=1}^N\sum _{j=1}^{N}P(i_{t+1}=q_j|o_1,\cdots ,o_t,i_t=q_i)P(o_1,\cdots ,o_t,i_t=q_i)P(o_{t+2},\cdots ,o_T|i_{t+1}=q_j)b_j(o_{t+1})& \\ & =\sum _{i=1}^N\sum _{j=1}^{N}P(i_{t+1}=q_j|i_t=q_i)P(o_1,\cdots ,o_t,i_t=q_i)P(o_{t+2},\cdots ,o_T|i_{t+1}=q_j)b_j(o_{t+1})& 齐次马尔可夫假设\\ & =\sum _{i=1}^N\sum _{j=1}^N{a}_{ij}{\alpha }_t(i){\beta }_{t+1}(j)b_j(o_{t+1})& \\ & =\sum _{i=1}^N\sum _{j=1}^N{\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j),t=1,2,\cdots ,T-1& \end{array}$$

这里利用了前面介绍的D-划分：
$$P(o_{t+2},\cdots ,o_T|o_1,\cdots ,o_t,i_t=q_i,i_{t+1}=q_j)=P(o_{t+2},\cdots ,o_T|i_{t+1}=q_j)$$
根据概率图就很直观了：

节点$o_{t+2},\cdots ,o_T$中的任何一个节点到节点$o_1,\cdots ,o_t$或$i_t$都要经过$i_{t+1}$，而$i_{t+1}$被观测到(它在条件中)，所有的这种路径都是头到尾的，因此该条件独立性质成立。

一些概率与期望值的计算

利用前向概率和后向概率，可以得到关于单个状态和两个状态概率的计算公式。

1.给定模型$\lambda$和观测$O$，在时刻$t$处于状态$q_i$的概率。记
$${\gamma }_t(i)=P(i_t=q_i|O,\lambda )\text{\hspace{1em}(10.23)}$$
可以通过前向概率后向概率计算
$${\gamma }_t(i)=P(i_t=q_i|O,\lambda )=\frac{P(i_t=q_i,O|\lambda )}{P(O|\lambda )}$$

同样，为了简便，省去已知参数$\lambda$。

$$\begin{array}{rlr}{\gamma }_t(i)& =P(i_t=q_i|O)& \\ & =\frac{P(i_t=q_i,O)}{P(O)}& \\ & =\frac{1}{P(O)}P(i_t=q_i)P(O|i_t=q_i)& \\ & =\frac{P(i_t=q_i)}{P(O)}P(o_1,\cdots ,o_t|i_t=q_i)P(o_{t+1},\cdots ,o_T|i_t=q_i)& D-划分\\ & =\frac{1}{P(O)}P(o_1,\cdots ,o_t,i_t=q_i)P(o_{t+1},\cdots ,o_T|i_t=q_i)& \\ & =\frac{{\alpha }_t(i){\beta }_t(i)}{P(O)}& \\ & =\frac{{\alpha }_t(i){\beta }_t(i)}{{\sum }_{j=1}^{N}P(O,i_t=q_j)}& \\ & =\frac{{\alpha }_t(i){\beta }_t(i)}{{\sum }_{j=1}^N{\alpha }_t(j){\beta }_t(j)}& \end{array}\text{\hspace{1em}(10.24)}$$
这里关键部分是$P(O|i_t=q_i)=P(o_1,\cdots ,o_t|i_t=q_i)P(o_{t+1},\cdots ,o_T|i_t=q_i)$。

这个是怎么得到的呢？还是利用了D-划分规则：

可能还是没那么直观，再展开就好了：
$$\begin{array}{rlr}P(O|i_t=q_i)& =P(o_1,\cdots ,o_T|i_t=q_i)& \\ & =P(o_1,\cdots ,o_t|i_t=q_i)P(o_{t+1},\cdots ,o_T|i_t=q_i,o_1,\cdots ,o_t)& \\ & =P(o_1,\cdots ,o_t|i_t=q_i)P(o_{t+1},\cdots ,o_T|i_t=q_i)& d-划分\end{array}$$
这里$P(o_{t+1},\cdots ,o_T|i_t=q_i,o_1,\cdots ,o_t)=P(o_{t+1},\cdots ,o_T|i_t=q_i)$利用了D-划分。

节点$o_{t+1},\cdots ,o_T$ 中的任何一个节点到节点 $o_1,\cdots ,o_t$中任意节点都要经过$i_t$，而$i_t$被观测到，所有的这种路径都是头到尾的，因此该条件独立性质成立。

2.给定模型$\lambda$和观测$O$，在时刻$t$处于状态$q_i$且在时刻$t+1$处于状态$q_j$的概率。记
$${\xi }_t(i,j)=P(i_t=q_i,i_{t+1}=q_j|O,\lambda )\text{\hspace{1em}(10.25)}$$
可以通过前向后向概率计算：
$$\begin{array}{rl}{\xi }_t(i,j)& =\frac{P(i_t=q_i,i_{t+1}=q_j,O)}{P(O)}\\ & =\frac{P(i_t=q_i,i_{t+1}=q_j,O)}{{\sum }_{i=1}^N{\sum }_{j=1}^{N}P(i_t=q_i,i_{t+1}=q_j,O)}\end{array}$$
而$P(i_t=q_i,i_{t+1}=q_j,O)={\alpha }_t(i)a_{ij}{b}_j(o_{t+1}){\beta }_{t+1}(j)$我们以经在公式$(10.22)$中证明过，这里直接拿来使用。

所以
$${\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)}\text{\hspace{1em}(10.26)}$$