【统计学习方法死磕系列】隐马尔可夫模型

写在前面

隐马尔可夫模型（Hidden Markov Model, HMM）是概率图模型中最经典的动态模型之一。它广泛应用于语音识别、自然语言处理（分词、词性标注、命名实体识别）、生物信息学（基因序列分析）、金融时序分析等领域。

HMM 的精妙之处在于：它用两层结构（隐藏的状态序列 + 可观测的输出序列）来描述一个动态系统。我们只能观察到输出，而隐藏的状态需要从观测中推断。这种”隐变量”建模思想是概率图模型、EM 算法乃至现代变分推断的基础。

本文将以数学严格性死磕 HMM 的三个经典问题（评估、解码、学习），并结合前向-后向算法、Viterbi 算法、Baum-Welch 算法（EM）进行全面推导。建议读者在阅读过程中用笔推演每个关键公式。

1. 马尔可夫链基础

在深入 HMM 之前，我们需要先建立马尔可夫链（Markov Chain）的数学基础。

1.1 马尔可夫性质

定义 1（马尔可夫性质，Markov Property）：一个随机过程 ${X_t, t = 1, 2, \ldots}$ 具有马尔可夫性质，当且仅当：

$$P(X_{t+1} = s_{t+1} | X_t = s_t, X_{t-1} = s_{t-1}, \ldots, X_1 = s_1) = P(X_{t+1} = s_{t+1} | X_t = s_t)$$

也就是说，给定当前状态，未来状态与过去状态条件独立。”未来只依赖现在，不依赖过去。”

这个性质极其强大，因为它意味着整个系统的动态可以完全由一步转移概率来描述。

$m$ 阶马尔可夫性质：$P(X_{t+1} | X_t, X_{t-1}, \ldots) = P(X_{t+1} | X_t, X_{t-1}, \ldots, X_{t-m+1})$——未来只依赖最近的 $m$ 个状态。本文只讨论一阶情况（$m=1$），这是 HMM 的标准假设。

1.2 离散时间马尔可夫链

考虑有限状态空间 $\mathcal{S} = {s_1, s_2, \ldots, s_N}$ 上的离散时间马尔可夫链。

定义 2（状态转移矩阵）：

$$A = [a_{ij}]_{N \times N}, \quad a_{ij} = P(X_{t+1} = s_j | X_t = s_i)$$

转移矩阵 $A$ 是一个随机矩阵（Stochastic Matrix）：

$a_{ij} \geq 0$（非负性）
$\sum_{j=1}^{N} a_{ij} = 1$（行和为 1，从状态 $i$ 出发的概率分布归一化）

定义 3（初始状态概率分布）：

$$\pi = (\pi_1, \pi_2, \ldots, \pi_N), \quad \pi_i = P(X_1 = s_i)$$

$\sum_{i=1}^{N} \pi_i = 1$。

有了 $\pi$ 和 $A$，一个马尔可夫链就被完全定义了。任意时刻的状态分布可以通过迭代计算：

$$P(X_t) = \pi A^{t-1}$$

其中 $[A^{t-1}]_{ij}$ 表示从状态 $i$ 出发，经过 $t-1$ 步转移到状态 $j$ 的概率。

1.3 Chapman-Kolmogorov 方程

对于 $m$ 步转移概率：

$$a_{ij}^{(m)} = P(X_{t+m} = s_j | X_t = s_i)$$

C-K 方程给出了多步转移的递归关系：

$$a_{ij}^{(m+n)} = \sum_{k=1}^{N} a_{ik}^{(m)} a_{kj}^{(n)}$$

在矩阵形式下，$A^{(m)} = A^m$。这正是矩阵乘法的定义，证明了转移概率的矩阵表示是自洽的。

1.4 马尔可夫链的分类与性质

定义 4（不可约性 Irreducibility）：如果从任意状态 $i$ 出发，存在一条路径以正概率到达任意状态 $j$（即 $\exists m > 0, a_{ij}^{(m)} > 0$），则称该马尔可夫链是不可约的。直观理解：整个状态空间是”连通”的，没有”死胡同”。

定义 5（周期性 Periodicity）：状态 $i$ 的周期 $d(i) = \gcd{m \geq 1 : a_{ii}^{(m)} > 0}$。若 $d(i) = 1$，状态 $i$ 是非周期的。若所有状态都是非周期的，称链是非周期的。

定义 6（常返性 Recurrence）：状态 $i$ 是常返的，若从 $i$ 出发以概率 1 返回 $i$。在有限状态空间中，不可约链的所有状态都是正常返的。

定义 7（遍历性 Ergodicity）：一个不可约、非周期的马尔可夫链是遍历的。遍历性保证了唯一平稳分布的存在性和收敛性。

1.5 平稳分布与遍历定理

定义 4（平稳分布）：若存在概率分布 $\mu = (\mu_1, \ldots, \mu_N)$ 满足：

$$\mu = \mu A \quad \text{即} \quad \mu_j = \sum_{i=1}^{N} \mu_i a_{ij}$$

则 $\mu$ 是马尔可夫链的平稳分布（Stationary Distribution）。一旦到达平稳分布，链在任何后续时刻都保持此分布。

遍历定理：对于不可约（irreducible，从任何状态可达任何状态）、非周期（aperiodic）的有限状态马尔可夫链，平稳分布 $\mu$ 存在且唯一，且从任意初始分布出发，$t$ 步后的分布都收敛到 $\mu$。

平稳分布的计算：求解线性方程组 $\mu A = \mu, \sum \mu_i = 1$。这是一个 $N+1$ 个方程、$N$ 个未知数的超定系统，去掉第一个方程（行相关）中的一个即可求解。

数值示例：考虑一个两状态 Markov 链

$$A = \begin{bmatrix} 0.7 & 0.3 \\ 0.4 & 0.6 \end{bmatrix}$$

求解 $\mu A = \mu$：

$$\begin{cases} 0.7\mu_1 + 0.4\mu_2 = \mu_1 \\ 0.3\mu_1 + 0.6\mu_2 = \mu_2 \\ \mu_1 + \mu_2 = 1 \end{cases}$$

从第二式：$0.3\mu_1 = 0.4\mu_2 \Rightarrow \mu_1 = \frac{4}{3}\mu_2$

结合 $\mu_1 + \mu_2 = 1$：$\frac{4}{3}\mu_2 + \mu_2 = 1 \Rightarrow \mu_2 = \frac{3}{7}, \mu_1 = \frac{4}{7}$

验证 $\mu A = (\frac{4}{7}, \frac{3}{7}) \begin{bmatrix} 0.7 & 0.3 \ 0.4 & 0.6 \end{bmatrix} = (\frac{2.8+1.2}{7}, \frac{1.2+1.8}{7}) = (\frac{4}{7}, \frac{3}{7}) = \mu \checkmark$

平稳分布的极限解释：从任意初始分布 $\pi$ 出发，经过足够多步后，状态分布收敛到 $\mu$。例如：

$\pi^{(0)} = (1, 0)$（从状态 1 开始）：
$\pi^{(1)} = (0.7, 0.3)$, $\pi^{(2)} = (0.61, 0.39)$, $\pi^{(5)} \approx (0.5719, 0.4281)$, $\pi^{(10)} \approx (0.57143, 0.42857) \approx \mu$
$\pi^{(0)} = (0, 1)$（从状态 2 开始）：
$\pi^{(1)} = (0.4, 0.6)$, $\pi^{(5)} \approx (0.5713, 0.4287)$, $\pi^{(10)} \approx (0.57143, 0.42857) \approx \mu$

两种初始分布在约 10 步后都收敛到同一个平稳分布——这是遍历马尔可夫链的本质特征：记忆在迭代中衰减。

马尔可夫链的混合时间（Mixing Time）：从任意初始分布出发，$t$ 步内总变差距离 $d_{TV}(t) = \frac{1}{2}|\pi^{(t)} - \mu|_1 < \varepsilon$ 所需的最小步数。对于上面的矩阵，混合时间（$\varepsilon = 0.01$）约为 6 步。混合时间是衡量马尔可夫链”多快遗忘初始状态”的关键指标，在 MCMC 采样中有重要的实际意义。

2. 隐马尔可夫模型的定义

2.1 模型三要素

HMM 包含两层随机过程：

隐藏的马尔可夫链：状态序列 $I = (i_1, i_2, \ldots, i_T)$，$i_t \in {1, 2, \ldots, N}$。这是不可观测的。
观测序列：$O = (o_1, o_2, \ldots, o_T)$，$o_t \in {v_1, v_2, \ldots, v_M}$。这是我们能看到的。

HMM 由三个参数组完全确定：$\lambda = (\pi, A, B)$

参数	含义	维度	数学定义
$\pi$	初始状态概率向量	$1 \times N$	$\pi_i = P(i_1 = s_i)$
$A$	状态转移概率矩阵	$N \times N$	$a_{ij} = P(i_{t+1} = s_j \mid i_t = s_i)$
$B$	观测概率矩阵（发射概率）	$N \times M$	$b_j(k) = P(o_t = v_k \mid i_t = s_j)$

其中 $N$ 是状态数，$M$ 是观测符号数，$T$ 是序列长度。

2.2 HMM 的两个基本假设

假设 1（齐次一阶马尔可夫性）：当前状态只依赖于前一个状态。

$$P(i_t | i_{t-1}, i_{t-2}, \ldots, i_1) = P(i_t | i_{t-1})$$

假设 2（观测独立性）：当前观测只依赖于当前状态。

$$P(o_t | i_T, o_T, i_{T-1}, o_{T-1}, \ldots, i_1, o_1) = P(o_t | i_t)$$

这两个假设极大地简化了模型。虽然现实中它们可能不完全成立，但在大量实际应用中，HMM 仍然表现出色。当观测独立性假设不成立时，条件随机场（CRF）提供了更灵活的替代方案（见本系列下一篇）。

2.3 HMM 的三个经典问题

问题	名称	目标	算法
问题 1	概率计算（Evaluation）	给定 $\lambda$ 和 $O$，计算 $P(O	\lambda)$
问题 2	学习（Learning）	给定 $O$，估计最优 $\lambda$	Baum-Welch (EM)
问题 3	解码（Decoding）	给定 $\lambda$ 和 $O$，找最可能的 $I$	Viterbi 算法

这三个问题覆盖了 HMM 应用的全部场景，分别对应了概率推理、参数学习和序列标注。

3. 问题一：概率计算

3.1 直接计算（不可行）

给定 $\lambda = (\pi, A, B)$ 和观测序列 $O = (o_1, \ldots, o_T)$，计算 $P(O|\lambda)$。

对任意一个状态序列 $I = (i_1, \ldots, i_T)$，产生 $O$ 的联合概率为：

$$P(O, I | \lambda) = \pi_{i_1} b_{i_1}(o_1) \cdot a_{i_1 i_2} b_{i_2}(o_2) \cdots a_{i_{T-1} i_T} b_{i_T}(o_T)$$

$P(O|\lambda) = \sum_I P(O, I|\lambda)$，对所有可能的状态序列求和。有 $N^T$ 个可能的状态序列，每个序列的联合概率需要 $O(T)$ 计算，总计算量为 $O(T \cdot N^T)$——指数级，完全不可行（$N=5, T=100 \Rightarrow 5^{100}$ 种序列）。

3.2 前向算法（Forward Algorithm）

前向算法利用动态规划（DP）将复杂度从 $O(T N^T)$ 降到 $O(T N^2)$——这是 HMM 最关键的计算技巧之一。

定义 5（前向概率）：

$$\alpha_t(i) = P(o_1, o_2, \ldots, o_t, i_t = s_i | \lambda)$$

即截止到时刻 $t$，部分观测序列为 $o_1, \ldots, o_t$ 且时刻 $t$ 的状态为 $s_i$ 的概率。

递推关系：

初始化（$t = 1$）：

$$\alpha_1(i) = \pi_i \cdot b_i(o_1), \quad i = 1, \ldots, N$$

递推（$t = 1, 2, \ldots, T-1$）：

$$\alpha_{t+1}(j) = \left[ \sum_{i=1}^{N} \alpha_t(i) \cdot a_{ij} \right] \cdot b_j(o_{t+1}), \quad j = 1, \ldots, N$$

直觉解释：到达时刻 $t+1$ 状态 $j$ 并看到观测 $o_{t+1}$，等于从 $t$ 时刻的所有可能状态 $i$ 转移过来（$\alpha_t(i) \cdot a_{ij}$），再乘以发射 $o_{t+1}$ 的概率 $b_j(o_{t+1})$。

终止：

$$P(O | \lambda) = \sum_{i=1}^{N} \alpha_T(i)$$

前向算法伪代码：

Algorithm: Forward
Input: λ = (π, A, B), O = (o_1, ..., o_T)
Output: P(O|λ) 和前向概率 α

1. // 初始化
2. for i = 1 to N:
3.     α_1(i) = π_i · b_i(o_1)
4. // 递推
5. for t = 1 to T-1:
6.     for j = 1 to N:
7.         α_{t+1}(j) = [Σ_{i=1}^{N} α_t(i) · a_{ij}] · b_j(o_{t+1})
8. // 终止
9. P = Σ_{i=1}^{N} α_T(i)
10. return P, α

3.3 后向算法（Backward Algorithm）

前向算法从前往后递推。后向算法从后往前递推，为问题三（学习）做准备。

定义 6（后向概率）：

$$\beta_t(i) = P(o_{t+1}, o_{t+2}, \ldots, o_T | i_t = s_i, \lambda)$$

即给定时刻 $t$ 状态为 $s_i$ 的条件下，从 $t+1$ 到 $T$ 的部分观测序列的概率。

递推关系：

初始化（$t = T$）：

$$\beta_T(i) = 1, \quad i = 1, \ldots, N$$

设定为 1 是因为，在 $t = T$ 时”之后的观测”为空——空的概率为 1。

递推（$t = T-1, T-2, \ldots, 1$）：

$$\beta_t(i) = \sum_{j=1}^{N} a_{ij} \cdot b_j(o_{t+1}) \cdot \beta_{t+1}(j)$$

直觉解释：从时刻 $t$ 的状态 $i$ 出发，后面发生的事 = 跳到所有可能的下一状态 $j$（$a_{ij}$），观测到 $o_{t+1}$（$b_j(o_{t+1})$），然后从 $j$ 继续（$\beta_{t+1}(j)$）。

后向算法伪代码：

Algorithm: Backward
Input: λ = (π, A, B), O = (o_1, ..., o_T)
Output: 后向概率 β

1. // 初始化
2. for i = 1 to N:
3.     β_T(i) = 1
4. // 递推
5. for t = T-1 down to 1:
6.     for i = 1 to N:
7.         β_t(i) = Σ_{j=1}^{N} a_{ij} · b_j(o_{t+1}) · β_{t+1}(j)
8. return β

注意：前向和后向都可以单独计算出 $P(O|\lambda)$：

$$P(O|\lambda) = \sum_{i=1}^{N} \alpha_t(i) \beta_t(i), \quad \forall t \in \{1, \ldots, T\}$$

这个性质非常有用：可以在任意时刻 $t$ 用前后向概率的乘积和来计算总概率。这提供了多种校验和数值稳定的手段。

3.4 数值稳定性：Scaling 技巧

前向概率 $\alpha_t(i)$ 和后向概率 $\beta_t(i)$ 中包含了 $t$ 个概率的连乘积。当 $T$ 较大时（例如 $T > 100$），这些值会指数级趋近于 0（下溢），超出浮点数的表示范围。

解决方案：对前向概率做重归一化（Scaling）

引入尺度因子：

$$c_t = \frac{1}{\sum_{i=1}^{N} \hat{\alpha}_t(i)}$$

其中 $\hat{\alpha}_t(i)$ 是未归一化的前向概率。实际存储归一化后的值：

$$\tilde{\alpha}_t(i) = c_t \cdot \hat{\alpha}_t(i)$$

最终的对数似然可以通过累加尺度因子恢复：

$$\log P(O|\lambda) = -\sum_{t=1}^{T} \log c_t$$

这就是在概率空间中取对数来做计算（log-space）——将连乘变成累加，彻底避免下溢。

Log-Sum-Exp 技巧：在需要累加多个小数时，使用：

$$\log \sum_i \exp(x_i) = x_{\max} + \log \sum_i \exp(x_i - x_{\max})$$

其中 $x_{\max} = \max_i x_i$。先减去最大值再进行指数运算，防止 $\exp(x_i)$ 溢出。

4. 问题三：解码问题（Viterbi 算法）

4.1 最可能的状态序列

给定 $\lambda$ 和 $O$，求最可能的状态序列：

$$I^* = \arg\max_I P(I | O, \lambda)$$

等价于最大化联合概率（因为 $P(O|\lambda)$ 对不同的 $I$ 是常数）：

$$I^* = \arg\max_I P(I, O | \lambda)$$

4.2 Viterbi 算法的 DP 结构

Viterbi 算法使用与前向算法类似的 DP 框架，但将求和（$\sum$）替换为求最大（$\max$），并保存回溯路径。

定义 7（Viterbi 变量）：

$$\delta_t(i) = \max_{i_1, i_2, \ldots, i_{t-1}} P(i_1, i_2, \ldots, i_{t-1}, i_t = s_i, o_1, o_2, \ldots, o_t | \lambda)$$

即沿某条路径到达时刻 $t$ 状态 $i$，并观察到 $o_1, \ldots, o_t$ 的所有路径中的最大概率。

递推关系：

初始化：

$$\delta_1(i) = \pi_i \cdot b_i(o_1), \quad i = 1, \ldots, N$$
$$\psi_1(i) = 0, \quad i = 1, \ldots, N$$

$ \psi_t(i) $ 保存 $ t $ 时刻状态 $ i $ 对应的上一时刻最优状态。

递推（$t = 1, 2, \ldots, T-1$）：

$$\delta_{t+1}(j) = \left[ \max_{1 \leq i \leq N} \delta_t(i) \cdot a_{ij} \right] \cdot b_j(o_{t+1})$$

$$\psi_{t+1}(j) = \arg\max_{1 \leq i \leq N} \delta_t(i) \cdot a_{ij}$$

终止：

$$P^* = \max_{1 \leq i \leq N} \delta_T(i)$$
$$i_T^* = \arg\max_{1 \leq i \leq N} \delta_T(i)$$

最优路径回溯（$t = T-1, T-2, \ldots, 1$）：

$$i_t^* = \psi_{t+1}(i_{t+1}^*)$$

4.3 Viterbi 算法伪代码

Algorithm: Viterbi
Input: λ = (π, A, B), O = (o_1, ..., o_T)
Output: 最优状态序列 I* = (i_1*, ..., i_T*), 最优路径概率 P*

1. // 初始化
2. for i = 1 to N:
3.     δ_1(i) = π_i · b_i(o_1)
4.     ψ_1(i) = 0
5. // 递推
6. for t = 1 to T-1:
7.     for j = 1 to N:
8.         δ_{t+1}(j) = max_i [δ_t(i) · a_{ij}] · b_j(o_{t+1})
9.         ψ_{t+1}(j) = argmax_i [δ_t(i) · a_{ij}]
10. // 终止
11. P* = max_i δ_T(i)
12. i_T* = argmax_i δ_T(i)
13. // 回溯
14. for t = T-1 down to 1:
15.     i_t* = ψ_{t+1}(i_{t+1}*)
16. return I* = (i_1*, ..., i_T*), P*

计算复杂度：$O(T \cdot N^2)$（与前向算法相同）。

4.4 Viterbi 与前向算法的区别

两者框架相同，但一个关键操作不同：

	前向算法	Viterbi 算法
聚合操作	$\sum$	$\max$
目标	计算总概率	找最优路径
额外存储	无（仅 $\alpha$）	$\psi_t(j)$（回溯指针）
计算结果	$P(O\vert\lambda)$	$I^$, $P^$

为什么 Viterbi 找到的”最可能状态序列”与逐时刻取最可能状态不同？

逐时刻取 $\arg\max_i P(i_t = s_i | O, \lambda)$（即 $\arg\max_i \gamma_t(i)$）得到的是每个时刻独立看最可能的状态，这并不能保证整个序列在联合分布下是最优的——因为状态之间有转换概率约束，两个连续时刻各自最可能的状态之间的转换概率可能为 0。

例如，一个人可以”明天是晴天”和”后天是晴天”各自都很可能，但如果晴天到晴天之间有”必定下雨”的转换约束，那”晴→晴”的序列整体概率可能就是 0 了。Viterbi 正确地考虑了这种全局一致性。

4.5 Viterbi 的数值示例

假设一个简单的词性标注 HMM（两个隐藏状态：N 名词、V 动词；三个观测：fish, swim, eat）：

$$\pi = [0.6, 0.4]$$

$$A = \begin{bmatrix} P(N|N) & P(V|N) \\ P(N|V) & P(V|V) \end{bmatrix} = \begin{bmatrix} 0.7 & 0.3 \\ 0.4 & 0.6 \end{bmatrix}$$

$$B = \begin{bmatrix} P(\text{fish}|N) & P(\text{swim}|N) & P(\text{eat}|N) \\ P(\text{fish}|V) & P(\text{swim}|V) & P(\text{eat}|V) \end{bmatrix} = \begin{bmatrix} 0.4 & 0.1 & 0.5 \\ 0.3 & 0.6 & 0.1 \end{bmatrix}$$

观测序列 $O = (\text{fish}, \text{swim})$。求最可能的状态序列。

$t = 1$（初始化）：

$\delta_1(N) = \pi_N \cdot P(\text{fish}|N) = 0.6 \times 0.4 = 0.24$
$\delta_1(V) = \pi_V \cdot P(\text{fish}|V) = 0.4 \times 0.3 = 0.12$

$t = 2$（递推）：

$\delta_2(N) = \max[\delta_1(N) \cdot a_{NN}, \delta_1(V) \cdot a_{VN}] \cdot P(\text{swim}|N)$
$= \max[0.24 \times 0.7, 0.12 \times 0.4] \times 0.1 = \max[0.168, 0.048] \times 0.1 = 0.0168$
$\psi_2(N) = N$（回溯到 N）
$\delta_2(V) = \max[\delta_1(N) \cdot a_{NV}, \delta_1(V) \cdot a_{VV}] \cdot P(\text{swim}|V)$
$= \max[0.24 \times 0.3, 0.12 \times 0.6] \times 0.6 = \max[0.072, 0.072] \times 0.6 = 0.0432$
$\psi_2(V) = N$ 或 $V$（任选一个，都行）

终止：

$P^* = \max[\delta_2(N), \delta_2(V)] = 0.0432$，$i_2^* = V$

回溯：

$i_1^* = \psi_2(V) = N$（或 V，取决于实现细节）

最可能的标注为：fish/N, swim/V 或 fish/V, swim/V（取决于 $\psi_2(V)$ 的选择）。计算两个完整序列的概率：

$P(N, V, \text{fish}, \text{swim}) = 0.6 \times 0.4 \times 0.3 \times 0.6 = 0.0432$
$P(V, V, \text{fish}, \text{swim}) = 0.4 \times 0.3 \times 0.6 \times 0.6 = 0.0432$

两者一样！说明在这个小例子中两种标注方案等可能。实践中，当出现平票时算法的实现细节决定了取舍。

5. 问题二：学习问题（Baum-Welch 算法）

5.1 问题的形式化

给定观测序列 $O = (o_1, \ldots, o_T)$，但不知道状态序列。我们要估计模型参数 $\lambda = (\pi, A, B)$，使得 $P(O|\lambda)$ 最大化：

$$\lambda^* = \arg\max_\lambda P(O|\lambda)$$

这是一个典型的”带有隐变量的极大似然估计”问题。天然适合 EM 算法。

与监督学习不同（监督学习知道状态序列，可以直接统计频率），HMM 的学习是无监督的——我们只有观测序列，状态序列（隐变量）是我们需要推断的。

5.2 EM 算法的核心思想

EM（Expectation-Maximization）交替执行两步：

E 步（Expectation）：基于当前参数 $\lambda^{\text{old}}$，计算隐变量的后验期望
M 步（Maximization）：基于 E 步得到的期望，重新估计参数 $\lambda^{\text{new}}$

对于 HMM，EM 的实例化就是 Baum-Welch 算法。

5.3 两个关键概率

在 E 步中，我们需要计算两个关键概率：

定义 8（单个状态的后验概率 $\gamma$）：

$$\gamma_t(i) = P(i_t = s_i | O, \lambda)$$

利用前后向概率：

$$\gamma_t(i) = \frac{\alpha_t(i) \beta_t(i)}{P(O|\lambda)} = \frac{\alpha_t(i) \beta_t(i)}{\sum_{j=1}^{N} \alpha_t(j) \beta_t(j)}$$

定义 9（相邻状态转移的后验概率 $\xi$）：

$$\xi_t(i, j) = P(i_t = s_i, i_{t+1} = s_j | O, \lambda)$$

利用前后向概率：

$$\xi_t(i, j) = \frac{\alpha_t(i) a_{ij} b_j(o_{t+1}) \beta_{t+1}(j)}{P(O|\lambda)} = \frac{\alpha_t(i) a_{ij} b_j(o_{t+1}) \beta_{t+1}(j)}{\sum_{k=1}^{N} \sum_{l=1}^{N} \alpha_t(k) a_{kl} b_l(o_{t+1}) \beta_{t+1}(l)}$$

$\gamma$ 和 $\xi$ 的关系：

$$\gamma_t(i) = \sum_{j=1}^{N} \xi_t(i, j)$$

(当前在 $i$ 的概率 = 从 $i$ 出发到任意 $j$ 的概率之和)

5.4 M 步：参数重估计公式

基于 E 步计算出的 $\gamma$ 和 $\xi$，我们可以重新估计参数。这里的推导使用了计数思想：

$\pi_i$ 的重估计：初始状态为 $s_i$ 的期望次数

$$\hat{\pi}_i = \gamma_1(i)$$

$a_{ij}$ 的重估计：从状态 $s_i$ 转移到 $s_j$ 的期望次数 ÷ 从状态 $s_i$ 出发的期望次数

$$\hat{a}_{ij} = \frac{\sum_{t=1}^{T-1} \xi_t(i, j)}{\sum_{t=1}^{T-1} \gamma_t(i)}$$

$b_j(k)$ 的重估计：在状态 $s_j$ 下观测为 $v_k$ 的期望次数 ÷ 处于状态 $s_j$ 的期望次数

$$\hat{b}_j(k) = \frac{\sum_{t=1, o_t = v_k}^{T} \gamma_t(j)}{\sum_{t=1}^{T} \gamma_t(j)}$$

其中分子表示”当观测是 $v_k$ 时处于状态 $j$”的时刻的 $\gamma_t(j)$ 之和。

5.5 Baum-Welch 算法伪代码

Algorithm: Baum-Welch (HMM 的 EM)
Input: 观测序列 O = (o_1, ..., o_T), 收敛容差 ε
Output: λ = (π, A, B)

1. 随机初始化 λ = (π, A, B)  // 确保 π 和 A 的行和为 1
2. repeat:
3.     // ===== E-step =====
4.     使用 Forward(λ, O) 计算 α
5.     使用 Backward(λ, O) 计算 β
6.     计算 P(O|λ) = Σ_i α_T(i)  // 当前似然
7.     for t = 1 to T:
8.         计算 γ_t(i) = α_t(i)β_t(i) / P(O|λ)
9.     for t = 1 to T-1:
10.        计算 ξ_t(i,j) = α_t(i)a_{ij}b_j(o_{t+1})β_{t+1}(j) / P(O|λ)
11.    // ===== M-step =====
12.    π_i = γ_1(i)
13.    a_{ij} = Σ_{t=1}^{T-1} ξ_t(i,j) / Σ_{t=1}^{T-1} γ_t(i)
14.    b_j(k) = Σ_{t: o_t=v_k} γ_t(j) / Σ_{t=1}^{T} γ_t(j)
15. until |P(O|λ_new) - P(O|λ_old)| < ε
16. return λ

5.6 EM 算法的收敛性

EM 算法保证每次迭代后 $P(O|\lambda^{\text{new}}) \geq P(O|\lambda^{\text{old}})$，即似然单调不减。但 EM 只能保证收敛到局部最优而非全局最优——这意味着 Baum-Welch 的结果依赖于初始化。

初始化策略：

随机初始化多次，取最终似然最大的
使用先验知识（如先均匀分布，再逐步调整）
对于有标签数据，使用监督统计作为初始值，再进行无监督精炼
使用分段 K-means 初始化（Viterbi 训练的一种变体）

5.7 多序列学习

当有 $S$ 条独立的观测序列 ${O^{(1)}, O^{(2)}, \ldots, O^{(S)}}$ 时，似然是各序列似然之积：

$$P(\text{All observations} | \lambda) = \prod_{s=1}^{S} P(O^{(s)} | \lambda)$$

M 步的重估计公式需要聚合所有序列的统计量：

$$\hat{a}_{ij} = \frac{\sum_{s=1}^{S} \sum_{t=1}^{T_s-1} \xi_t^{(s)}(i, j)}{\sum_{s=1}^{S} \sum_{t=1}^{T_s-1} \gamma_t^{(s)}(i)}$$

$$\hat{b}_j(k) = \frac{\sum_{s=1}^{S} \sum_{t=1, o_t^{(s)} = v_k}^{T_s} \gamma_t^{(s)}(j)}{\sum_{s=1}^{S} \sum_{t=1}^{T_s} \gamma_t^{(s)}(j)}$$

6. 扩展 HMM

6.1 高斯 HMM（Gaussian HMM）

当观测是连续值而非离散符号时，$B$ 不再是离散概率分布矩阵，而是连续概率密度函数（通常使用高斯分布或高斯混合模型）。

对状态 $j$，假设观测服从多元高斯分布：

$$b_j(o) = \mathcal{N}(o; \mu_j, \Sigma_j) = \frac{1}{(2\pi)^{d/2} |\Sigma_j|^{1/2}} \exp\left(-\frac{1}{2}(o - \mu_j)^T \Sigma_j^{-1} (o - \mu_j)\right)$$

参数估计（M 步）变为：

$$\hat{\mu}_j = \frac{\sum_{t=1}^{T} \gamma_t(j) \cdot o_t}{\sum_{t=1}^{T} \gamma_t(j)}$$

$$\hat{\Sigma}_j = \frac{\sum_{t=1}^{T} \gamma_t(j) \cdot (o_t - \hat{\mu}_j)(o_t - \hat{\mu}_j)^T}{\sum_{t=1}^{T} \gamma_t(j)}$$

注意这里 $\mu_j$ 的更新公式是加权平均（权重为 $\gamma_t(j)$），这与 K-means 的软聚类更新非常类似。

Gaussian HMM 广泛应用于：

金融市场中的波动率状态建模
传感器信号处理
语音识别中的声学特征建模

6.2 高斯混合 HMM（GMM-HMM）

更进一步，每个状态可以用高斯混合模型（GMM）来建模更复杂的观测分布：

$$b_j(o) = \sum_{m=1}^{M} c_{jm} \mathcal{N}(o; \mu_{jm}, \Sigma_{jm})$$

其中 $c_{jm}$ 是状态 $j$ 下第 $m$ 个混合成分的权重，$\sum_m c_{jm} = 1$。

这引入了额外的隐变量（混合成分），形成”隐状态中的隐变量”结构。GMM-HMM 在传统的语音识别系统中曾是标准配置（在深度学习出现之前）。

6.3 自回归 HMM（AR-HMM）

当观测之间存在时间依赖时（违反观测独立性假设），可以使用自回归 HMM：

$$b_j(o_t | o_{t-1}, \ldots, o_{t-p}) = \mathcal{N}(o_t; \sum_{r=1}^{p} A_j^{(r)} o_{t-r}, \Sigma_j)$$

每个状态有自己的自回归系数矩阵 $A_j^{(r)}$。

6.4 输入-输出 HMM（IO-HMM）

IO-HMM 在转移概率和/或发射概率中引入外生输入变量 $u_t$（特征）：

$$a_{ij}(u_t) = \frac{\exp(w_{ij} \cdot u_t)}{\sum_k \exp(w_{ik} \cdot u_t)}$$

这使转移概率依赖于外部特征——一个直接的例子是：在序列标注任务中，特征可以包括当前词的上下文、词性等，而不只是 HMM 自身的状态转移。IO-HMM 可以看作是从 HMM 到 CRF 的过渡形式。

7. HMM 的统计性质与局限性

7.1 HMM 的指数族表示

HMM 的完整数据（状态 $I$ + 观测 $O$）的联合分布属于指数族。其对数似然可以写为参数与充分统计量的内积形式。这使得 M 步有闭式解（如前所示），也是 EM 算法在 HMM 上如此自然的深层原因。

7.2 HMM 的局限性：为什么需要 CRF

观测独立性假设过强：HMM 假设 $o_t$ 只依赖于 $i_t$。在序列标注任务中，这意味着相邻词的标注之间只能通过状态转移概率来交互，但不能直接利用观测之间的相关性（例如”形容词后跟名词”这一模式无法被直接建模）。
生成式 vs 判别式：HMM 是生成式模型（model $P(X,Y)$），而序列标注本质上是一个判别式任务（model $P(Y|X)$）。生成式模型需要为 $P(X)$ 建模——这是我们并不关心的。
无法使用重叠特征：HMM 中每个观测只能由当前状态”发射”。CRF 的势函数可以包含任意的、重叠的、手工设计的特征函数——这通常是序列标注任务中获得高精度的关键。

CRF（条件随机场）通过直接建模 $P(Y|X)$（判别式）并允许任意特征函数，解决了上述问题。但这并不意味着 HMM 被淘汰——HMM 作为生成式模型，在数据量较小、对模型的概率解释要求较强的场景中仍然有价值。且理解 HMM 是理解 CRF 的必要前提。

8. Baum-Welch 算法的深入分析

8.1 EM 算法的 Jensen 不等式推导

Baum-Welch 是 EM 算法的具体实例化。EM 的核心是构造对数似然的一个下界（Q 函数），然后迭代优化这个下界。

对于 HMM，对数似然为：

$$\log P(O|\lambda) = \log \sum_I P(O, I|\lambda) = \log \sum_I Q(I) \frac{P(O, I|\lambda)}{Q(I)}$$

由 Jensen 不等式（对 concave 函数 $\log$）：

$$\log \sum_I Q(I) \frac{P(O, I|\lambda)}{Q(I)} \geq \sum_I Q(I) \log \frac{P(O, I|\lambda)}{Q(I)} = \mathcal{F}(Q, \lambda)$$

这个下界 $\mathcal{F}(Q, \lambda)$ 称为变分自由能或 ELBO（Evidence Lower BOund）。

EM 算法的两个步骤可以等价地理解为：

E 步：固定 $\lambda^{\text{old}}$，最大化 $\mathcal{F}(Q, \lambda^{\text{old}})$ 关于 $Q$。最优解为 $Q(I) = P(I|O, \lambda^{\text{old}})$，即隐变量的后验分布。
M 步：固定 $Q$，最大化 $\mathcal{F}(Q, \lambda)$ 关于 $\lambda$。等价于最大化 $Q(\lambda, \lambda^{\text{old}}) = \sum_I P(I|O, \lambda^{\text{old}}) \log P(O, I|\lambda)$。

8.2 EM 收敛性的几何解释

在参数空间中，EM 算法可以理解为：在当前参数点处，构造一个在切线方向和曲率上匹配真实对数似然的辅助函数（下界），然后直接跳到该辅助函数的最高点。由于辅助函数是下界，跳到其最高点必然至少不低于当前点。与前向梯度下降不同，EM 的”步长”是由辅助函数的最优解自然决定的，不需要调学习率。

8.3 多个局部最优的应对

Baum-Welch 对初始化敏感。常见策略：随机初始化 $K=10-50$ 组参数，每组运行 Baum-Welch 至收敛，选最终似然最高的那组。对于 $N$ 较大的模型，可先用 $K$-means 聚类观测向量初始化发射概率，再用均匀分布初始化转移概率。

9. HMM 中的参数平滑

9.1 数据稀疏问题

在参数估计中，如果某个转移或发射事件在训练数据中从未出现，对应的概率被估计为 0。这导致测试时一旦遇到该事件，整个序列的概率变为 0——模型完全拒绝了一个”可能”但”未见过”的序列。

9.2 拉普拉斯平滑（Laplace Smoothing）

最简单的解决方案：给每个计数加一个小的正数：

$$a_{ij}^{\text{smoothed}} = \frac{C(i \to j) + \alpha}{\sum_k (C(i \to k) + \alpha)} = \frac{C(i \to j) + \alpha}{C(i) + \alpha N}$$

其中 $\alpha > 0$ 是平滑参数。$\alpha = 1$ 即为经典 Laplace 平滑。类似地发射概率：

$$b_j(k)^{\text{smoothed}} = \frac{C(j \to v_k) + \beta}{C(j) + \beta M}$$

9.3 其他平滑方法

Dirichlet 先验（贝叶斯平滑）：在 Baum-Welch 的 M 步中使用伪计数（pseudo-counts），等价于对 $A$ 和 $B$ 的行使用 Dirichlet 先验。

Jelinek-Mercer 平滑：$\hat{a}{ij} = \lambda a{ij}^{\text{ML}} + (1-\lambda) \frac{1}{N}$，其中 $\lambda \in [0,1]$ 控制对训练数据的信任度。

实践中推荐使用 Dirichlet 正则化（在 M 步分子分母加一个小的常数，如 0.01）。

10. HMM 与深度学习的关系

10.1 从 HMM 到 RNN/LSTM

传统 HMM 有两处限制：(1) 转移概率 $a_{ij}$ 是静态的（不依赖观测内容或更长历史）；(2) 发射概率 $b_j(o_t)$ 仅依赖当前状态。RNN/LSTM 在一个统一的神经网络架构中同时绕过这两个限制——隐藏状态通过可微递归计算，当前输出可依赖整个历史。从概率图视角看，HMM 是马尔可夫链（线性链 + 隐变量依赖），RNN 是确定性递归 + 非线性变换。

10.2 Bi-LSTM-CRF 的混合架构

这是序列标注的经典 SOTA 架构（在 Transformer 出现之前）：

Bi-LSTM 层：对每个词提取上下文相关的特征表示，捕获长程依赖
CRF 层：基于 LSTM 特征学习标注之间的转移约束

CRF 层确保标签之间的全局一致性（如 I-PER 前必须有 B-PER）——这与 HMM 的 Viterbi 解码非常相似，但发射得分来自 LSTM。Bi-LSTM-CRF 可视为 HMM 的判别式版本（CRF）+ 深度学习特征提取（LSTM）。

10.3 HMM 为什么未被淘汰

小样本场景：HMM 作为生成式模型（利用了 $P(X)$），在训练数据极少时可能优于判别式方法
可解释性要求：HMM 的状态有明确语义（如金融中”牛市/熊市”），深度网络隐藏状态难以解释
专家知识融入：HMM 的转移矩阵和发射矩阵可直接融入领域专家知识
计算资源受限：HMM 在嵌入式设备中 Viterbi 推理极快
生成任务：需要从模型中采样生成序列时（如蛋白质序列设计），HMM 更自然

11. 数值示例：完整的前向-后向-Viterbi-Baum-Welch 走查

11.1 问题设定

假设有一个简单的天气 HMM：

隐藏状态：Sunny（S，晴天）和 Rainy（R，雨天），$N=2$
观测：Walk（W，散步）、Shop（P，购物）、Clean（C，打扫），$M=3$

参数：

$$\pi = \begin{bmatrix} 0.6 & 0.4 \end{bmatrix}$$

$$A = \begin{bmatrix} 0.7 & 0.3 \\ 0.4 & 0.6 \end{bmatrix}$$

$$B = \begin{bmatrix} 0.6 & 0.3 & 0.1 \\ 0.1 & 0.4 & 0.5 \end{bmatrix} \quad \text{(行=S,R; 列=W,P,C)}$$

观测序列 $O = (W, P, C)$，$T=3$。

11.2 前向算法计算

**$t=1$**：

$\alpha_1(S) = \pi_S \cdot P(W|S) = 0.6 \times 0.6 = 0.36$
$\alpha_1(R) = \pi_R \cdot P(W|R) = 0.4 \times 0.1 = 0.04$

**$t=2$**：

$\alpha_2(S) = [\alpha_1(S) \cdot a_{SS} + \alpha_1(R) \cdot a_{RS}] \cdot P(P|S)$
$= [0.36 \times 0.7 + 0.04 \times 0.4] \times 0.3 = [0.252 + 0.016] \times 0.3 = 0.0804$
$\alpha_2(R) = [0.36 \times 0.3 + 0.04 \times 0.6] \times 0.4 = [0.108 + 0.024] \times 0.4 = 0.0528$

**$t=3$**：

$\alpha_3(S) = [0.0804 \times 0.7 + 0.0528 \times 0.4] \times 0.1 = [0.05628 + 0.02112] \times 0.1 = 0.00774$
$\alpha_3(R) = [0.0804 \times 0.3 + 0.0528 \times 0.6] \times 0.5 = [0.02412 + 0.03168] \times 0.5 = 0.0279$

$P(O|\lambda) = \alpha_3(S) + \alpha_3(R) = 0.00774 + 0.0279 = 0.03564$

11.3 后向算法计算

**$t=3$**：

$\beta_3(S) = 1, \beta_3(R) = 1$

**$t=2$**：

$\beta_2(S) = a_{SS} \cdot P(C|S) \cdot \beta_3(S) + a_{SR} \cdot P(C|R) \cdot \beta_3(R)$
$= 0.7 \times 0.1 \times 1 + 0.3 \times 0.5 \times 1 = 0.07 + 0.15 = 0.22$
$\beta_2(R) = 0.4 \times 0.1 \times 1 + 0.6 \times 0.5 \times 1 = 0.04 + 0.30 = 0.34$

**$t=1$**：

$\beta_1(S) = 0.7 \times 0.3 \times 0.22 + 0.3 \times 0.4 \times 0.34 = 0.0462 + 0.0408 = 0.087$
$\beta_1(R) = 0.4 \times 0.3 \times 0.22 + 0.6 \times 0.4 \times 0.34 = 0.0264 + 0.0816 = 0.108$

验证 $P(O|\lambda) = \pi_S \cdot P(W|S) \cdot \beta_1(S) + \pi_R \cdot P(W|R) \cdot \beta_1(R) = 0.6 \times 0.6 \times 0.087 + 0.4 \times 0.1 \times 0.108 = 0.03132 + 0.00432 = 0.03564 \checkmark$

11.4 $\gamma$ 和 $\xi$ 的计算

**$\gamma_t(i)$**：

以 $t=2$ 为例：

$\gamma_2(S) = \frac{\alpha_2(S) \cdot \beta_2(S)}{P(O)} = \frac{0.0804 \times 0.22}{0.03564} = \frac{0.017688}{0.03564} \approx 0.496$
$\gamma_2(R) = \frac{0.0528 \times 0.34}{0.03564} = \frac{0.017952}{0.03564} \approx 0.504$

解释：在 $t=2$ 看到 W,P,C 这些观测后，天气是晴天的后验概率约 49.6%，是雨天的后验概率约 50.4%——几乎等可能。

**$\xi_t(i,j)$**（以 $t=1, i=S, j=S$ 为例）：

$$\xi_1(S, S) = \frac{\alpha_1(S) \cdot a_{SS} \cdot P(P|S) \cdot \beta_2(S)}{P(O)} = \frac{0.36 \times 0.7 \times 0.3 \times 0.22}{0.03564} = \frac{0.016632}{0.03564} \approx 0.467$$

解释：给定三个观测，在 $t=1$ 是晴天、$t=2$ 还是晴天的后验概率约 46.7%。

11.5 Viterbi 解码

**$t=1$**：

$\delta_1(S) = 0.36, \psi_1(S) = 0$
$\delta_1(R) = 0.04, \psi_1(R) = 0$

**$t=2$**：

$\delta_2(S) = \max[0.36 \times 0.7, 0.04 \times 0.4] \times 0.3 = \max[0.252, 0.016] \times 0.3 = 0.0756$, $\psi_2(S)=S$
$\delta_2(R) = \max[0.36 \times 0.3, 0.04 \times 0.6] \times 0.4 = \max[0.108, 0.024] \times 0.4 = 0.0432$, $\psi_2(R)=S$

**$t=3$**：

$\delta_3(S) = \max[0.0756 \times 0.7, 0.0432 \times 0.4] \times 0.1 = \max[0.05292, 0.01728] \times 0.1 = 0.005292$, $\psi_3(S)=S$
$\delta_3(R) = \max[0.0756 \times 0.3, 0.0432 \times 0.6] \times 0.5 = \max[0.02268, 0.02592] \times 0.5 = 0.01296$, $\psi_3(R)=R$

回溯：$i_3^* = R$, $i_2^* = \psi_3(R) = R$, $i_1^* = \psi_2(R) = S$

最优状态序列：(Sunny, Rainy, Rainy) — 对应观测 (Walk, Shop, Clean)。

解释：第一天晴天去散步，然后变雨天，继续购物和打扫。这与直觉吻合——雨天更可能在家购物和打扫，晴天更可能外出散步。但第二天开始下雨后转为室内活动。

11.6 对比：逐时刻最优 vs Viterbi 全局最优

我们计算出 $\gamma_t(i)$ 以对每一时刻独立选最优状态：

$\gamma_1(S) = \frac{\alpha_1(S)\beta_1(S)}{P(O)} = \frac{0.36 \times 0.087}{0.03564} \approx 0.878$, $\gamma_1(R) \approx 0.122$ → 选 S
$\gamma_2(S) \approx 0.496$, $\gamma_2(R) \approx 0.504$ → 选 R
$\gamma_3(S) = \frac{0.00774 \times 1}{0.03564} \approx 0.217$, $\gamma_3(R) \approx 0.783$ → 选 R

逐时刻最优序列也是 (S, R, R)，与 Viterbi 一致——说明在这个例子中两种方法恰巧一致。但这个一致性不总是成立的。

反例：假设修改转移矩阵使 $a_{SR} = 0$（从晴天不能直接变雨天）但 $a_{RR} = 1.0$（雨天永远是雨天）。如果 $\gamma_1$ 选择 S 但 $\gamma_2$ 选择 R，Viterbi 会拒绝这个序列（因为 $a_{SR} = 0$ 使该路径概率为 0），转而选择次优但全局相容的序列。而逐时刻贪心法则会输出一个概率为 0 的不合法序列。这就是 Viterbi 的全局最优与逐时刻最优的核心分歧。

12. HMM 扩展应用详解

12.1 中文分词中的 HMM

中文分词是序列标注问题（BEMS 标注：Begin, End, Middle, Single）。观测是字符序列，隐藏状态是 BEMS 标签。

观测序列：”我爱北京天安门” → “我/爱/北京/天安门”

HMM 分词流程：

训练阶段：在有标注的分词语料上统计初始状态概率、状态转移概率（如 P(E|B) 应该很低）、发射概率（如 P(‘京’|E) 应该较高）
解码阶段：对任意句子，运行 Viterbi 算法找到最优的 BEMS 标签序列
后处理：根据 BEMS 标签切词

现代中文分词系统（如 jieba）在 HMM 基础上增加词典匹配和规则来提升准确率。

12.2 语音识别中的 GMM-HMM

在深度学习统治 ASR 之前，GMM-HMM 是语音识别的标准框架：

HMM 的状态对应音素（phone）的子段（通常一个音素 3-5 个状态）
GMM 建模每个状态下声学特征（MFCC 等）的分布
语言模型（N-gram）提供状态转移的先验

整个系统是 HMM（时序结构）+ GMM（声学建模）+ N-gram（语言约束）的组合。Viterbi 解码在音素、词、句子的三级网络中搜索最优路径。

12.3 生物信息学中的 Profile HMM

Profile HMM 用于蛋白质/基因家族的序列比对。它在标准 HMM 基础上引入了三个状态类型：

Match 状态：匹配位置（发射氨基酸）
Insert 状态：插入片段（额外序列）
Delete 状态：缺失（静默状态，不发射）

这使得它可以建模序列进化中的 indel（插入/缺失）事件。HMMER 是使用 Profile HMM 的流行工具，在 Pfam、InterPro 等蛋白质家族数据库中发挥关键作用。

13. 面试高频问答

Q1: 前向算法和后向算法的关系是什么？为什么两者都能计算 $P(O|\lambda)$？

A: 前向概率 $\alpha_t(i) = P(o_1, \ldots, o_t, i_t = s_i | \lambda)$ 涵盖”过去到现在的观测和现在状态”。后向概率 $\beta_t(i) = P(o_{t+1}, \ldots, o_T | i_t = s_i, \lambda)$ 涵盖”现在状态给定下未来观测”。两者之积 $\alpha_t(i) \beta_t(i) = P(O, i_t = s_i | \lambda)$，求和即得 $P(O|\lambda)$。这个等式对任意 $t$ 都成立。

Q2: Viterbi 算法中的 $\max$ 操作可以换成 $\sum$ 吗？如果可以，结果是什么？

A: 如果把 Viterbi 算法中的 $\max$ 全部换成 $\sum$，得到的恰好就是前向算法。两者使用完全相同的 DP 网格结构，只是聚合算子不同：$\max$ 找最优单条路径概率，$\sum$ 累积所有路径的总概率。Viterbi 是前向算法在 max-plus 半环上的版本。

Q3: Baum-Welch 算法的 E 步为什么要计算 $\xi_t(i,j)$？

A: $\xi_t(i,j) = P(i_t = s_i, i_{t+1} = s_j | O, \lambda)$ 是给定观测后相邻状态转移的后验概率。M 步重估计转移概率 $a_{ij}$ 需要”从 $s_i$ 到 $s_j$ 的期望次数”（$\sum_t \xi_t(i,j)$）除以”处于 $s_i$ 的期望次数”（$\sum_t \gamma_t(i)$）。这是 EM 中”用后验期望替代缺失数据”的具体体现。

Q4: HMM 为什么使用 EM 算法？为什么不能直接解析求解？

A: HMM 的似然函数 $P(O|\lambda) = \sum_I P(O, I|\lambda)$ 是 log-sum 形式：$L(\lambda) = \log \sum_I P(O, I|\lambda)$，参数耦合在一起无法解析求导。EM 算法通过 Jensen 不等式将 $\log \sum$ 交换为 $\sum \log$，构造下界（Q 函数）来间接优化，M 步有闭式解。

Q5: 前向算法的数值下溢问题如何解决？

A: 前向概率 $\alpha_t(i)$ 随着 $T$ 增大呈指数级趋于 0。解决方案是 Scaling（重归一化）：每一步将 $\alpha_t(i)$ 归一化为概率分布，记录尺度因子 $c_t$。最终 $\log P(O|\lambda) = -\sum_t \log c_t$，在 log-space 中完成所有计算。Baum-Welch 中 $\gamma_t(i)$ 和 $\xi_t(i,j)$ 的分子分母 scale factor 可互相抵消。

Q6: 为什么 HMM 是生成式模型，而 CRF 是判别式模型？

A: HMM 对联合分布 $P(O, I)$ 建模（生成式）：定义状态序列如何生成（$\pi, A$），以及每个状态下观测如何生成（$B$），可以从模型中采样完整的数据对。CRF 直接对条件分布 $P(I|O)$ 建模（判别式）：不建模 $P(O)$，可自由使用 $O$ 中的任意特征（全局、重叠特征），无需担心如何生成观测。生成式模型在小样本下可能更好（建模 $P(X)$ 是一种正则化）；判别式模型在有足够标注数据时通常精度更高。

Q7: Viterbi 算法找到的”逐时刻最可能状态”和”全局最可能状态序列”有什么不同？为什么前者可能不合法？

A: 逐时刻取 $\arg\max_i \gamma_t(i)$（即给定观测下每个时刻独立最优的状态）不一定构成联合概率最大的完整序列。原因：两个连续时刻各自最可能状态之间的转移概率 $a_{ij}$ 可能为 0 或极低。例如 $\gamma_1=(0.6, 0.4)$ 选择状态 1，$\gamma_2=(0.5, 0.5)$ 选择状态 2，但如果 $a_{12}=0$（从状态 1 无法直接跳到状态 2），这个序列的联合概率就是 0。Viterbi 正确地考虑了状态间的全局相容性，找到的是联合最可能的完整序列。

Q8: 如何处理观测序列长度不同的问题？

A: HMM 天然支持变长序列——参数 $\pi, A, B$ 与序列长度 $T$ 无关。训练时对每条序列独立做前向-后向（用其自身的 $T_s$），将所有序列的 $\gamma_t^{(s)}$ 和 $\xi_t^{(s)}$ 聚合计入 M 步。这也适用于多个独立序列的 Baum-Welch 训练。

Q9: 描述 Baum-Welch 算法中 E 步和 M 步的具体计算过程。

E 步：使用当前参数 $\lambda^{\text{old}}$ 运行前向-后向算法，计算 $\alpha_t(i), \beta_t(i)$ 以及 $P(O|\lambda^{\text{old}})$。然后计算后验概率 $\gamma_t(i) = \alpha_t(i)\beta_t(i)/P(O|\lambda)$ 和 $\xi_t(i,j) = \alpha_t(i)a_{ij}b_j(o_{t+1})\beta_{t+1}(j)/P(O|\lambda)$。
M 步：用 E 步的期望统计量重新估计参数。$\hat{\pi}i = \gamma_1(i)$，$\hat{a}{ij} = \sum_t \xi_t(i,j) / \sum_t \gamma_t(i)$，$\hat{b}j(k) = \sum{t: o_t=v_k} \gamma_t(j) / \sum_t \gamma_t(j)$。这些公式的直觉是”期望计数除以期望总数”，是完整数据 MLE 的”软”版本（用后验概率代替了硬分配）。

Q10: 如果 HMM 中有 $N=100$ 个隐藏状态，观测序列 $T=1000$，运行一次 Viterbi 解码需要多少次基本运算？

A: Viterbi 的主要计算量在递推步骤。对每个时刻 $t = 1, \ldots, T-1$，对每个当前状态 $j$ 计算 $\delta_{t+1}(j) = \max_i [\delta_t(i) \cdot a_{ij}] \cdot b_j(o_{t+1})$，这需要遍历所有 $N$ 个前一状态 $i$ 来求 max。因此每步需要 $O(N^2)$ 次乘法+比较。总计 $(T-1) \times N^2 \approx 1000 \times 10000 = 10^7$ 次基本运算。加上回溯阶段需要 $T-1$ 次查表，总计约 10M 次运算——在现代 CPU 上不到 1 毫秒即可完成。这体现了 DP 的强大：指数级搜索空间被压缩到了二次复杂度。

14. 总结

HMM 是时序概率建模的基石。三个经典问题与三个核心算法构成了一个完整的理论体系：

问题	算法	类型	复杂度
概率计算	前向-后向	DP / 精确	$O(T N^2)$
解码	Viterbi	DP / 精确	$O(T N^2)$
学习	Baum-Welch (EM)	迭代 / 局部最优	$O(ITER \cdot T N^2)$

HMM 的三个关键假设——马尔可夫性、齐次性、观测独立性——使计算变得可行，但也限制了模型的表达能力。理解这些假设何时合理、何时需要被松弛，是选择 HMM vs CRF vs RNN/LSTM/Transformer 的判断依据。

HMM 的思想遗产远不止模型本身：

动态规划（前向、Viterbi）在序列建模中的广泛应用
EM 算法的无监督学习范式
隐变量模型的概率推断框架
从 HMM 到 CRF、从 HMM 到 HMM+LSTM 的发展脉络

都深刻影响了后续的统计学习和深度学习模型。

HMM 三大算法复杂度汇总：

算法	问题	时间复杂度	空间复杂度	递归方向
Forward	概率计算	$O(T N^2)$	$O(TN)$ 或优化为 $O(N)$	前向递推
Backward	概率计算	$O(T N^2)$	$O(TN)$ 或优化为 $O(N)$	后向递推
Viterbi	解码	$O(T N^2)$	$O(TN)$（存回溯指针）	前向递推 + 回溯
Baum-Welch	学习	$O(\text{iter} \cdot T N^2)$	$O(TN)$（存储 $\alpha, \beta, \gamma, \xi$）	前向+后向 + 参数重估

HMM 与相关模型的对比：

模型	类型	观测独立性	特征灵活性	训练	状态空间
HMM	生成式	是	低（预定义发射概率）	EM (Baum-Welch)	离散
CRF	判别式	否	高（任意特征）	梯度下降/L-BFGS	离散
HMM+SVM	混合	是	中	分段训练	离散
MEMM	判别式	否	高	MLE	离散（有标签偏置）
RNN/LSTM	判别式	否	极高	BPTT / SGD	连续（隐藏向量）
Transformer	判别式	否	极高	Adam/SGD	连续（自注意力）

在表格数据或小样本序列建模中，HMM 和 CRF 仍然是最稳健的选择；在大规模 NLP 任务中，Transformer 已经取代了 LSTM-CRF 的 SOTA 地位。

扩展阅读推荐：

Rabiner (1989), A Tutorial on Hidden Markov Models and Selected Applications in Speech Recognition — HMM 的经典教程，引用量超过 30000
Bishop (2006), Pattern Recognition and Machine Learning, Chapter 13 — HMM 的 PRML 视角
Murphy (2012), Machine Learning: A Probabilistic Perspective, Chapter 17 — HMM 在现代概率图框架下的全面讲解
L. R. Rabiner & B. H. Juang (1986), An Introduction to Hidden Markov Models — 另一篇经典综述

HMM 实践工具推荐：

hmmlearn（Python）：pip install hmmlearn，支持 Gaussian HMM、GMM-HMM、Multinomial HMM，API 设计类似 sklearn
pomegranate（Python）：更通用的概率模型库，支持 HMM、贝叶斯网络、混合模型
seqHMM（R）：用于序列分析和纵向数据的 HMM 库
GHMM（C）：GNU 的通用 HMM 工具箱，支持自定义发射分布

Q11: HMM 中的马尔可夫假设（一阶）如果被违反了会怎样？如何检测和缓解？

A: 如果真实数据不满足一阶马尔可夫性（即 $P(i_t | i_{t-1}, i_{t-2}) \neq P(i_t | i_{t-1})$），模型仍有偏——它会把高阶依赖”压缩”到一阶转移概率中，但这会导致模型不够精确。检测方法：(1) 在验证集上检查残差是否仍有自相关性；(2) 比较一阶和二阶 HMM 的似然（似然比检验）；(3) 在训练数据上统计二阶转移频率并与一阶模型的预测对比。缓解方法：(1) 使用 $m$ 阶 HMM（代价是参数从 $O(N^2)$ 增加到 $O(N^{m+1})$）；(2) 使用分段平稳假设，将时间序列切分为多个近似平稳的段；(3) 使用 RNN/LSTM 替代 HMM 的离散状态转移建模。

Q12: 多序列 Baum-Welch 训练中，如果各序列的质量差异很大怎么办？

A: 可以通过为每条序列分配权重 $w^{(s)}$ 来反映其置信度或重要性。在 M 步中，将 $\gamma_t^{(s)}(i)$ 替换为 $w^{(s)} \cdot \gamma_t^{(s)}(i)$（加权期望计数）。这等价于在似然函数中使用加权对数似然 $\sum_s w^{(s)} \log P(O^{(s)}|\lambda)$。常见加权策略：(1) 按序列长度加权（长序列可能更可靠）；(2) 按序列的标注质量/置信度加权（如果部分序列有弱标注或自动标注）；(3) 对于异常序列（如含大量未见观测的序列）降低权重。

Q13: 在深度学习时代，为什么还要学习 HMM？

概率图模型的思维范式：理解 HMM 的推断（前向-后向）、解码（Viterbi）和学习（EM）三件套，是理解更复杂概率模型（CRF, LDA, HDP-HMM, 状态空间模型, 深度生成模型 VAE/扩散模型）的基础。概率推断的动态规划范式（sum-product 和 max-product）贯穿了所有图模型。
小数据 + 强先验场景：在很多科学领域（如生物信息学、金融计量、语音学），可用的标注数据远少于 NLP，但领域知识可以自然编码为 HMM 的参数结构。此时 HMM 显著优于深度网络。
解释性要求：在医疗诊断、金融监管等场景，模型分析需要一个可解释的”状态”概念。HMM 的隐藏状态天然对应可解释的系统模式。
算法思想的迁移：Viterbi 算法拓展为 Beam Search（现代 NLP 解码的标配），前向-后向算法拓展为 seq2seq 的注意力权重计算（本质上是软对齐），EM 思想深刻影响了 VAE 和扩散模型的训练。

本系列下一篇文章（也是本系列最后一篇）将深入探讨条件随机场（CRF）——从无向图模型到线性链 CRF，完成从 HMM 到 CRF 的认知跨越，敬请期待。

HMM 关键公式速查：

公式	表达式	用途
前向概率	$\alpha_{t+1}(j) = [\sum_i \alpha_t(i) a_{ij}] b_j(o_{t+1})$	概率计算
后向概率	$\beta_t(i) = \sum_j a_{ij} b_j(o_{t+1}) \beta_{t+1}(j)$	概率计算
似然	$P(O\vert\lambda) = \sum_i \alpha_T(i) = \sum_i \pi_i b_i(o_1) \beta_1(i)$	概率计算
$\gamma$	$\gamma_t(i) = \frac{\alpha_t(i) \beta_t(i)}{P(O\vert\lambda)}$	E 步
$\xi$	$\xi_t(i,j) = \frac{\alpha_t(i) a_{ij} b_j(o_{t+1}) \beta_{t+1}(j)}{P(O\vert\lambda)}$	E 步
Viterbi 递推	$\delta_{t+1}(j) = [\max_i \delta_t(i) a_{ij}] b_j(o_{t+1})$	解码
M 步-$A$	$\hat{a}_{ij} = \frac{\sum_t \xi_t(i,j)}{\sum_t \gamma_t(i)}$	学习
M 步-$B$	$\hat{b}j(k) = \frac{\sum{t: o_t=v_k} \gamma_t(j)}{\sum_t \gamma_t(j)}$	学习
Scaling	$\log P(O\vert\lambda) = -\sum_t \log c_t$	数值稳定

HMM 学习资源推荐：

《统计自然语言处理》（宗成庆）：HMM 在中文分词和词性标注中的应用
《Speech and Language Processing》（Jurafsky & Martin, 3rd ed.）：NLP 视角的 HMM 教程
《Pattern Recognition and Machine Learning》（Bishop）第 13 章：概率图模型框架下的 HMM
《Biological Sequence Analysis》（Durbin et al.）：HMM 在生物信息学中的经典应用
hmmlearn 文档：Python 中 Gaussian HMM 和 Multinomial HMM 的最佳入门

最后的思考：HMM 的价值不在于它今天的 SOTA 地位（它早已不是），而在于它给序列建模提供了一套完整的思维工具——状态空间建模、动态规划推断、EM 参数学习。这套工具在后来的 CRF、LSTM、seq2seq、Transformer 中一再重现。理解了 HMM 的三问题三算法，你就掌握了序列建模的”语法”——具体用哪种”词汇”（HMM/CRF/LSTM/Transformer）去表达，都只是选择而已。

HMM 是统计学习领域少有的”数学优雅 + 计算可行 + 应用广泛”三位一体的模型。从 1960 年代 Baum 等人的开创性工作，到 1980 年代 Rabiner 将其系统化并应用于语音识别，再到今天在生物信息学和金融时序中的持续使用——HMM 的半个多世纪的历史证明了它的基础性地位。在深度学习时代重读 HMM，你会发现动态规划、消息传递、对数空间计算、EM 变分下界这些概念，在注意力和 Transformer 中一次次重现——只是换了个名字。

这也是本系列的倒数第二篇。下一篇——条件随机场——将把 HMM 的”局部归一化概率空间”升级为 CRF 的”全局归一化得分空间”，完成从生成式到判别式的序列建模范式迁移。

感谢阅读。 如果对你有所启发，欢迎分享和讨论。

统计学习不是记忆公式，是理解每个公式背后的优化目标和概率假设。
HMM 是这种理解的绝佳起点——三个问题、三个算法、一个统一的 DP 框架。