8.4 粒子滤波 (Particle Filtering)

回想一下，对于贝叶斯网络，当运行精确推理计算成本太高时，使用我们讨论过的采样技术之一是有效近似我们想要的所需概率分布的可行替代方案。隐马尔可夫模型也有同样的缺点——使用前向算法运行精确推理所需的时间与随机变量定义域中值的数量成正比。在我们当前的天气问题公式中，天气只能取 2 个值，\(W_i \in \{sun, rain\}\)，这是可以接受的，但假设我们想运行推理来计算某一天实际温度的分布，精确到十分之一度。

贝叶斯网络采样的隐马尔可夫模型模拟称为粒子滤波 (particle filtering)，它涉及模拟一组粒子通过状态图的运动，以近似相关随机变量的概率（信念）分布。这解决了与前向算法相同的问题：它给了我们 \(P(X_N \mid e_{1:N})\) 的近似值。

我们不再存储将每个状态映射到其信念概率的完整概率表，而是存储 \(n\) 个粒子 (particles) 的列表，其中每个粒子处于我们时间依赖随机变量定义域中的 \(d\) 个可能状态之一。通常，\(n\) 明显小于 \(d\)（符号表示为 \(n << d\)），但仍然足够大以产生有意义的近似值；否则，粒子滤波的性能优势将变得微不足道。粒子只是该算法中样本的名称。

我们对粒子在任何给定时间步处于任何给定状态的信念完全取决于我们的模拟中该时间步处于该状态的粒子数量。例如，假设我们确实想模拟某一天 \(i\) 的温度 \(T\) 的信念分布，并为简单起见假设此温度只能取范围 \([10, 20]\) 内的整数值（\(d = 11\) 个可能状态）。进一步假设我们有 \(n = 10\) 个粒子，它们在模拟的时间步 \(i\) 取以下值：

\[[15, 12, 12, 10, 18, 14, 12, 11, 11, 10]\]

通过计算我们的粒子列表中出现的每个温度的计数并除以粒子总数，我们可以生成我们想要的时间 \(i\) 的温度经验分布：

\(T_i\)	10	11	12	13	14	15	16	17	18	19	20
\(B(T_i)\)	0.2	0.2	0.3	0	0.1	0.1	0	0	0.1	0	0

现在我们已经看到了如何从粒子列表恢复信念分布，剩下的就是讨论如何为我们选择的时间步生成这样的列表。

8.4.1 粒子滤波模拟 (Particle Filtering Simulation)

粒子滤波模拟从粒子初始化开始，这可以非常灵活地完成——我们可以随机、均匀或从某个初始分布中采样粒子。一旦我们采样了初始粒子列表，模拟就采用与前向算法类似的形式，在每个时间步进行时间流逝更新，然后进行观测更新：

时间流逝更新 (Time Elapse Update) — 根据转移模型更新每个粒子的值。对于处于状态 \(t_i\) 的粒子，从 \(P(T_{i+1} \mid t_i)\) 给出的概率分布中采样更新值。注意时间流逝更新与贝叶斯网络先验采样的相似性，因为任何给定状态下粒子的频率反映了转移概率。
观测更新 (Observation Update) — 在粒子滤波的观测更新期间，我们使用传感器模型 \(P(F_i \mid T_i)\) 根据观测到的证据和粒子状态所指示的概率对每个粒子进行加权。具体来说，对于处于状态 \(t_i\) 且传感器读数为 \(f_i\) 的粒子，分配权重 \(P(f_i \mid t_i)\)。观测更新的算法如下：
1. 如上所述计算所有粒子的权重。
2. 计算每个状态的总权重。
3. 如果所有状态的所有权重之和为 0，则重新初始化所有粒子。
4. 否则，归一化状态上的总权重分布，并从该分布中重新采样你的粒子列表。

注意观测更新与似然加权的相似性，我们再次根据我们的证据降低样本的权重。

让我们看看我们是否可以通过示例稍微更好地理解这个过程。使用温度作为时间依赖随机变量为我们的天气场景定义一个转移模型如下：对于特定的温度状态，你可以保持在同一状态或转移到一度之外的状态，在范围 \([10, 20]\) 内。在可能的合成状态中，转移到最接近 15 的状态的概率为 80%，其余合成状态在其自身之间均匀分配剩余的 20% 概率。

我们的温度粒子列表如下：

\[[15, 12, 12, 10, 18, 14, 12, 11, 11, 10]\]

为了对该粒子列表中的第一个粒子（处于状态 \(T_i = 15\)）执行时间流逝更新，我们需要相应的转移模型：

\(T_{i+1}\)	14	15	16
\(P(T_{i+1} \mid T_i = 15)\)	0.1	0.8	0.1

实际上，我们为 \(T_{i+1}\) 定义域中的每个值分配不同的值范围，以便这些范围一起完全跨越区间 \([0, 1)\) 且不重叠。对于上述转移模型，范围如下：

\(T_{i+1} = 14\) 的范围是 \(0 \leq r < 0.1\)。
\(T_{i+1} = 15\) 的范围是 \(0.1 \leq r < 0.9\)。
\(T_{i+1} = 16\) 的范围是 \(0.9 \leq r < 1\)。

为了重新采样我们处于状态 \(T_i = 15\) 的粒子，我们只需生成一个范围 \([0, 1)\) 内的随机数，看看它落在哪个范围内。因此，如果我们的随机数是 \(r = 0.467\)，那么 \(T_i = 15\) 处的粒子保持在 \(T_{i+1} = 15\)，因为 \(0.1 \leq r < 0.9\)。现在考虑区间 \([0, 1)\) 中的以下 10 个随机数列表：

\[[0.467, 0.452, 0.583, 0.604, 0.748, 0.932, 0.609, 0.372, 0.402, 0.026]\]

如果我们使用这 10 个值作为重新采样我们 10 个粒子的随机值，那么完全时间流逝更新后的新粒子列表应该如下所示：

\[[15, 13, 13, 11, 17, 15, 13, 12, 12, 10]\]

自己验证一下！更新后的粒子列表产生了相应的更新后的信念分布 \(B(T_{i+1})\)：

\(T_i\)	10	11	12	13	14	15	16	17	18	19	20
\(B(T_{i+1})\)	0.1	0.1	0.2	0.3	0	0.2	0	0.1	0	0	0

将我们更新后的信念分布 \(B(T_{i+1})\) 与我们的初始信念分布 \(B(T_i)\) 进行比较，我们可以看到作为一个总体趋势，粒子倾向于收敛到 \(T = 15\) 的温度。

接下来，让我们执行观测更新，假设我们的传感器模型 \(P(F_i | T_i)\) 指出正确预测 \(f_i = t_i\) 的概率为 80%，预测其他 10 个状态中的任何一个的概率均匀为 2%。假设预测为 \(F_{i+1} = 13\)，我们 10 个粒子的权重如下：

Particle	\(p_1\)	\(p_2\)	\(p_3\)	\(p_4\)	\(p_5\)	\(p_6\)	\(p_7\)	\(p_8\)	\(p_9\)	\(p_{10}\)
State	15	13	13	11	17	15	13	12	12	10
Weight	0.02	0.8	0.8	0.02	0.02	0.02	0.8	0.02	0.02	0.02

然后我们按状态聚合权重：

State	10	11	12	13	15	17
Weight	0.02	0.02	0.04	2.4	0.04	0.02

将所有权重的值相加得出总和 2.54，我们可以通过将每个条目除以此总和来归一化我们的权重表以生成概率分布：

State	10	11	12	13	15	17
Weight	0.02	0.02	0.04	2.4	0.04	0.02
Normalized Weight	0.0079	0.0079	0.0157	0.9449	0.0157	0.0079

最后一步是从这个概率分布中重新采样，使用我们在时间流逝更新期间用于重新采样的相同技术。假设我们在范围 \([0, 1)\) 内生成 10 个随机数，其值如下：

\[[0.315, 0.829, 0.304, 0.368, 0.459, 0.891, 0.282, 0.980, 0.898, 0.341]\]

这产生了一个重新采样的粒子列表如下：

\[[13, 13, 13, 13, 13, 13, 13, 15, 13, 13]\]

具有相应的最终新信念分布：

\(T_i\)	10	11	12	13	14	15	16	17	18	19	20
\(B(T_{i+1})\)	0	0	0	0.9	0	0.1	0	0	0	0	0

观察到我们的传感器模型编码了我们的天气预测以 80% 的概率非常准确，并且我们的新粒子列表与此一致，因为大多数粒子被重新采样为 \(T_{i+1} = 13\)。