隐马尔科夫模型（二）
回忆上次所提到的，隐马尔科夫模型引出的三个基本问题：
1)概率计算问题：给定模型参数$\lambda=(A,B,\pi)$ 以及某个观测序列$O=(0_{1},...,0_{T})$, 如何高效地计算出概率$P(O\mid \lambda)$?
2)学习问题：已知某个观测序列$O=(o_{1},...,o_{T})$, 如何得到参数$\lambda=(A,B,\pi)$使得似然函数$P(O\mid\lambda)$极大，或者至少尽量大?
3)预测问题：预测问题，也成为解码问题(Decoding Problem)，也就是已知观测序列$O=(o_{1},...,o_{T})$, 以及参数$\lambda=(A,B,\pi)$, 求不可观测序列$I$使得概率: \begin{equation}P(I\mid O,\lambda)\end{equation}极大。
这一次就总结一下后面两个问题的算法。
1.学习算法：Baum-Welch算法
将Baum-Welch算法就是将EM算法应用于HMM的参数学习。这时候我们假设某带未知参数$\lambda=(A,B,\pi)$有一可观测序列$O=(o_{1},...,o_{T})$，现在我们开始计算E步。
在$E$步中，假设我们第(s)次迭代已经有某个参数$\lambda^{(s)}=(A^{(s)},B^{(s)},\pi^{(s)})$, 现在我们要计算的是条件期望：
\begin{equation}\begin{split}&E_{I}(\log P(O,I\mid \lambda)\mid O,\lambda^{(s)})\newline =&\sum_{I}P(I\mid O,\lambda^{(s)})\log P(O,I\mid \lambda)\newline =&\sum_{I}\frac{P(O,I\mid \lambda^{(s)})}{P(O\mid \lambda^{(s)})}\log P(O,I\mid \lambda)\end{split}\end{equation}
由于$P(O\mid \lambda^{(s)})$为常数，求上述期望的极小值等价于求$P(O\mid \lambda^{(s)})\cdot E_{I}(\log P(O,I\mid \lambda)\mid O,\lambda^{(s)})$之极小值，这时为了方便起见我们令函数$Q(\lambda,\lambda^{(s)})\triangleq P(O\mid \lambda^{(s)})\cdot E_{I}(\log P(O,I\mid \lambda)\mid O,\lambda^{(s)})=\sum_{I}P(O,I\mid \lambda^{(s)})\log P(O,I\mid \lambda)$, 而不是直接将$Q$定义为条件期望$E_{I}(\log P(O,I\mid \lambda)\mid O,\lambda^{(s)})$。
注意到，我们直接计算概率：
$$P(O,I\mid \lambda)=\prod_{t=1}^{T}B_{o_{t}i_{t}}\cdot\prod_{t=1}^{T-1}A_{i_{t+1}i_{t}}\cdot \pi_{i_{1}},$$
从而：
\begin{equation}\log P(O,I\mid\lambda)=\sum_{t=1}^{T}\log B_{o_{t}i_{t}}+\sum_{t}^{T-1}\log A_{i_{t+1}i_{t}}+\log\pi_{i_{1}},\end{equation}
结合定义我们立即得到：
\begin{equation}Q(\lambda,\lambda^{(s)})=\sum_{I}\sum_{t=1}^{T}P(O,I\mid\lambda^{(s)})\log B_{o_{t}i_{t}}+\sum_{I}\sum_{t=1}^{T-1}P(O,I\mid\lambda^{(s)})+\sum_{I}P(O,I\mid \lambda^{(s)})\log\pi_{i_{1}}.\end{equation}
计算上式第一项：
\begin{equation}\begin{split}&\sum_{I}\sum_{t=1}^{T}P(O,I\mid \lambda^{(s)})\log B_{o_{t}i_{t}}\newline =&\sum_{t=1}^{T}\sum_{k=1}^{K}\sum_{I:i_{t}=k}P(O,I\mid \lambda^{(s)})\log B_{o_{t}k}\newline =&\sum_{t=1}^{T}\sum_{k=1}^{K}P(O,i_{t}=k\mid\lambda^{(s)})\log B_{o_{t}k}.\end{split}\end{equation}
类似地，第二，第三项分别为$\sum_{k,l}\sum_{t=1}^{T-1}\lbrace P(O,i_{t}=k,i_{t+1}=l\mid\lambda^{(s)})\rbrace\log A_{lk},$ $\sum_{k=1}^{K}P(O,i_{1}=k\mid\lambda^{(s)})\log\pi_{k},$
所以我们自然有：
\begin{equation}\begin{split}Q(\lambda,\lambda^{(s)})=&\sum_{t=1}^{T}\sum_{k=1}^{K}P(O,i_{t}=k\mid\lambda^{(s)})\log B_{o_{t}k}+\sum_{k,l}\sum_{t=1}^{T-1}\lbrace P(O,i_{t}=k,i_{t+1}=l\mid\lambda^{(s)})\rbrace\log A_{lk}+\newline &\sum_{k=1}^{K}P(O,i_{1}=k\mid\lambda^{(s)})\log\pi_{k},\end{split}\end{equation}
现在E步完成。
M步中，我们很容易得到当$Q(\lambda,\lambda^{(s)})$最大的时候，有：
\begin{equation}A_{lk}=\frac{\sum_{t=1}^{T-1}P(O,i_{t}=k,i_{t+1}=l\mid \lambda^{(s)})}{\sum_{t=1}^{T-1}P(O,i_{t}=k\mid \lambda^{(s)})}\end{equation}
\begin{equation}B_{\alpha k}=\frac{\sum_{t:o_{t}=\alpha}P(O,i_{t}=k\mid\lambda^{(s)})}{\sum_{t=1}^{T}P(O,i_{t}=k\mid\lambda^{(s)})}\end{equation}
\begin{equation}\pi_{i}=\frac{P(O,i_{1}=i\mid \lambda^{(s)})}{P(O\mid \lambda^{(s)})}\end{equation}
算法一（Baum-Welch算法）
输入：观测序列$O=(o_{1},...,o_{T})$；
输出：隐马尔科夫模型参数$\lambda=(A,B,\lambda)$
Step1. 初始化参数$\lambda^{(0)}=(A^{(0)},B^{(0)},\lambda^{(0)})$;
Step2. 对$s=1,...$不断执行如下操作，直到停止条件满足：
1) 对$t=1,...,T$,$k=1,...K$递归计算：
$\alpha_{t,k}^{(s)}=B_{o_{1}k}^{(s)}\pi_{k}^{(s)}$ 当$t=1$,
$\alpha_{t,k}^{(s)}=\sum_{l=1}^{K}B_{o_{t}l}^{(s)}A_{lk}^{(s)}\alpha_{t-1,k}^{(s)}$当$t>1$
2) 对$t=T,T-1,...,1$,$k=1,...,K$递归计算：
$\beta_{t,k}^{(s)}=1$ 当$t=T$,
$\beta_{t,k}^{(s)}=\sum_{l=1}^{K}\beta_{t+1,l}^{(s)}B_{o_{t+1}l}^{(s)}A_{lk}^{(s)}$
3) 计算概率：
对$t=1,...,T-1$, $k,l=1,...,K$:
$$P(O,i_{t}=k,i_{t+1}=l\mid \lambda^{(s)})=\beta_{t+1,l}^{(s)}B_{o_{t+1}l}^{(s)}A_{lk}^{(s)}\alpha_{t,k}^{(s)}$$
对$t=1,...,T$,$k=1,...,K$:
$$P(O,i_{t}=k\mid\lambda^{(s)})=\beta_{t,k}^{(s)}\alpha_{t,k}^{(s)}$$
计算：
$$P(O\mid\lambda^{(s)})=\sum_{k=1}^{K}\alpha_{T,k}^{(s)}$$
4) 计算：
\begin{equation}A_{lk}^{(s+1)}=\frac{\sum_{t=1}^{T-1}P(O,i_{t}=k,i_{t+1}=l\mid \lambda^{(s)})}{\sum_{t=1}^{T-1}P(O,i_{t}=k\mid \lambda^{(s)})}\end{equation}
\begin{equation}B_{\alpha k}^{(s+1)}=\frac{\sum_{t:o_{t}=\alpha}P(O,i_{t}=k\mid\lambda^{(s)})}{\sum_{t=1}^{T}P(O,i_{t}=k\mid\lambda^{(s)})}\end{equation}
\begin{equation}\pi_{i}^{(s+1)}=\frac{P(O,i_{1}=i\mid \lambda^{(s)})}{P(O\mid \lambda^{(s)})}\end{equation}
Step3 输出最终得到的参数$\lambda^{(S)}=(A^{(S)},B^{(S)},\pi^{(S)})$, 其中$S$为学习步长。
2.预测问题：维特比算法
维特比算法(Viterbi Algorithm)就是用动态规划法求最短路径$I$使得$P(I\mid O,\lambda)$极大(等价于$P(O,I)$极大)。首先我们观察一个事实是，如果不可观测序列（也称状态路径）$I=(i_{1},...,i_{T})$是一个最优路径的话，那么对于任意$t=1,...,T$，前$t$个状态组成的路径$I_{t}=(i_{1},...,i_{t})$也必定是一个固定$t$时刻的状态下的最优路径，具体来说就是：
\begin{equation}(i_{1},...,i_{t})=\mathop{\text{argmax}}\limits_{\lbrace(i^{\prime}_{1},...,i^{\prime}_{t})\mid i^{\prime}_{s}\in\lbrace 1,...,K\rbrace, i_{t}^{\prime}=i_{t}\rbrace}P(O,i_{1}^{\prime},...,i_{t}^{\prime},i_{t+1},...,i_{T}\mid\lambda).\end{equation}
注意到对于满足$t_{t}^{\prime}=i_{t}$之序列我们又有： \begin{equation}P(O,i_{1}^{\prime},...,i_{t}^{\prime},i_{t+1},...,i_{T}\mid\lambda)=P(i_{t+1},...,i_{T},o_{t+1},...,o_{T}\mid i_{t})P(o_{1},...,o_{t},i_{1}^{\prime},...,i_{t}^{\prime}\mid\lambda),\end{equation}
所以(13)等价于：
\begin{equation}(i_{1},...,i_{t})=\mathop{\text{argmax}}\limits_{\lbrace(i^{\prime}_{1},...,i^{\prime}_{t})\mid i^{\prime}_{s}\in\lbrace 1,...,K\rbrace, i_{t}^{\prime}=i_{t}\rbrace}P(o_{1},...,o_{t},i_{1}^{\prime},...,i_{t}^{\prime}\mid\lambda).\end{equation}
(15)可以启发我们设计一个动态规划算法如下：
我们现在想递归地求：
\begin{equation}\delta_{t,i}\triangleq \mathop{\max}\limits_{(i_{1},...,i_{t-1})}P(o_{1},...,o_{t},i_{1},...,i_{t-1},i_{t}=i\mid\lambda)\end{equation}
以及： \begin{equation}\psi_{t,i}\triangleq (\mathop{\text{argmax}}\limits_{(i_{1},...,i_{t-1})}P(o_{1},...,o_{t},i_{1},...,i_{t-1},i_{t}=i\mid\lambda))_{t-1},\end{equation}
注意的这里$\mathop{\text{argmax}}\limits_{(i_{1},...,i_{t-1})}P(o_{1},...,o_{t},i_{1},...,i_{t-1},i_{t}=i\mid\lambda)$是一个长度为$t-1$的状态序列（路径），而上面公式的下标表示取该序列的最后一个状态，也就是路径终点。上面的定义式中$\delta_{t,i}$表示前$t$步最优概率值，当固定第$t$步状态为$i$的情况下；$\psi_{t,i}$表示该最优解的倒数第二步的状态，记录该状态以便于回溯。
我们很容易得到递推关系：
\begin{equation}\delta_{t+1,i}=\max_{j} B_{o_{t+1}i}A_{ij}\delta_{t,j}\end{equation}
\begin{equation}\psi_{t+1,i}=\mathop{\text{argmax}}\limits_{j}A_{ij}\delta_{t,j}\end{equation}
如果以此由上式递推得到$\delta$,$\psi$, 我们立即可以得到最优路径$I^{\ast}=(i^{\ast}_{1},...,i^{\ast}_{T})$：
\begin{equation}i^{\ast}_{T}=\mathop{\text{argmax}}\limits_{j}\delta_{T,j},\end{equation}
\begin{equation}i^{\ast}_{t}=\psi_{t+1,i^{\ast}_{t+1}}(t=T-1,...,1),\end{equation}
综上我们总结一下：
算法二（维特比算法）：
输入：观测序列$O=(o_{1},...,o_{T})$, 参数$\lambda=(A,B,\pi)$;
输出：预测最优路径$I^{\ast}=(i^{\ast}_{1},...,i^{\ast}_{T})$.
Step1. 初始化：$\delta_{1,i}=B_{o_{1},i}\pi_{i}$, $\psi_{1,i}=0$
Step2. 对$t=2,...T$, $i=1,...,K$递推得到：
$$\delta_{t,i}=\max_{j}B_{o_{t}i}A_{ij}\delta_{t-1,j},$$
$$\psi_{t,i}=\mathop{\text{argmax}}\limits_{j}A_{ij}\delta_{t-1,j}$$
Step3. 回溯得到最优路径：
$$i^{\ast}_{T}=\mathop{\text{argmax}}\limits_{i}\delta_{T,i}$$,
对$t=T-1,...,1$:
$$i^{\ast}_{t}=\psi_{t+1,i^{\ast}_{t+1}}$$,
输出$I^{\ast}=(i^{\ast}_{1},...,i^{\ast}_{T})$。
