HMMs (Forward Algorithm, Viterbi Algorithm)(15)

上一讲我们学习了如何对一个状态完全可观察的、随时间演化的系统进行建模和预测。然而，在现实世界中，我们往往无法直接观察到系统的真实状态，只能通过一些间接的、带有噪声的观测来进行推断。

这一讲的核心就是解决这个问题：当系统的真实状态是“隐藏”的时，我们如何通过一系列的观测证据来推断这个隐藏状态？这个强大的工具就是隐马尔可夫模型 (Hidden Markov Models - HMMs)。

马尔可夫模型 (Markov Model)：我们能直接看到每天的天气（状态）。隐马尔可夫模型 (Hidden Markov Model)：我们看不到真实的天气（状态是隐藏的），但我们每天能看到天气预报（证据/观测）。我们需要根据天气预报来推断真实的天气。

隐马尔可夫模型 (Hidden Markov Models)

HMM 在马尔可夫模型的基础上，增加了一层可观测的证据变量 (evidence variable)****

状态变量 (State Variable)：W_i，代表在时间 i 的隐藏状态（例如，真实天气）。
证据变量 (Evidence Variable)：F_i，代表在时间 i 的可观测证据（例如，天气预报）。

因此, 一个平稳的 HMM 可以由三个部分紧凑地表示：

初始分布：P(W₀)，和马尔可夫模型一样。
转移模型：P(W_i + 1|W_i)，和马尔可夫模型一样。
传感器模型 (Sensor Model)：P(F_i|W_i)，这是新增的部分。它描述了在给定真实状态 Wᵢ 的情况下，观测到证据 Fᵢ 的概率。

另一方面, 对于HMM模型, 由于我们往往无法直接得到隐藏状态, 所以更关注两个概念:

置信度 (Belief) B(Wₜ): 在观测到截至时刻 t 的所有证据后，对状态 Wₜ 的概率分布。即, B(Wₜ) = P(Wₜ|f₁,…,fₜ)
预测置信度 (Predicted Belief) B’(Wₜ): 在观测到截至时刻 t-1 的证据后，对状态 Wₜ 的概率分布。即, B’(Wₜ) = P(Wₜ|f₁,…,fₜ₋₁)

更多时候, 我们使用 λ = (A, B, π) 代表整个HMM, 其中: - A：状态转移概率矩阵, 表示从一个隐藏状态转移到另一个隐藏状态的概率。 - 如果有 N 个隐藏状态，A 是一个 N×N 的矩阵，其中 A[i][j] 表示从状态 i 转移到状态 j 的概率。 - B：观测概率矩阵（或称发射概率矩阵）: 表示在某个隐藏状态下，生成某个观测值的概率。 - 如果有 M 个可能的观测值，B 是一个 N×M 的矩阵，其中 B[i][k] 表示在隐藏状态 i 下观察到观测值 k 的概率。 - π：初始状态概率向量: 表示系统在时间 t=1 时处于每个隐藏状态的概率。 - π[i] 是系统一开始处于状态 i 的概率。

前向算法

在隐马尔可夫模型中，我们拥有一系列随时间变化的、不可直接观测的隐藏状态（例如，每天的真实天气 W_t），以及在每个时间点上可以观测到的证据（例如，天气预报 F_t）。

前向算法的核心目标是计算在观测到所有历史证据后，当前隐藏状态的概率分布。这个分布就是我们前面提到的置信度（Belief）。换句话说，算法要计算的是 P(W_t|f₁, f₂, ..., f_t)，即在已知从第1天到第t天的所有天气预报的情况下，第t天真实天气的概率分布。

前向算法通过一个迭代过程，从 B(W_i) 计算出 B(W_i + 1)。这个过程可以分解为两个关键步骤：时间流逝更新(Time Elapse Update) 和 观测更新 (Observation Update)

第一步：时间流逝更新 (Time Elapse Update)

目标：从当前时刻 i 的置信度 B(W_i) 推导出下一时刻 i+1 的预测置信度 B^′(W_i + 1)。这一步相当于在没有看到新证据之前，对未来状态进行预测。

推导过程如下：

从定义出发：我们想要求解 B^′(W_i + 1) = P(W_i + 1|f₁, ..., f_i)。

引入上一时刻的状态 W_i：为了将 B(W_i) 引入计算，我们使用边缘化（Marginalization），对 W_i 的所有可能取值 w_i 进行求和: B^′(W_i + 1) = ∑_{w_i}P(W_i + 1, w_i|f₁, ..., f_i)
这个公式的意义是，“在已知历史证据的情况下，下一状态为 W_i + 1 的概率”等于”在同样证据下，当前状态为某个具体值 w_i 且下一状态为 W_i + 1 的所有可能情况的概率之和”。

应用条件概率的链式法则：将联合概率 P(W_i + 1, w_i|...) 分解: B^′(W_i + 1) = ∑_{w_i}P(W_i + 1|w_i, f₁, ..., f_i)P(w_i|f₁, ..., f_i)

这是条件概率的基本性质，将联合概率分解为两个条件的乘积。
利用HMM的条件独立性假设：HMM假设未来的状态 W_i + 1 仅依赖于当前状态 W_i，而与过去的所有证据 f₁, ..., f_i 无关。这被称为马尔可夫性。
因此，P(W_i + 1|w_i, f₁, ..., f_i) 可以简化为 P(W_i + 1|w_i)，这就是转移模型（Transition Model）。

识别置信度：同时，我们发现 P(w_i|f₁, ..., f_i) 正是当前置信度 B(w_i) 的定义。

得出最终公式：将简化后的项代入，得到时间更新的最终形式: B^′(W_i + 1) = ∑_{w_i}P(W_i + 1|w_i)B(w_i)
这个公式表明，对下一状态的预测是所有可能当前状态的置信度 B(w_i)，按照它们各自转移到下一状态的概率 P(W_i + 1|w_i) 进行加权求和的结果。

第二步：观测更新 (Observation Update)

目标：在获得 i+1 时刻的新证据 f_i + 1 后，将预测置信度 B^′(W_i + 1) 更新为真实置信度 B(W_i + 1)。

推导过程如下：

从定义出发：我们想要求解 B(W_i + 1) = P(W_i + 1|f₁, ..., f_i, f_i + 1)。

应用贝叶斯法则（或条件概率定义）：
我们将求解目标中的条件 f_i + 1 移到分子中，形成一个联合概率。

引入归一化技巧：分母 P(f_i + 1|f₁, ..., f_i) 是一个常数，因为它不随 W_i + 1 的变化而改变。为了简化计算，我们可以暂时忽略它，用比例关系 ∝ 来表示: B(W_i + 1) ∝ P(W_i + 1, f_i + 1|f₁, ..., f_i)

这意味着 B(W_i + 1) 的分布形状由分子决定，我们可以在所有计算完成后再进行归一化，使其概率和为1。
接着再次应用链式法则：对分子进行分解: B(W_i + 1) ∝ P(f_i + 1|W_i + 1, f₁, ..., f_i)P(W_i + 1|f₁, ..., f_i)

利用HMM的条件独立性假设：HMM假设当前证据 f_i + 1 仅依赖于当前状态 W_i + 1，而与过去的状态和证据无关。

因此，P(f_i + 1|W_i + 1, f₁, ..., f_i) 可以简化为 P(f_i + 1|W_i + 1)，这被称为传感器模型（Sensor Model）。

识别预测置信度：同时，我们发现 P(W_i + 1|f₁, ..., f_i) 正是我们在第一步中计算出的预测置信度 B^′(W_i + 1)。

得出最终公式：将简化后的项代入，得到观测更新的最终形式, 即: B(W_i + 1) ∝ P(f_i + 1|W_i + 1)B^′(W_i + 1)
这个公式的意义是，结合新证据后的置信度，与”预测置信度”和”新证据在该状态下的出现概率”的乘积成正比。如果观测到的证据在某个状态下很可能发生，那么这个状态的置信度就会相应提高。

前向算法的核心

通过将上述两个步骤结合起来，我们就得到了一个从 B(W_i) 计算 B(W_i + 1) 的完整迭代公式：

B(W_i + 1) ∝ P(f_i + 1|W_i + 1)B^′(W_i + 1)

B(W_i + 1) ∝ P(f_i + 1|W_i + 1)∑_{w_i}P(W_i + 1|w_i)B(w_i)

这个公式就是前向算法的核心。从初始置信度 B(W₀) 开始，我们可以通过这个公式逐步递推，计算出任何时刻 t 的置信度分布 B(W_t)。在计算出每个 B(W_i + 1) 的非归一化值后，只需将所有可能状态的概率值相加，再用每个值除以这个总和，即可完成归一化，得到最终的概率分布。

示例分析

我们首先需要定义模型的三个核心部分：初始分布、转移模型和传感器模型。

初始置信度 B(W₀)表示第0天的天气概率分布。

W₀	B(W₀)
sun	0.8
rain	0.2

转移模型 P(W_i + 1|W_i)表示从今天的天气转移到明天天气的概率。

W_i	W_i + 1	P(W_i + 1\|W_i)
sun	sun	0.6
sun	rain	0.4
rain	sun	0.1
rain	rain	0.9

传感器模型 P(F_i|W_i)表示在某种真实天气下，观测到特定天气预报的概率。

W_i	F_i	P(F_i\|W_i)
sun	good	0.8
sun	bad	0.2
rain	good	0.3
rain	bad	0.7

我们的目标是计算第1天的置信度 B(W₁)，假设我们观测到第1天的天气预报是”good” (f₁ = good)。

整个计算过程遵循前向算法的两个核心步骤：时间流逝更新和观测更新。

第一步是时间流逝更新（预测）: 根据第0天的置信度 B(W₀)，预测第1天的天气概率分布 B^′(W₁)，此时我们尚未考虑第1天的天气预报。

公式为：B^′(W₁) = ∑_w₀P(W₁|w₀)B(w₀)

第1天是晴天的预测概率，等于”第0天是晴天且第1天也是晴天”的概率，加上”第0天是雨天但第1天是晴天”的概率。

B^′(W₁ = sun) = P(W₁ = sun|W₀ = sun)B(W₀ = sun) + P(W₁ = sun|W₀ = rain)B(W₀ = rain) = (0.6 · 0.8) + (0.1 · 0.2) = 0.48 + 0.02 = 0.5

B^′(W₁ = rain) = P(W₁ = rain|W₀ = sun)B(W₀ = sun) + P(W₁ = rain|W₀ = rain)B(W₀ = rain) = (0.4 · 0.8) + (0.9 · 0.2) = 0.32 + 0.18 = 0.5

经过时间更新后，我们得到的预测置信度为：

W₁	B^′(W₁)
sun	0.5
rain	0.5

这表明，在不考虑任何新证据的情况下，我们预测第1天晴天和雨天的概率各占一半。

第二步是观测更新（修正）: 现在我们引入观测证据 f₁ = good，用它来修正我们的预测，得到最终的置信度 B(W₁)。

公式为：B(W₁) ∝ P(f₁|W₁)B^′(W₁)

B(W₁ = sun) ∝ P(F₁ = good|W₁ = sun)B^′(W₁ = sun) ∝ 0.8 · 0.5 = 0.4

B(W₁ = rain) ∝ P(F₁ = good|W₁ = rain)B^′(W₁ = rain) ∝ 0.3 · 0.5 = 0.15

在第二步之后, 我们得到了未经归一化的置信度值：sun为0.4，rain为0.15。

第三步是归一化, 将上述结果转换为合法的概率分布（即各项概率之和为1）。

计算总和：0.4 + 0.15 = 0.55

归一化： B(W₁ = sun) = 0.4/0.55 = 8/11 B(W₁ = rain) = 0.15/0.55 = 3/11

经过完整的计算，我们得到第1天的最终置信度分布 B(W₁) 如下：

W₁	B(W₁)
sun	8/11 (≈ 0.727)
rain	3/11 (≈ 0.273)

最开始，我们对第1天的天气预测是晴雨概率各占50% (B^′(W₁))。然而，当我们观测到”预报为good”这个证据后，我们对”第1天是晴天”的信心大大增加，从50% (1/2) 上升到了约73% (8/11)。这完全符合直觉，因为晴天时出现”good”预报的概率（0.8）远高于雨天时（0.3）。这个例子清晰地展示了前向算法是如何融合历史信息和当前证据来动态更新对隐藏状态的信念的。

维特比算法（Viterbi Algorithm）

前向算法回答的是“在时刻 N，最可能的状态是什么？”。而 Viterbi 算法回答的是一个不同但同样重要的问题：“给定从时刻1到 N 的所有观测，最可能的状态序列是什么？”

换句话说，算法旨在找到一个完整的隐藏状态序列（一条路径），使得这个序列在给定观测证据的情况下出现的概率最大。

公式化表达即为寻找最可能的隐藏状态序列 x_1 : N^*。 x_1 : N^* = arg max_{x_1 : N}P(x_1 : N|e_1 : N)

想象一下，你每天都会收到一份天气预报（观测证据），但你无法直接知道每天的真实天气（隐藏状态）。几天过后，你手上有一系列的预报记录，比如“（好，好，坏）”。你现在想推断出这几天最可能发生的真实天气序列是什么。是“（晴，晴，雨）”还是“（晴，雨，雨）”？

维特比算法采用动态规划来解决这个问题。其核心思想可以概括为一句话：

“如果要找在时刻 t 到达状态 A 的最优路径，我们只需要知道在时刻 t-1 到达所有可能前置状态（A, B, C, …）的最优路径即可。”

我们不需要关心到达状态 B 或 C 的完整历史路径是什么，只需要知道它们各自由最优路径到达自己的概率是多少。然后我们就可以计算从 A 到 A, 从 B 到 A 和从 C 到 A 的新得分，选择其中更高分的那条作为到达 A 的最优路径。

基本步骤

假设我们有一个隐藏马尔可夫模型，它由以下要素组成：

隐藏状态集合 S = {s₁, s₂, ..., s_N}：我们无法直接观察到的状态。
观测状态集合 O = {o₁, o₂, ..., o_M}：我们能够直接观察到的结果。
初始状态概率 π = {π₁, π₂, ..., π_N}：模型在时间 t=1 时处于每个隐藏状态的概率。
状态转移概率矩阵 A = {a_ij}：从隐藏状态 s_i 转移到 s_j 的概率。
发射概率矩阵 B = {b_j(k)}：在隐藏状态 s_j 下观察到观测 o_k 的概率。

维特比算法的目标是找到最有可能的隐藏状态序列 Q = {q₁, q₂, ..., q_T}，使得 P(Q|O) 的值最大，其中 O = {O₁, O₂, ..., O_T} 是已知的观测序列。

维特比算法分为四个主要步骤：初始化、递归计算、终结和回溯。

步骤一：初始化

在第一个时间步 t=1，我们计算到达每个隐藏状态 s_i 的最大概率。这个概率是该状态的初始概率 π_i 与在该状态下观测到第一个观测值 O₁ 的发射概率 b_i(O₁) 的乘积。

操作：对于每个隐藏状态 s_i，计算维特比变量 δ₁(i): δ₁(i) = π_i ⋅ b_i(O₁)

这一步是为动态规划的起点做准备。我们计算了所有可能的路径中，在时间 t=1 结束于状态 s_i 的最大概率。

步骤二：递归计算

从时间步 t=2 到 T，我们对每个时间步和每个隐藏状态进行迭代计算。对于每个状态 s_j，我们需要考虑所有从上一个时间步 t-1 的状态 s_i 转移到 s_j 的可能性，并选择其中概率最大的路径。

操作：对于每个时间步 t (2≤t≤T) 和每个隐藏状态 s_j，计算维特比变量 δ_t(j) 和回溯指针 ψ_t(j)。

数学公式: δ_t(j) = max_i[δ_t − 1(i) ⋅ a_ij] ⋅ b_j(O_t) ψ_t(j) = arg max_i[δ_t − 1(i) ⋅ a_ij]

δ_t(j) 的计算：max_i[δ_t − 1(i) ⋅ a_ij] 寻找了从前一个时间步的任意状态 s_i 到达当前状态 s_j 的最可能路径。然后将其乘以在当前状态 s_j 下观测到 O_t 的概率 b_j(O_t)，得到在时间 t 结束于状态 s_j 的整体最大概率。

ψ_t(j) 的计算：它记录了在时间 t 达到状态 s_j 的最可能路径是从时间 t-1 的哪个状态 s_i 转移过来的。这个回溯指针是最终重建路径的关键。

步骤三：终结在处理完所有时间步后，我们找到最后一个时间步 T 中，概率最大的那个状态，这将是我们的最优路径的终点。

操作：找到在时间 t=T 时，所有维特比变量 δ_T(i) 中的最大值。

数学公式：

P^* = max_i[δ_T(i)] q_T^* = arg max_i[δ_T(i)]

解释：P^* 是整个观测序列下，最有可能隐藏状态序列的概率。q_T^* 是该最优序列的最后一个状态。

步骤四：回溯

这是重建最优路径的关键步骤。我们从最后一个状态 q_T^* 开始，利用之前记录的回溯指针 ψ_t(j) 逆向推导出整个序列。

操作：从 t=T 向 t=1 循环，通过回溯指针确定前一个状态。

数学公式：

q_t − 1^* = ψ_t(q_t^*)

这一步利用了动态规划的特性。由于我们在每一步都记录了最优的前驱状态，我们只需要从终点倒推，就能唯一地确定一条由局部最优解构成的全局最优路径。

示例分析

让我们用一个经典的”天气-心情”例子来演示维特比算法。

隐藏状态：S = {晴天, 雨天}
观测序列：O = {开心, 难过, 开心}

参数：

初始概率 π：π_晴天 = 0.6, π_雨天 = 0.4
转移概率 A：
- $晴天晴天$
- $晴天雨天$
- $雨天晴天$
- $雨天雨天$
发射概率 B：
- P(开心|晴天) = 0.8
- P(难过|晴天) = 0.2
- P(开心|雨天) = 0.3
- P(难过|雨天) = 0.7

计算过程如下：

时间步 t = 1 (观测：开心): 初始概率乘发射概率来实现初始化

δ₁(晴天) = π_晴天 ⋅ b_晴天(开心) = 0.6 ⋅ 0.8 = 0.48

δ₁(雨天) = π_雨天 ⋅ b_雨天(开心) = 0.4 ⋅ 0.3 = 0.12

这样, 我们计算了第一天是晴天（概率0.48）和雨天（概率0.12）的各自最大概率路径。

时间步 t = 2 (观测：难过): 计算 δ₂ 和 ψ₂

δ₂(晴天) = max [δ₁(晴天) ⋅ a_{晴天, 晴天}, δ₁(雨天) ⋅ a_{雨天, 晴天}] ⋅ b_晴天(难过)

= max [0.48 ⋅ 0.7, 0.12 ⋅ 0.4] ⋅ 0.2 = max [0.336, 0.048] ⋅ 0.2 = 0.336 ⋅ 0.2 = 0.0672

ψ₂(晴天) = 晴天（因为 0.336 > 0.048）

δ₂(雨天) = max [δ₁(晴天) ⋅ a_{晴天, 雨天}, δ₁(雨天) ⋅ a_{雨天, 雨天}] ⋅ b_雨天(难过)

= max [0.48 ⋅ 0.3, 0.12 ⋅ 0.6] ⋅ 0.7 = max [0.144, 0.072] ⋅ 0.7 = 0.144 ⋅ 0.7 = 0.1008

ψ₂(雨天) = 晴天（因为 0.144 > 0.072）

在第二天，我们分别计算了结束于”晴天”和”雨天”的两种可能路径的最大概率。我们发现，要到达”晴天”状态，最可能的前一天状态是”晴天”；要到达”雨天”状态，最可能的前一天状态也是”晴天”。

时间步 t = 3 (观测：开心): 计算 δ₃ 和 ψ₃

δ₃(晴天) = max [δ₂(晴天) ⋅ a_{晴天, 晴天}, δ₂(雨天) ⋅ a_{雨天, 晴天}] ⋅ b_晴天(开心)

= max [0.0672 ⋅ 0.7, 0.1008 ⋅ 0.4] ⋅ 0.8 = max [0.04704, 0.04032] ⋅ 0.8 = 0.04704 ⋅ 0.8 = 0.037632

ψ₃(晴天) = 晴天

δ₃(雨天) = max [δ₂(晴天) ⋅ a_{晴天, 雨天}, δ₂(雨天) ⋅ a_{雨天, 雨天}] ⋅ b_雨天(开心)

= max [0.0672 ⋅ 0.3, 0.1008 ⋅ 0.6] ⋅ 0.3 = max [0.02016, 0.06048] ⋅ 0.3 = 0.06048 ⋅ 0.3 = 0.018144

ψ₃(雨天) = 雨天

终结和回溯

比较 δ₃(晴天) 和 δ₃(雨天), 有0.037632 > 0.018144

因此，最优路径的最后一个状态是晴天。

最后我们回溯路径即可得到结果：

q₃^* = 晴天

q₂^* = ψ₃(q₃^*) = ψ₃(晴天) = 晴天

q₁^* = ψ₂(q₂^*) = ψ₂(晴天) = 晴天

因此, 最有可能的隐藏状态序列是：{晴天, 晴天, 晴天}。

应用领域

维特比算法因其高效和准确性，在许多领域有广泛应用，例如：

语音识别：根据声学模型（隐藏状态）和声音信号（观测）找到最可能的单词序列。
自然语言处理：
- 用于词性标注（Part-of-Speech Tagging），根据单词序列（观测）找到最可能的词性序列（隐藏状态）。
- 用于命名实体识别（Named Entity Recognition），根据文本序列（观测）找到最可能的实体序列（隐藏状态）。
生物信息学：用于基因识别，根据 DNA 序列（观测）找到编码区和非编码区（隐藏状态）。
无线通信：在编码理论中用于信道解码。