隐马尔科夫模型
生成模型与判别模型
在分类任务中，给定样本 $x$ ，需要对 $p(y|x)$ 建模，建模的方法分为生成方法与判别方法，对应的分别为生成模型与判别模型，生成模型直接对 $p(x,y)$ 建模，然后通过贝叶斯公式得到类别 $y$ 的后验分布即可：
\[p(y|x) = \frac{p(x,y)}{p(x)}\]
使用生成模型建模可以体现不同类型数据各自的特点，从统计的角度表示数据的分布情况，能够反映出同类数据本身的相似度，而不关心到底决策边界在哪里。最后找到后验概率最大的值对应的 $y$ 即为 $x$ 对应的标记。
判别模型则直接对 $p(y|x)$ 建模，或者学习决策函数 $y = f(x)$ ，目的是寻找异类数据的差异，通过决策平面来区分异类数据，两种建模方法的区别是：
- 生成模型可以还原出联合分布，而判别模型则不可以；
- 判别模型直接学习决策面，所以实践中一般准确率高于生成模型；
- 生成模型收敛速度较快，当数据有隐变量时只能用生成模型；
- 生成模型反应同类数据相似度，而判别模型只关心异类数据的差别；
典型的生成模型：Navie Bayes, HMM, Lattent Dirichlet Allocation ,PLSA, GMM；
典型的判别模型：Logistic Regression, SVM, KNN, Decision Tree, MaxEnt Model, Boosting, CRF.
隐马尔科夫模型
隐马尔科夫模型 HMM 是用于标注问题的监督学习模型，它是一种生成模型，描述一个隐藏的马尔科夫链随机生成状态序列，再由各个状态状态生成可观测序列的过程，HMM 的示意图如下所示：
下面把 HMM 形式化为数学模型，假设 Q 是所有可能的隐含状态集合， V 是所有可能的观测集合：
\[Q = \left \{q_1,q_2,…,q_N \right\}, \ \ V = \left \{v_1,v_2,…,v_M \right\}\]
一般 N $\ne$ M ，比如中文分词中采用 4 标记法则 $Q = \left \{ B,M,E,S \right \}$， 即 N = 4， V 则是语料库中字的数目。用 H 表示长度为 T 的隐态序列， O是对应的观测序列：
\[H = (h_1,h_2,…,h_T), \ \ O=(o_1,o_2,…,o_T)\]
用 A 表示状态转移矩阵，即状态转移概率构成一个 N $\times$ N 的矩阵：
\[A = [a_{ij}]_{N \times N}\]
其中 A 中的每个元素代表当前时刻状态 为 $q_i$ 下一时刻转到状态 $q_j$ 的概率：
\[a_{ij} = P(h_{i+1} = q_j | h_{t} = q_i)\]
B 是观测概率矩阵：
\[B = [b_j(k)]_{N \times M}\]
B 中的每个元素代表了在 t 时刻处于状态 $q_j$ 的条件下生成观测态 $v_k$ 的概率。
\[ b_j(k) = P(o_t = v_k | h_t = q_j ) \]
$\pi$ 是初始状态概率向量：
\[ \pi = \left \{ \pi_i \right \}_{i=1}^N , \ \ \pi_i = P(h_1 = q_i) \]
综上，初始状态概率向量 $\pi$ ，状态转移概率矩阵 A 与观测概率矩阵 B 共同构成了隐马尔科夫模型 $\lambda$：
\[\lambda = \left \{A,B ,\pi \right \} \]
HMM 中 $\pi$ 与 A 决定了马尔科夫链， B 决定了观测态的生成，模型有两个比较强的假设：
齐次马尔科夫假设：马尔科夫链中在任意时刻 t 的状态依赖于前一时刻的状态，与观测态和时刻 t 无关，即：
\[P(h_t | h_{t-1},o_{t-1},…,h_1,o_1) = P(h_t|h_{t-1})\]
观测独立假设：任意时刻的观测状态只依赖于当前的隐态，与其他状态无关：
\[P(o_t | h_T,o_T,…,h_1,o_1) = P(o_t|h_t)\]
HMM 三个基本问题
HMM 模型可以解决三个基本问题，分别如下所示：
概率计算问题：给定模型 $\lambda = \left \{A,B ,\pi \right \}$ ，和观测序列 $O=(o_1,o_2,…,o_T)$, 计算模型 $\lambda$ 中观测序列 O 出现的概率 $P(O|\lambda)$。
学习问题：已知观测序列 $O=(o_1,o_2,…,o_T)$， 估计模型参数 $ \lambda = \left \{A,B ,\pi \right \}$， 使得$P(O|\lambda)$ 在该模型下出现的概率最大， 如果同时给定隐态序列 H ，则直接 MLE 即可估计参数，否则采用 EM 算法即可。
预测问题：给定模型 $ \lambda = \left \{A,B ,\pi \right \}$ 和预测序列 $O=(o_1,o_2,…,o_T)$ ，求给定观测序列的条件下最有可能的隐态序列 $H = (h_1,h_2,…,h_T)$ ，即对应的 $P(H|O)$ 最大的序列即可。
概率计算问题
给定模型 $\lambda = \left \{A,B ,\pi \right \}$ ，和观测序列 $O=(o_1,o_2,…,o_T)$，最显而易见的做法是根据边缘概率与条件概率公式求得观测序列的边缘分布：
\[P(O|\lambda) = \sum_H P(H,O|\lambda) = \sum_H P(O|H,\lambda)P(H|\lambda)\]
根据上述公式，只需枚举所有可能的状态序列然后求和即可，但是该方法计算量非常大，是 $O(TN^T)$ 阶的，所以该方法只是理论上可行，一种比较简便的做法便是前向算法，在前向算法中有一个前向概率的概念，对于当前模型 $\lambda$ ，定义在时刻 t 的观测序列为 $o_1,o_2,…,o_t$ ，且当前状态为 $q_i$ 的概率为前向概率，形式如下：
\[a_t(i) = P(o_1,o_2,…,o_t,h_t = q_i|\lambda)\]
前向概率可以递推求解，首先给定初值：
\[a_1(i) = \pi_ib_i(o_1)\]
用 $a_t(j)$ 表示时刻 t 观测序列为 $o_1,o_2,…,o_t$ 且在时刻 t 处于状态 $q_j$ 的前向概率，那么 $a_t(j)a_{ji}$ 则代表了在时刻 t 处于状态 $q_j$ 而在 时刻 $t+1$ 处于 $q_i$ 的联合概率，对时刻 $t$ 所有可能的 N 个状态求和 就是到时刻 t 观测为 $o_1,o_2,…,o_t$ 并在时刻 t 处于状态 $q_j$ 且时刻 $t+1$ 到达状态 $q_i$ 的联合概率，最后乘以观测概率 $b_i(o_{t+1})$ 得到 $a_{t+1}(i)$ 即可：
\[a_{t+1}(i) = \left \{ \sum_{j=1}^N a_t(j)a_{ji} \right \} b_i(o_{t+1}) \tag{*}\]
$(*)$ 式完整的计算过程如下：
\begin{aligned}
a_{t+1}(i)
&= P(o_1,..,o_{t+1},h_{t+1} = q_i) = P(o_1,..,o_{t+1}|h_{t+1} = q_i)P(h_{t+1} = q_i) \\ & \ \\
&= P(o_1,..,o_t|h_{t+1} = q_i)P(o_{t+1}|h_{t+1} = q_i)P(h_{t+1} = q_i) = P(o_1,..,o_t,h_{t+1} = q_i)P(o_{t+1}|h_{t+1} = q_i)\\
&= \left \{\sum_{j=1}^NP(o_1,..,o_t,h_t = q_j,h_{t+1} = q_i)\right\} P(o_{t+1}|h_{t+1} = q_i)\\
&= \left \{\sum_{j=1}^NP(o_1,..,o_t,h_{t+1} = q_i|h_t = q_j)P(h_t = q_j)\right\}P(o_{t+1}|h_{t+1} = q_i)\\
&= \left \{\sum_{j=1}^NP(o_1,..,o_t,h_t = q_j)P(h_{t+1} = q_i|h_t = q_j)\right\}P(o_{t+1}|h_{t+1} = q_i)\\
&= \left \{\sum_{j=1}^Na_t(j)a_{ji} \right\}b_i(o_{t+1})
\end{aligned}
有了以上递推公式，就可以一直计算下去，一直到时刻 T ：
\[a_T(i) = P(o_1,o_2,…,o_T,h_T = q_i|\lambda)\]
因而得到：一直递推下去，最后一列的部分概率是所有可能的路径的概率和，所以就是这个观察序列 O 在给定 $\lambda$ 下的概率了：
\[P(O|\lambda) = \sum_{i=1}^N a_T(i) \]
相较于枚举所有序列，前向算法递归的利用之前时刻的概率值，避免了重复计算，在 $t = 1$ 时刻，计算 $a_1(i)$ 的N 个值，接下来每次计算 $a_{t+1}(i)$ 时都会利用前一时刻的 N 个值，计算的时间复杂度为 $O(N^2T)$ ，如下如所示：
前向算法就可以解决观测序列 O 出现的概率 $P(O|\lambda)$ ，与前向算法相对应的还有一个后向算法定义时刻 t 的状态为 $q_i$ 的条件下，从 $t+1$ 到 T 的观测序列为 $o_{t+1},o_{t+2},…,o_T$ 的概率为后向概率，记做：
\[\beta_t(i) = P(o_{t+1},o_{t+2},…,o_T|h_t = q_i,\lambda)\]
首先对于所有状态 $q_i$ 定义：
\[ \beta_T(i) = 1 ,i=1,2,…,N \]
这里，要计算 $\beta_t(i)$，只需考虑在 $t+1$ 时刻所有可能的 N 个状态 $q_j$ 的转移概率，即 $a_{ij}$ ，以及在此状态下观测状态为 $o _{t+1}$ 的概率 $b_j(o_{t+1})$ ，不断从后往前递推即可：
\[\beta_t(i) = \sum_{j=1}^Na_{ij}b_j(o_{t+1})\beta_{t+1}(j)\]
最后类似于前向算法的求和：
\[P(O|\lambda) = \sum_{i=1}^N\pi_i b_i(o_1)\beta_1(i)\]
最后给出一个前向概率与后向概率的图，如下所示，左边为前向概率，右边为后向概率：
前向概率与后向概率结合，可以得到一系列重要的结论：
结论1： 给定 $\lambda$ 和 观测序列 O ，则在时刻 t 处于状态 $q_i$ 的概率为：
\[\gamma_t(i) = P(h_t =q_i|O,\lambda) = \frac{P(h_t = q_i,O,\lambda)}{P(O|\lambda)}\]
由前向概率 $a_t(i)$ 与后向概率 $\beta_t(i)$ 可知：
\[ a_t(i)\beta_t(i) = P(h_t = q_i,O|\lambda) \]
于是可得：
\[\gamma_t(i) = \frac{a_t(i) \beta_t(i)}{P(O|\lambda)}=\frac{a_t(i)\beta_t(i)}{\sum_{j=1}^Na_t(j)\beta_t(j)}\]
结论2： 给定 $\lambda$ 和 观测序列 O ，则在时刻 t 处于状态 $q_i$ ，t+1 时刻状态为 $q_j$ 的概率；
\[\xi _t(i,j) = P(h_t = q_i ,h_{t+1} = q_j|O，\lambda)\]
用前后向概率来表示：
\[P(h_t = q_i,h_{t+1}=q_j,O,|\lambda) = a_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)\]
因此可得
\begin{aligned}
\xi_t(i,j) &= \frac{P(h_t = q_i ,h_{t+1} = q_j|O,\lambda)}{P(O|\lambda)}\\
&= \frac{P(h_t = q_i ,h_{t+1} = q_j|O,\lambda)}{\sum_{i=1}^N\sum_{j=1}^N P(h_t = q_i ,h_{t+1} = q_j|O,\lambda) }\\
&=\frac{a_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)}{\sum_{i=1}^N\sum_{j=1}^Na_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j) }\\
\end{aligned}
结论3：前向后向算法是已知模型和观测序列求观测序列出现概率的算法，也是用于模型参数学习的 Baum-Welch 算法的循环中的一个步骤，有几个关键的量需要计算，将 $\gamma_t(i)$ 与 $\xi_t(i,j)$ 对各个时刻 t 求和，可以得到几个期望值：
1）在观测序列 O 下状态 i 出现的期望为：$\sum_{t=1}^T \gamma_t(i)$；
2) 在观测序列 O 下由状态 i 转移的期望值为：$ \sum_{t=1}^{T-1}\gamma_t(i) $；
3）在观测序列 O 下由状态 i 转移到状态 j 的期望为：$\sum_{t=1}^{T-1}\xi_t(i,j)$.
参数学习
HMM 的参数学习算法分为监督学习与非监督学习算法，给定数据集 $\left \{ (H_i,O_i)\right\}_{i=1}^S$，监督学习的方式直接采用 MLE 来估计模型的参数 $\lambda = \left \{A,B ,\pi \right \}$ ，设样本在时刻 t 处于状态 i 而在 t+1 转移到状态 j 的频数为$A_{ij}$，而样本中状态为 j 并且观测为 k 的频数记做 $B_{jk}$ ，则 $\hat{a}_{ij}$ 与 $\hat{b}_j(k)$ 的估计值分别为：
\[\hat{a}_{ij} = \frac{A_{ij}}{\sum_{j=1}^NA_{ij}} , \ \ \ \ \hat{b}_j(k) =\frac{B_{jk}}{\sum_{k=1}^MB_{jk}}, \ \ i,j=1,...,N\]
初始状态 $\pi_i$ 的估计值 $\hat{\pi}_i$ 为S 个样本中初始状态为 $q_i$ 的频率，这便得到了MLE 的模型 $\hat{\lambda}$。对于无标记数据，只会有一些观测数据 $\left\{ O_1,O_2,…,O_S\right\}$， 而没有对应的隐态序列，这时可以考虑 O 为可见变量， H 为隐藏变量，这就变成了一个含有隐变量的概率模型：
\[P(O|\lambda) = \sum_HP(O|H,\lambda)P(H|\lambda)\]
其参数学习可以用 Baum-Welch 算法来实现 ，Baum-Welch 是一种特殊形式的 EM 算法来实现，其 Q 函数的形式与 EM 稍有差别，首先根据：
\[P(H|O,\lambda^{old})= \frac{P(O,H|\lambda^{old}) }{P(O|\lambda^{old})}\]
在 Baum-Welch 算法里，省去了 对 $\lambda$ 来说的常数因子，$1/P(O|\lambda^{old})$ ，因此其 Q 函数为：
\[Q(\lambda|\lambda^{old}) = \sum_H \log P(O,H|\lambda) P(O,H|\lambda^{old})\]
完全数据的联合分布为：
\[P(O,H|\lambda) = \pi_{h_1}b_{h_1}(o_1)a_{h_1h_2}b_{h_2}(o_2),…,a_{h_{T-1}h_T}b_{h_T}(o_T)\]
于是 Q 函数可以写作：
\begin{aligned}
Q(\lambda,\lambda^{old})
= \sum_H \log \pi_{h_1}P(O,H|\lambda^{old}) +\sum_H\left \{\sum_{t=1}^{T-1}\log a_{h_th_{t+1}} \right \}P(O,I|\lambda^{old})+\sum_H\left \{\sum_{t=1}^{T}\log b_{h_t}(o_t) \right \}P(O,I|\lambda^{old})\end{aligned}
对于这个复杂的 Q 函数，需逐项求解，第一项可写作：
\[\sum_H \log \pi_{h_1}P(O,H|\lambda^{old}) = \sum_{i=1}^N \log \pi_i P(O,h_1 = q_i| \lambda^{old})\]
这里 $\pi_i$ 满足约束 $\sum_{i=1}^N\pi_i = 1$ ，这便是一个等式约束的拉格朗日优化问题，其拉格朗日函数为：
\[L(\mu) = \sum_{i=1}^N \log_{\pi_i} P(O,h_1 = q_i|\lambda) + \mu \left( \sum_{i=1}^N\pi_i -1 \right)\]
对 $\pi_i$ 求导，令导数得 0 ，进而可得：
\begin{aligned}
\frac{\partial L(\mu)}{\partial \pi_i} = 0 \Rightarrow & \ P(O,h_1 = q_i|\lambda^{old}) + \mu \pi_i = 0 \\
\Rightarrow & \ \mu \sum_{i=1}^N \pi_i =-\sum_{i=1}^N P(O,h_1 = q_i|\lambda^{old}) \\
\Rightarrow & \ \mu =- P(O|\lambda^{old})
\end{aligned}
带入 $L(\mu)$ 即可得：
\[\pi_i = \frac{P(O,h_t = q_i|\lambda)}{P(O|\lambda^{old})}\]
Q 函数的第二项可以写作：
\[\sum_H\left \{\sum_{t=1}^{T-1}\log a_{h_th_{t+1}} \right \}P(O,I|\lambda^{old}) = \sum_{i=1}^N\sum_{j=1}^N\sum_{t=1}^{T-1} \log a_{ij} P(O,h_t = q_i,h_{t+1} = q_j| \lambda^{old})\]
类似的根据约束 $\sum_{j=1}^N a_{ij} = 1$ ，构建拉格朗日函数，即可求得：
\[a_{ij} = \frac{\sum_{t=1}^{T-1} P(O,h_t = q_i,h_{t+1} = q_j|\lambda^{old})}{\sum_{t=1}^{T-1} P(O,h_t = q_i|\lambda^{old})}\]
Q 函数的第三项可以写作：
\[\sum_H\left \{\sum_{t=1}^{T}\log b_{h_t}(o_t) \right \}P(O,I|\lambda^{old})= \sum_{j=1}^N\sum_{t=1}^T \log b_j(o_t)P(O,h_t = q_j|\lambda^{old})\]
同样应用拉格朗日乘子法，约束条件为 $\sum_{k=1}^Mb_j(k) = 1$ ，这里只有 $o_t = v_k$ 时， $b_j(o_t)$ 对 $b_j(k)$ 的偏导数才不为 0 ，用 $I(o_t = v_k)$ 表示即可：
\[b_j(k) = \frac{\sum_{t=1}^T P(O,h_t = q_j|\lambda^{old})I(o_t = v_k)}{\sum_{t=1}^T P(O,h_t = q_j|\lambda^{old})}\]
这里用之前前后向的概率来表示各个 HMM 中的参数 $\lambda ^{new}$ ：
\[a_{ij}^{new} = \frac{\sum_{t=1}^{T-1}\xi_t(i,j) }{\sum_{t=1}^{T-1}\gamma_t(i)} ,\ \ b_j(k)^{new} = \frac{\sum_{t=1,o_t = v_k}^T\gamma_t(j) }{\sum_{t=1}^T\gamma_t(j)} ,\ \ \pi_i^{new} = \gamma_1(i) \]
这便是 Baum-Welch 算法，跟 EM 算法思想一模一样，首先初始化参数 $\lambda^{old} = \left\{ A^{old},B^{old} ,\pi^{old} \right\}$ ，按照以上步骤递推计算 $\lambda^{new}$ 即可，迭代终止便得到了最终模型的参数.
预测算法
HMM 的预测算法是这样的，给定模型 $ \lambda = \left \{A,B ,\pi \right \}$ 和预测序列 $O=(o_1,o_2,…,o_T)$ ，求给定观测序列的条件下最有可能的隐态序列 $H = (h_1,h_2,…,h_T)$ ，最直接的做法就是根据前后向算法得到近似结果：
\[\gamma_t(i) = \frac{a_t(i)\beta_t(i)}{P(O|\lambda) } = \frac{a_t(i)\beta_t(i)}{\sum_{j=1}^Na_t(j)\beta_t(j)}\]
则在每个时刻最可能的状态时：
\[h_t^* = \arg\max_{1 \le i \le N}[\gamma_t(i)]\]
从而得到状态序列 $H^* = \left\{h_1^* ,h_2^* ,…,h_T^* \right\}$，但是这种方式只是针对单个节点计算，甚至得到的序列中存在两个相邻的状态转移概率为 0 的情况，所以通常 HMM 的解码算法采用的是维特比算法，该算法是针对篱笆网络的有向图（Lattice )的最短路径问题而提出的。篱笆网络就是类似于深度学习那样的有向网络，维特比算法是一种动态规划算法，凡是使用隐含马尔可夫模型描述的问题都可以用它来解码，包括数字通信、语音识别、机器翻译、拼音转汉字、分词等问题。算法的精髓就是，既然知道时刻 $t$ 所有隐状态 $h_t = q_i, \ i = 1,2,…,N$ 的最短路径，那么 t+1 时刻的最短路径就等于时刻 t 的第 j 个节点的最短路径加上该状态到第 t+1 列各个状态的距离的最小值。这里路径便为对应的隐藏状态序列，每个状态都有一个最可能的路径，比如对于有三个隐态的序列，在 $T=3$ 时刻的三个状态都有一个如下的最可能的路径，找到三个路径中中最终节点的代价最小值（代表全局最优路径），回溯即可找到整条路径：
和前向算法中的部分概率不一样，这里的概率只是单条最优路径的概率，而不是所有路径的概率和，接下来定义在时刻 $t$ 状态为 $q_i$ 的所有单个路径 $(h_1,h_2,…,h_t)$ 中概率最大的:
\[\delta_t(i) = \max_{h_1,h_2,…,h_{t-1}}P(h_t = q_i,h_{t-1},…,h_1,o_t,…,o_1|\lambda)\]
由此可得变量 $\delta$ 的递推公式：
\begin{aligned}
\delta_{t+1}(i) &= \max_{h_1,...,h_t}P(h_{t+1} = q_i,h_t,...,h_1,o_{t+1},...,o_1|\lambda)\\
&= \max_{1\le j\le N}\left \{ \delta_t(j)a_{ji} \right \}b_i(o_{t+1}), \ \ i=1,2,...,N,\ t=1,2,...,T-1
\end{aligned}
考虑到要计算 t 时刻的部分概率，只需要知道 t-1 时刻的部分概率，所以只需要记录那个导致了 t 时刻最大部分概率的的状态，定义在时刻 $t$ 的状态为 $q_i$ 的所有单个路径 $\left\{h_1 ,h_2 ,…,h_{t-1},h_t = q_i \right\}$ 中概率最大的路径的第 $t-1$ 个节点为：
\[ \Psi_i(t) = \arg \max_{1 \le j \le N} \left \{ \delta _{t-1} (j) a_{ji}\right \} , \ \ \ i = 1,2,…,N\]
综上，最后总结一下 维特比算法：
输入：给定模型 $\lambda = \left\{A,B,\pi \right\}$ 与观测序列 $O =\left\{o_1 ,o_2 ,…,o_t \right\}$;
输出：最优路径 $H ^* = \left\{h_1^* ,h_2^* ,…,h_t^* \right\}$.
1. 初始化： $\delta_1(i) = \pi_ib_i(o_1), \ \ \ \ \Psi_1(i) = 0 , \ \ \ i = 1,2,…,N$；
2.对 $t = 2,3,…,T$ ：
\begin{aligned}
\delta_t(i) &= \max_{1\le j\le N}\left \{ \delta_{t-1}(j)a_{ji} \right \}b_i(o_t), &\ \ \ i=1,2,...,N\\
\Psi_i(t)& = \arg \max_{1 \le j \le N} \left \{ \delta _{t-1} (j) a_{ji}\right \}, &\ \ \ i=1,2,...,N
\end{aligned}3. 终止时，得到发生概率的最大值：
\[P^* = \max_{1\le i\le N}\delta_T(i)\]
\[h^*_T = \arg \max_{1\le i\le N}\delta_T(i)\]
4. 回溯找到整条最优路径，对 $t = T-1,T-2,…1$ :
\[h_i^ * = \Psi_{t+1}(h^*_{t+1})\]
这便是维特比的过程，因为没时间也没有写实现，但是实现相对简单，Hancks 那里有代码，时间都去哪了,还没好好感受年轻就老了…
参考：
http://www.cnblogs.com/xiapeiyong/archive/2009/05/15/1457888.html
统计学习方法,PRML
http://freemind.pluskid.org/machine-learning/discriminative-modeling-vs-generative-modeling/
