是否有一个 Monte Carlo/MCMC 采样器可以处理后验分布的孤立局部最大值?

机器算法验证 贝叶斯 马尔可夫链蒙特卡罗 蒙特卡洛 遍历的
2022-03-05 04:29:42

我目前正在使用贝叶斯方法来估计由多个 ODE 组成的模型的参数。因为我有 15 个参数要估计,所以我的采样空间是 15 维的,我搜索的后验分布似乎有许多局部最大值,这些最大值被非常低概率的大区域隔离开来。

这导致了我的蒙特卡洛链的混合问题,因为一条链“跳出”一个局部最大值并意外撞到另一个最大值之一的可能性很小。

在这个领域似乎有很多研究,因为很容易找到处理这个问题的论文(见下文),但很难找到实际的实现。我只找到了与分子动力学相关的软件包,但没有找到贝叶斯推理。是否有能够处理孤立的局部最大值的(MC)MC 采样器的实现?

我被迫使用 Matlab,因为我的 ODE 模型就是用这种方法编写的,因此非常欢迎有关 Matlab 的建议;-)。但是,如果有其他语言的“杀手级应用程序”,也许我可以说服我的 PI 切换 ;-)。

我目前正在使用由Haario、Laine 等人编写的延迟拒绝/自适应蒙特卡洛采样器。,这也是迄今为止我能找到的唯一一个比标准 Metropolis-Hastings 算法更复杂的采样器


值得注意的方法似乎是:

编辑于 2017 年 3 月 7 日更新了我在此期间学到的知识

多个起点不同的相似链

链间适配。使用多个独立链生成的池化样本的经验协方差矩阵来更新链的提议分布的协方差矩阵。(1)

不同回火的多条链条

回火: 某种“温度”似乎会改变后部的景观,使链条更可能混合。(我还没有深入研究这个)(1)回火的目的是为了使后验概率分布形成的(高维)概率景观变平。它通常通过将后验概率取为的幂来完成,其中后验景观在时变平(3,第 298 页)。这意味着,不是计算状态 ,而是给定数据计算缓和后验概率1/TT>1p(θD)θD

p(θD)1/T(p(Dθ)p(θ))1/T

越高,概率景观中的峰越平坦、越宽。因此,较高的值导致采样器从一个局部最大值切换到另一个局部最大值的概率更高。但是,不是在 时搜索的后验分布。因此,必须使用该分布的样本链来启用从之后的采样。TTp(θD)1/TT1p(θD)

来自原始的、未经调和的后验分布的样本,给定来自该分布的调和版本的样本,可以通过几种方法获得:

  • Metropolis 耦合 MCMC同时运行多个链,每个链具有不同但恒定的值。概率性地切换两条链的状态。的样本进行下游估计;其他链只是确保对所有峰进行采样。参考。(4) 有一个并行算法,并引用了一篇会议文章和一本教科书的想法 (5,6)TT=1

  • 小世界 MCMC。采样器在两个提议之间切换。大多数情况下,使用方差小的提案分布,很少使用方差大的提案。这两个提议之间的选择是随机的。具有大方差的建议也可以从另一个链中提取,该链只会产生非常大的跳跃,以粗略的方式尽可能多地采样样本空间。(2,7)

哈密​​顿蒙特卡洛 (HMC)

我对此了解不多,但JAGS的 No-U-Turn 采样器 (NUTS) 似乎使用它。参见参考文献。(8)。Alex Rogozhnikov 创建了一个关于该主题的可视化教程。


参考:

(1) Craiu 等人,2009:向你的邻居学习:平行链和区域自适应 MCMC。 J Am Stat Assoc 104:488,第 1454-1466 页。http://www.jstor.org/stable/40592353

(2) Guam et al., 2012: Small World MCMC with tempering: E​​rgocity and Spectrum gap。https://arxiv.org/abs/1211.4675仅在 arXiv 上

(3):布鲁克斯等人。(2011)。马尔可夫链蒙特卡罗手册。CRC出版社。

(4):Altekar 等人。(2004):平行都市耦合马尔可夫链蒙特卡罗用于贝叶斯系统发育推断。生物信息学20 (3) 2004,第 407–415 页,http: //dx.doi.org/10.1093/bioinformatics/btg427

(5):Geyer CJ (1991) 马尔可夫链蒙特卡罗最大似然。在:Keramidas(编辑),计算科学与统计:第 23 届接口研讨会论文集界面基金会,费尔法克斯站,第 156-163 页。

(6):吉尔克斯 WR 和罗伯茨 GO (1996)。改进 MCMC 的策略。在: Gilks​​ WR、Richardson S 和 Spiegelhalter (eds) Markov chain Monte Carlo in Practice查普曼和霍尔,p。89–114。

(7):关永等。小世界中的马尔可夫链蒙特卡罗。统计与计算(2006) 16 (2),第 193-202 页。http://dx.doi.org/10.1007/s11222-006-6966-6

(8): Hoffmann M and Gelman A (2014): The No-U-Turn Sampler: Adaptively Setting Path Lengths in Hamiltonian Monte Carlo。 机器学习研究杂志,15,第 1351-1381 页。https://arxiv.org/abs/1111.4246

1个回答

上述两种策略都不是特别适合多重最优。

更好的选择是差分进化 MCMC 和派生的 MCMC,例如 DREAM。这些算法与多个混合生成提案的 MCMC 链一起工作。如果您在每个最优值中至少有一条链,它们将能够在最优值之间有效地跳转。R 中的实现可在此处获得https://cran.r-project.org/web/packages/BayesianTools/index.html