FLARE MCMC：基于保真度的层自适应递归MCMC提案

# 基于保真度的自适应分层递归MCMC提议方法
来源：https://arxiv.org/html/2608.13774
## FLARE MCMC：基于保真度的自适应分层递归MCMC提议方法
**致谢：** 本文为作者接受的稿件，原发表于《SIAM不确定性量化期刊》。最终版本可通过以下链接获取：https://doi.org/10.1137/25M1795194
**Harini Venkatesan**  
邮箱：[[email protected]](mailto:[email protected])  
**Christopher Shelton**  
邮箱：[[email protected]](mailto:[email protected])  
**Mengxuan Wu**  
邮箱：[[email protected]](mailto:[email protected])  
**致谢单位：** 加州大学河滨分校计算机科学与工程系
**Ming-Feng Ho**  
邮箱：[[email protected]](mailto:[email protected])  
**致谢单位：** 密歇根大学莱因韦伯理论物理中心
**Simeon Bird**  
邮箱：[[email protected]](mailto:[email protected])  
**致谢单位：** 加州大学河滨分校物理与天文系
**Mengxuan Wu**<sup>1</sup>

#### 摘要
马尔可夫链蒙特卡洛（MCMC）方法仅需评估似然函数即可工作，这使其成为复杂模型推断的常用技术。然而，该方法可能存在混合速率慢的问题，需要生成大量样本才能获得良好估计，总体计算成本较高。FLARE MCMC是一种多保真度分层MCMC方法，它利用真实似然计算的低保真度近似来提升混合性能，从而实现整体更快的计算速度。这种低保真度似然函数在科学与工程应用中十分常见，特别是涉及仿真模拟且其分辨率或精度可调节的模型。本技术采用递归分层链与简单的层级调优机制，不要求似然函数具备特定形式或内部数学结构。实验表明，FLARE MCMC在包括水文学和宇宙学在内的不同科学领域中，能在相同计算时间内实现更大的有效样本量。

**关键词**：马尔可夫链蒙特卡洛、分层模型、贝叶斯推断、基于模拟的推断、水文学、宇宙学

**资金来源**：本研究由美国国家科学基金会（NSF）资助，资助编号IIS-2435579。SB获得NASA-80NSSC21K1840项目支持。
<sup>††</sup>页眉：FLARE MCMC：基于保真度的自适应分层递归MCMC提议方法 / H. Venkatesan, C. Shelton, M. Ho, S. Bird, and M. Wu

**MSC分类**：62F15, 62M05, 65C05, 65C40, 85A35

## 1 引言
马尔可夫链蒙特卡洛（MCMC）是科学与工程计算的重要工具。最常用于根据观测数据计算模型参数的后验分布。计算得到的分布（以样本形式表示）为校准和不确定性量化提供估计，可用于指导新科学实验的生成、阐明模型的可观测性以及解决科学理论问题。在众多MCMC算法中，Metropolis-Hastings MCMC（MH-MCMC）因其能够从几乎任何分布中采样而广受欢迎，且仅需评估给定参数设置下的模型似然函数。然而，这也正是其弱点：该方法缺乏额外的问题背景知识来有效指导采样过程。因此，其混合时间（生成有效新样本的速度）可能较慢，导致整体算法计算负担沉重。

哈密顿蒙特卡洛及其变体[17,47,36]等方法通过添加辅助动量变量来加速混合，允许更长的步长以减少连续样本间的相关性。但此类方法需要计算对数目标分布相对于参数的梯度，当分布通过冗长模拟代码评估时，计算梯度可能极其昂贵。例如，我们在实验结果中使用的宇宙学模拟旨在近似基于SDSS-III重子振荡光谱巡天（BOSS）数据[13,4]的星系功率谱条件后验密度，该模拟无法修改以生成梯度。由于正向宇宙学模拟的复杂性和不可微性，模型参数的梯度无法获得，因此无法应用自动微分或解析梯度形式，这阻碍了基于梯度的推断方法的使用。

FLARE MCMC通过利用同一问题的低保真度模型来加速MH-MCMC的混合时间。许多工程或科学计算模型可以运行在多个保真度级别。FLARE MCMC利用一组计算成本更低的后验计算，每个都是真实后验的近似。许多后验涉及求解偏微分方程、常微分方程或积分问题：对这些问题，粗化空间或时间网格可以产生更廉价的近似；对于包含约束或优化求解器的问题，降低求解器容差或最大迭代次数同样可以产生更廉价的近似。我们进一步展示了一个物理示例，其中可以通过减少代表性粒子数量来粗化底层模拟。通过递归应用MCMC链，我们可以利用较粗分辨率模型来指导较高分辨率的MCMC链。最终结果是为目标模型提供了收敛更快、每计算时间产生更多有效样本的采样器，即使考虑到应用低保真度计算所需的额外时间。

## 2 背景
马尔可夫链蒙特卡洛是一类通过构造易于模拟的马尔可夫链来从复杂目标分布中采样的算法，该马尔可夫链的平稳分布即为目标分布。通常，该目标分布是给定观测条件下一组参数的后验分布。设\(D\)为观测数据，\(\theta \in \Theta \subset \mathbb{R}^R\)为参数。假设参数服从先验分布\(p(\theta)\)，则感兴趣的目标后验分布\(\pi(\theta \mid D)\)通过贝叶斯定理获得：
\[
\pi(\theta \mid D) = \frac{\mathcal{L}(D \mid \theta) p(\theta)}{p(D)} \propto \mathcal{L}(D \mid \theta) p(\theta) \tag{1}
\]
其中\(\mathcal{L}(D \mid \theta)\)是数据的似然函数，在许多科学应用中需要长时间模拟才能评估。我们仅要求能评估\(\pi(\theta \mid D)\)至比例常数，因此可以安全忽略分母\(p(D)\)。后验\(\pi(\theta \mid D)\)作为条件分布这一特性对MCMC基本无关，因此我们直接用\(\pi(\theta)\)表示感兴趣的分布（若底层分布为后验，则等于\(\pi(\theta \mid D)\)，但也可以是\(\theta\)上的任意分布）。

### 2.1 Metropolis-Hastings MCMC
我们首先关注由Hastings（1970）[30]提出的用于马尔可夫链蒙特卡洛的Metropolis-Hastings方法（MH-MCMC）。第\((i+1)\)个样本\(\theta^{i+1}\)基于链中前一个样本\(\theta^i\)通过两步过程生成。首先，从提议分布\(q(\tilde{\theta}^i \mid \theta^i)\)生成提议的下一个状态\(\tilde{\theta}^i\)。然后，根据精心构造的接受概率，决定将\(\tilde{\theta}^i\)接受为\(\theta^{i+1}\)或拒绝。若接受，则\(\theta^{i+1} = \tilde{\theta}^i\)；否则\(\theta^{i+1} = \theta^i\)。通常使用以\(\theta^i\)为中心的正态分布作为提议分布\(q(\tilde{\theta}^i \mid \theta^i)\)，但在温和条件下（例如\(q(\tilde{\theta}^i \mid \theta^i)\)处处为正）几乎任何提议分布均可使用。选定\(q(\tilde{\theta}^i \mid \theta^i)\)后，转移\(\theta^i \to \tilde{\theta}^i\)的接受概率\(\mathcal{A}\)为：
\[
\mathcal{A}(\theta^i \to \tilde{\theta}^i) = \min(1, r(\theta^i \to \tilde{\theta}^i)) \tag{2}
\]
其中
\[
r(\theta^i \to \tilde{\theta}^i) = \frac{\pi(\tilde{\theta}^i)}{\pi(\theta^i)} \frac{q(\theta^i \mid \tilde{\theta}^i)}{q(\tilde{\theta}^i \mid \theta^i)} \tag{3}
\]
尽管标准Metropolis-Hastings MCMC算法可能是从后验分布采样的简便方法，但它需要足够样本才能有效近似后验分布。当链混合缓慢（由于次优提议分布）时，连续样本高度相关，需要采集更多样本才能获得代表真实分布的集合。当评估\(\pi(\theta^i)\)（公式3计算所需）计算代价高昂时，这一问题尤为突出。

### 2.2 相关工作
加速MCMC采样的目标与大量现有工作一致。这些方法包括：链耦合方法（如模拟退火）、旨在使用多保真度廉价近似降低目标估计器方差的方法（如MLMC），以及使用廉价模型构建MCMC提议的方法。与FLARE MCMC类似，模拟退火和耦合MCMC[57,46,5]等方法也使用多条链。通过评估过程的能量并调整模型温度来接受或拒绝样本。两条链在不同温度下并行运行，系统在不同温度间切换。可逆跳跃MCMC[24,3]也在不同维度的链间跳跃。虽然FLARE MCMC具有多链概念，但由于它解决的问题不同（旨在利用计算成本低数个数量级的模拟），其最终结构非常不同。

序列MCMC或粒子滤波[43,16]等方法使用大量称为粒子的样本集合来近似目标，并通过重要性采样在时间上推进这些粒子。但这些是滤波框架，不会收敛到平稳分布。因此，尽管结构看似相关，FLARE MCMC与这些方法存在本质区别。

一个极具影响力的研究方向聚焦于多层级蒙特卡洛（MLMC）框架中降低目标期望最终估计器的方差。受高维参数积分的多层级蒙特卡洛方法[31]及多层级蒙特卡洛路径模拟[22]启发，Hoang等人[34]提出了一种应用于贝叶斯逆问题的多层级MCMC方法。其核心思想是利用分辨率递增的计算模型层级，将高保真期望分解为伸缩和形式。该方法通过使用少量样本估计低方差差值项实现计算加速，而大部分计算资源花费在廉价低保真估计器上。作者还提供了严格的复杂度证明，展示了当在稀疏网格上运行迭代采样器并利用离散化误差的伸缩展开时，后验期望的收敛速度。该MLMC框架已在许多方向得到扩展：多层级序列蒙特卡洛采样器[7]和多层级粒子滤波[39]，结合前期工作[35,25,26]，将MLMC扩展到序列蒙特卡洛；Jasra等人[41]将MLMC扩展到部分观测扩散的静态参数估计问题；多重索引MCMC[29,40]处理了存在多种离散化方式的问题。这些MLMC方法在伸缩估计器中使用所有链的样本，针对特定感兴趣量的均方误差（MSE）。相比之下，FLARE MCMC仅使用最细链的样本（类似下述方法），并以链的混合时间而非MSE为目标。本文的理论分析聚焦于FLARE MCMC的遍历性和收敛速率，并不针对特定问题领域。

与FLARE MCMC更相似的是，多种先前方法表明，用通常具有高接受概率的近似替代提议，可以显著降低标准Metropolis-Hastings算法的计算成本。该思想最早由Christen和Fox[9,20]作为两阶段MCMC方法提出：使用廉价近似测试原始提议，以寻找链中更可能被接受的移动。换言之，在用更昂贵模型评估之前，先使用近似模型的似然函数接受候选样本。Efendiev、Hou和Luo[18]提出的基于粗尺度模拟的预条件MCMC中，采用两阶段测试粗模型（基于高保真多尺度有限体积模型）以减少细保真度产生的计算成本。但该方法仅执行单次廉价近似检查，并未利用其运行完整的MCMC子链。多层级马尔可夫链蒙特卡洛（MLMCMC）[14]通过在细层级实现计算效率：若粗提议被细层级拒绝，粗链独立于细链继续运行，而非递归地从细链当前样本启动下一个粗链。MLMCMC使用用户指定的、在似然计算内部共享跨层级的变量（例如，通过噪声模型与真实观测进行比较的预测观测值）。从粗近似中抽取的样本用于降低该内部变量的方差，从而从粗保真度获得更优的提议。Lykkegaard等人[44]提出了自适应多层级延迟接受（MLDA）...