软件X 14（2021）100664原始软件出版物CU-MSDSp：一种灵活的并行可逆跳马尔可夫链蒙特卡罗方法约翰·泰勒·查维斯三世a，b，艾米·路易丝·科克伦b，c，克里斯托弗·詹姆斯·厄尔斯a，daCenter for Applied Mathematics，Cornell University，Ithaca，NY 14850，United Statesb人口健康科学系，威斯康星大学麦迪逊分校，威斯康星州53726，美国威斯康星大学数学系，Madison，WI 53706，美国d康奈尔大学土木与环境工程学院，伊萨卡，纽约州14850，美国ar t i cl e i nf o文章历史记录：收到2020年2021年1月20日收到修订版2021年1月21日接受保留字：可逆跳马尔可夫链蒙特卡罗方法RJMCMC模型选择并行MCMCa b st ra ct可逆跳马尔可夫链蒙特卡罗（RJMCMC）是一种功能强大的贝叶斯跨维算法，用于执行模型选择，同时推断模型参数的分布。 目前的工作介绍CU-MSDSp作为一个开源和全自动的并行RJMCMC实现，旨在增加从业者的RJMCMC的可访问性。CU-MSDSp通过独立地形成马尔可夫链来近似每个模型参数的后验分布。然后，这些近似值用于估计模型空间的后验分布。这一可并行化的软件消除了设计跨维建议分布和雅可比矩阵的需要，同时确保与非并行RJMCMC算法相同的理论保证。最后，CU-MSDSp使从业者能够依靠他们以前的知识，固定维MCMC收敛评估和仿真设计。版权所有©2021作者。由爱思唯尔公司出版这是CC BY许可下的开放获取文章（http://creativecommons.org/licenses/by/4.0/）中找到。代码元数据当前代码版本v1用于此代码版本的代码/存储库的永久链接https://github.com/ElsevierSoftwareX/SOFTX-D-20-00033Code Ocean compute capsule法律代码许可证MIT许可证使用git的代码版本控制系统软件代码语言、工具和服务使用C、python、MPI、bash、C++、Stan编译要求、操作环境依赖性符合POSIX标准的平台上的GSL、BLAS库（如果可用）开发人员文档/手册链接问题支持电子邮件jtc278@cornell.edu1. 动机和意义模型选择是一个具有挑战性的问题，经常出现在应用科学和工程，当从一个系统的观察是使用参数化模型来解释。当进行这样的选择时，期望用于量化竞争模型的准确性的原则手段可逆跳马尔可夫链蒙特卡罗（RJMCMC）方法是一种贝叶斯框架，通过提供计算模型空间上的联合后验分布估计以及实例化它们的相关模型参数，来辅助模型选择。*通讯作者。电子邮件地址：jtc278@cornell.edu（John Taylor Chavis III）.https://doi.org/10.1016/j.softx.2021.100664RJMCMC自1995年推出以来一直受到欢迎[1]。自早期发表详尽的文献综述[2]以来的几年中，文献中出现了2500多篇引用GreenRJMCMC已被用于各种学科，如统计/概率，计算机科学，工程，生物学和计算机图形学[3，4]。虽然科学界对RJMCMC的兴趣是意料之中的，但有点出乎意料的是，由于其应用带来的令人印象深刻的理论好处，它的使用并不广泛。我们认为，这种状况是一个人工制品的综合因素，使实施RJMCMC复杂，并创造困难时，解释结果从RJMCMC分析。2352-7110/©2021作者。 由Elsevier B.V.出版。这是一篇开放获取的文章，使用CC BY许可证（http://creativecommons.org/licenses/by/4.0/）。可在ScienceDirect上获得目录列表SoftwareX期刊主页：www.elsevier.com/locate/softxJohn Taylor Chavis III，Amy Louise Cochran和Christopher James Earls软件X 14（2021）1006642||=∈|=<${ }×∈K{M∈K}K={}|=|||||尽管RJMCMC有能力导航存在于潜在不同维度（即，trans-dimensional）并随后收敛到正确的平稳联合后验分布，但是目前使用该方法存在显著的实际障碍，例如：（1）如何将复杂的数学公式转换成可靠的代码;以及（2）当将该方法应用于实际数据时，如何正确地评估收敛性和准确性。例如，个体模型马尔可夫链是否包含足够有代表性的样本？该方法是否适当考虑了所有竞争模型为了克服这些问题，本文介绍了我们的软件，康奈尔大学-模型选择相关搜索虽然标准RJMCMC算法将模型指数的更新与模型参数的更新混合在一起，但CU-MSDSp以不同的方式评估这种联合后验分布[1]。首先，独立地形成马尔可夫链以近似与每个候选模型相关联的模型参数的后验分布π（φkk，x这些链被假定为已经适当地收敛到它们的目标平稳分布，并且因此被称为“金标准”链。然后，将这些金标准参数链合并到单独的并行RJMCMC算法中，以近似后验分布π（kx）模型指数。 重要的是，后一种分布可以从金本位链中恢复，因为并行（CU-MSDSp）。CU-MSDSp是一个开源的，完全一个自动化的软件系统，提供了一个可并行化的RJMCMC实现。在高级别上，CU-MSDSpπ（k|x）=p（x）|k）的p（x）p（x）|k，φk）p（φk|（k）、（2）π（φk|k，x）p（x）首先通过C++ Stan的流行应用程序编程接口（API）独立地形成马尔可夫链来近似每个模型参数的后验分布然后，这些近似值用于估计模型指数的后验分布。我们的软件提供以下贡献：通过整合广泛使用的斯坦语，增加非专业人士对RJMCMC的访问一种直观的并行RJMCMC方法的软件，大大提高了标准公式的计算效率;一种RJMCMC方法，它消除了设计标准方法中所需的跨维建议分布和雅可比矩阵的需要，同时确保MCMC的理论要求，如详细平衡，保持不变;一种RJMCMC方法，该方法基于从业者本文系统地介绍了我们的软件CU-MSDSp。第2节介绍了CU-MSDSp并行化RJMCMC的理论基础。我们在第3节中给出了一个说明性示例，以突出CU-MSDSp的灵活性和实现功能。最后，在第5节中，我们将讨论CU-MSDSp的影响和未来发展方向。2. 软件描述CU-MSDSp的核心功能是采用一组合理的参数模型，以及来自感兴趣现象的观察结果，随后使用RJMCMC形成模型空间的联合后验，以及相关的模型内参数。具体来说，设x是观测数据的向量，k，k，1， 2， 3，表示可数可以描述x的候选模型的集合。每个模型索引为k，并由nk个参数的集合实例化，形成向量φkRnk。重要的是，每个模型的维数nk可能 最后，假设马尔可夫状态是一个元素θk（k，φk），它由模型索引和相关的模型参数组成，位于状态空间Θk（kRnk）内。RJMCMC的目标是估计模型指数和模型参数的联合后验分布，即。 π（θkx），通过贝叶斯定理的应用，可以用似然p（x）表示|k，φk）和联合先验p（k，φk）= p（φk|k）p（k）：p（x）|k，φk）p（φk|k）p（k）对于任意φk，π（φkk，x）> 0。最后，能够近似π（k x）和π（φkk，x）允许近似目标联合分布：π（θkx）π（φkk，x）π（k x）。通过分别近似模型指数和模型参数的后验分布，RJMCMC方法将马尔可夫链收敛性的总体评估简化为模型指数的马尔可夫链和模型参数的马尔可夫链的收敛性的特别地，用于评估固定维MCMC中的收敛性的任何工具这种简化有望帮助研究人员更好地识别从RJMCMC获得的结果何时收敛，因此对重复分析具有鲁棒性。我们相信，虽然标准MCMC收敛诊断是不完美的，但它们是从业者熟悉的，因此允许这些标准方法在复杂的RJMCMC环境中应用此外，我们观察到CU-MSDSp如何消除RJMCMC的一些典型的跨维问题，通过令人尴尬地并行化模型内参数和模型间概率估计，作为一系列独立的固定MCMC模拟。这种精简包括消除设计建议分布和雅可比矩阵的需要，以确保各种模型空间之间的同构的存在图1给出了我们提出的并行方法的示意图。RJMCMC的完整理论基础超出了本文的范围，但我们建议读者参考几部作品[1，3，6]。虽然图2的并行化可能看起来不起眼，但并行化与模型选择自动化的结合是我们软件提供的一个重要贡献。特别是，我们已经消除了设计跨维提案分布和雅可比矩阵以确保存在性的需要，的一个例子。这在许多模型（如图形模型和有限高斯混合模型[7，8]）的应用中变得非常重要。CU-MSDSp的一些理论基础也作为CU-MSDSp存储库中的补充文档提供它还指出，像标准MCMC，RJMCMC确实遭受的维数灾难，因此作为一个实际问题，建议限制应用程序的模型需要少于十个参数的实例。我们不是第一个实现并行RJMCMC算法。最值得注意的是YeThompson 2009的工作[6]。它们的并行实现需要使用我们提出的金标准链的等价物π（θk|x）=p（x）。（一个）后验模式或最大似然，以便估计模型对概率。我们的工作不同于这个早期的贡献分母p（x）是确保联合后验分布π（θk）的归一化项|x）积分为1。在两个重要方面：第一，我们提供了一个开放的，完全自动化的软件库，它是可推广到特定问题的=····John Taylor Chavis III，Amy Louise Cochran和Christopher James Earls软件X 14（2021）1006643Fig. 1. 并行化RJMCMC的示意图。每个模型的联合参数后验分布通过固定维MCMC方法（Stan）估计。然后，这些金标准链用于确定模型空间的后验分布所有链都使用MPI并行采样并且不需要用户对他们的问题进行“硬编码”其次，CU-MSDSp利用完整的黄金标准链进行模型选择。这种在金标准链上的爬行为我们提供了与标准RJMCMC算法相同的关于平稳分布和目标分布重合的理论保证;尽管CU-MSDSp需要使用核密度估计来估计后验分布（计算工作量的适度增加在应用中实现了更大的灵活性）。2.1. 软件构架CU-MSDSp是开源的，由C++，C，Bash和Python软件元素的混合构建。软件语言的选择是基于速度和效率的考虑。特别是，Bash脚本被选为CU-MSDSp接口语言，因为它的速度和低级别的功 能 。 在 更 高 的 层 次 上 ， CU-MSDSp 采 用 单 程 序 多 数 据（SPMD）编程模型，使用消息传递接口John Taylor Chavis III，Amy Louise Cochran和Christopher James Earls软件X 14（2021）1006644===∼=∼=-=-图二. CU-MSDSp的文件结构（MPI）的实现。MPI是用于高性能计算（HPC）的通信协议，其允许消息分组在分布式计算节点（例如，主节点“等级0”与所有其它等级（工作者）双向通信。图2突出显示了用户在安装时将看到的文件结构。在执行过程中（第一阶段），CU-MSDSp首先通过标准固定维MCMC实现[ 9 - 11 ]估计I链的2 T样本来充分探索K模型如上所述，这些链条被称为我们的黄金标准链条。固定维MCMC通过调用C++的Stan应用程序编程接口（API）在CU-MSDSp中实现[5]。Stan是一个最先进的、广泛使用的开源API，可以自动化贝叶斯MCMC采样。然而，Stan目前没有RJMCMC能力。我们建议用户参考[12]以获得有关如何安装和实现Stan的具体指南。从CU-MSDSp的角度来看，用户只需要通过Stan语言定义他们的模型和数据。用户还需要定义其模型的对数似然我们注意到，除了前面提到的并行化、多线程的使用可以在Stan [5]中使用。实际上，MCMC的任何加速都是对马尔可夫链潜在长期运行的改进[6，13]。一旦用户使用Stan定义了他们的模型，CU-MSDSp就会并行编译和运行每个模型所得.csv输出然后被用作CU-MSDSp第二阶段的黄金标准链：模型选择。第二阶段的目标是使用这些金标准链来近似模型的后验概率。这是通过MPI并行实现的我是工人。这些I工作人员每个实例化一个马尔可夫链来估计模型的概率，通过使用一个特殊的情况下，绿色bash脚本负责编写、编译和执行在每个模型之间正确跳转的C程序。各模型参数后验分布的数值估计 为了估计模型选择概率，通过格林的RJMCMC算法。目前，CU-MSDSp仅实现高斯核密度估计，但用户愿意实现他们选择的核密度（例如有界支持的核密度）。内核函数位于bodypRJMCMC. c文件中。2.2. 软件功能除了估计模型和模型参数的联合后验分布之外，CU-MSDSp系统架构还包括几个有用的功能，这些功能为用户在RJMCMC实施中提供了更大的灵活性，同时还为决策提供了有用的帮助。2.2.1. MCMC诊断CU-MSDSp版本1.0不提供自动MCMC收敛诊断。但是，CU-MSDSp产生的样本输出允许用户轻松地使用他们选择的另一个程序进行性能收敛诊断一个简单的选择是利用Stan中内置的收敛诊断2.2.2. 可视化工具用户可以通过CU-MSDSp中的工具可视化他们的结果特别是，用户可以在模型选择和模型参数空间中可视化每个模型参数的边缘后验分布的迹线图最后，提供了模型选择（第二阶段）的接受概率。这些可视化是支持贝叶斯统计分析的标准输出。特别是，我们的软件自动提供这些数字，并随后允许从业者专注于建模决策，如先验的选择，候选人生成分布等。.3. 说明性实例我 们 现 在 在 RJMCMC 文 献 中 的 一 个 典 型 示 例 上演示CU-MSDSp：确定计数数据的最佳模型。有趣的问题是：在英格兰超级足球比赛中，在特定的时间间隔内，总进球数是否过于分散？具体而言，数据是否可以使用Poison或Negative Binomial分布进行最佳描述？让x（x1，. . .，xN），表示来自2005-2007赛季期间的N1040场英格兰超级足球比赛的数据。进一步地，设一个以k为索引的模型，1代表进球数数据作为泊松随机变量X的独立同分布（iid）实现，平均值为λ：X泊松（λ）。 在这种情况下，马尔可夫状态由θ1（1，λ）表示。类似地，假设一个以k为索引的模型2代表相同目标计数数据作为具有均值λ和参数κ（影响方差）的重新参数化负二项随机变量Y的实现：XNegBin（λ，κ）。在这种情况下，马尔可夫状态是θ2（2，λ，κ）。我们引导感兴趣的读者到文献中，为所需的似然性的理论推导，并选择先验，对于这个例子问题[14]。在频率论中，可以使用泊松模型比负二项模型更受欢迎的零假设的似然比检验来解决这个问题似然比检验结果（χ2（1）0的情况。706页0的情况。407）将表明我们不能在5%的显著性水平上拒绝上述零假设。我们还注意到，其他技术，如过度分散的Cameron和Trivedi检验，可用于评估该特定示例的模型选择[15]。然而，我们选择这个典型的例子来突出CU-MSDSp的特征，并允许读者将我们的结果与RJMCMC文献[14]进行比较。我们通过将CU-MSDSp应用于JSON格式的观测值x来评估两个竞争的统计模型，该格式与Stan兼容。必要的JSON文件被放置在CU-MSDSp的数据目录中，使用标准名称data.json。一旦数据被正确定义，并以正确的John Taylor Chavis III，Amy Louise Cochran和Christopher James Earls软件X 14（2021）1006645目录，用户必须接下来定义竞争模型，以供RJMCMC考虑。这种模型规范使用Stan在其API中使用的数据结构来实现，并且随后采取单独文件的形式（每个计算模型一个），这些文件被放置在CU-MSDSp的models文件夹泊松模型的一个例子，写在斯坦兼容的格式，如下所示。再次提醒读者，模型的对数似然必须定义为Stan中的生成量（作为logModel）。/ / Poisson模型数据{我n不int maximum =0>maximum [N];}参数{实际下限=0>lambda;}模型{/ /Gamma优先于lambda形状和尺度lambda ~gamma（25.0，10。0）的范围内;（lambda）;}生成量{1}真正的logModel;logModel=Poisson_lpmf（y|λ）+gamma_lpdf（lambda|二十五0，10。0）的范围内;图三. CU-MSDSp生成的泊松与负二项分布的模型后验概率的可视化CU-MSDSp提供了这种可视化，允许用户快速评估他们的模型选择。}我们在具有T=50，000的I= 5链上运行算法样品，在50，000个初始“老化”取样后收集。这个例子是运行在Macbook Pro与2.4 GHz的四核心英特尔酷睿i5和16 GB RAM。只需要四个“处理器核心”就可以执行Stan API来解决这个问题。用户必须在setup.config文件中指定其计算环境。此配置文件还包含Stan的绝对路径。要从命令行执行CU-MSDSp，请键入以下内容：sh runParallel每个模型的黄金标准链的输出可以在新创建的modelSelection文件夹中找到从ModelSelection阶段输出的接受概率和模型空间自动生成的文件modelIndex.txt包含模型名称到CU-MSDSp使用的索引此外，CU-MSDSp将默认显示分析结果。这些可视化被放置在目录图片。例如，泊松分布的模型分布空间、接受概率、模型迹线图和对数似然图显示在图1和图2中。3-六、同样，我们注意到这些可视化是用于评估RJMCMC模型参数和模型索引结果的典型结果这些自动绘图使用户能够专注于解释其模拟输出的质量并微调其模型。4. 影响CU-MSDSp为研究人员和从业者提供了一种开源、可访问和便携的方法，当遇到可能用几个竞争模型描述的观测数据时，可以实现原则模型选择（附带模型参数估计）。CU-MSDSp中的RJMCMC实现遵循了作者对RJMCMC的扩展;因此，用户可以更加确信该方法的性能令人满意（即，收敛）。见图4。II期模型选择的CU-MSDSp生成的泊松与负二项分布的验收概率迹线图。这种可视化使从业者能够快速评估MCMC模拟的质量。目前，RJMCMC只有少数几种跨维收敛诊断[8，16尽管这些诊断具有令人信服的性质，但它们确实存在缺点（例如，不频繁的模型访问），这削弱了用户对其性能的信心。CU-MSDSp去除了RJMCMC的跨维性质。在RJMCMC收敛中，人们同时关注John Taylor Chavis III，Amy Louise Cochran和Christopher James Earls软件X 14（2021）1006646图五. 参数分布和迹线绘制CU-MSDSp生成的泊松分布的lambda参数。Stan输出的可视化允许用户评估模型内参数MCMC收敛性。见图6。由CU-MSDSp生成的泊松分布的对数似然迹图。这种标准的可视化在评估模型参数如何描述数据方面非常重要它们的模拟输出的模型间和模型内混合。CU-MSDSp消除了这一问题，并允许用户依靠其固定维MCMC直觉来评估单个模型参数的收敛性（CU-MSDSp的第一阶段），然后评估离散模型概率空间的收敛性这两个阶段不再具有跨维度的诅咒最后，CU-MSDSpCU-MSDSp还5. 结论可逆跳马尔可夫链蒙特卡罗（RJMCMC）是一个强大的工具，用于在竞争参数模型之间进行原则性区分，这些模型在解释来自感兴趣的系统/现象的观察时是合理的。尽管RJMCMC有许多积极的特性，但它并没有在科学界获得人们所期望的吸引力。我们认为，这源于该方法复杂的理论基础以及将该方法的“数学转化为代码”的难度。本文介绍了一个开源的通用RJMCMC软件套件CU-MSDSp来处理这些困难。CU-MSDSp利用作者尽管CU-MSDSp能够消除维度的诅咒，但该软件在收敛方面仍然受到任何MCMC方法的相同限制虽然理论上保证了所有因此，所有MCMC方法，无论是跨维还是固定维，都必须依赖于只能显示收敛失败的收敛诊断。然而，尽管有这些限制，我们希望CU-MSDSp将使专家和非专家都受益。竞合利益作者声明，他们没有已知的竞争性财务利益或个人关系，可能会影响本文报告的工作致谢JTCIII衷心感谢由美国国防部/研究与工程部（国防研究与工程部副部长）和美国国防教育计划（NDEP）/基础研究BA-1资助的科学数学与研究转型（SMART）作者还感谢Alfred P。斯隆基金会（MPHD）计划，美国（第70481号裁决）。引用[1] 绿 色 睡 衣 可 逆 跳 马 尔 可 夫 链 蒙 特 卡 罗 计 算 和 马 尔 可 夫 模 型 确 定 。Biometrika1995;82（4）：711-32.得双曲正切值.doi.org/10.1093/biomet/82.4.711网站。[2]绿色PJ，Hastie DI。可逆跳变MCMC。遗传学2009;155（3）：1391 - 403.[3]Sisson SA.超维马尔可夫链：十年的进展和未来展望。美国统计学会杂志2005;100（471）：1077-89。网址：//dx.doi.org/10.1198/016214505000000664网站。[4]放大图片作者：Jakob W，J.使用逆映射的可逆跳跃Metropolis光传输ACMTrans Graph 2017;37（1）.网址：//dx.doi.org/10.1145/3132704网站。[5]Carpenter B，Gelman A，Hoffman MD，Lee D，Goodrich B，BetancourtM，et al. Stan：A probabilistic programming language. J Stat Softw 2017;76（1）. http://dx.doi.org/10.18637/jss.v076.i01网站。John Taylor Chavis III，Amy Louise Cochran和Christopher James Earls软件X 14（2021）1006647[6] Ye J，Wallace A，Thompson J.变维信号分析的并行马尔可夫链蒙特卡罗计算。2009年第17届欧洲信号处理会议。IEEE; 2009，p.2673-7[7] Brooks SP，Gelman A.迭代模拟收敛监控的一般方法。J Comput Graph Stat1998;7（4）：434-55.http://dx.doi的网站。org/10.1080/10618600.1998.10474787。[8] Brooks SP，Giudici P.通过双向方差分析的马尔可夫链蒙特卡罗收敛性评估。JComputGraphStatistist2000;9（2）：266-85.http://dx.doi.org/10.1080/10618600.2000.10474880网站。[9] 杨文伟，李文伟.马尔可夫链蒙特卡罗的固定宽度输出分析。美国统计学会杂志2006;101（476）：1537-47。http://dx.doi.org/10.1198/016214506000000492网站。[10]龚·弗莱格尔·JM马尔可夫链蒙特卡罗模拟的相对固定宽度停止规则。StatististSinica2015;655-75.[11]Gong L，Flegal JM.高维马尔可夫链蒙特卡罗的一个实用序贯停止规则。JComputGraphStatist2016;25（3）：684-700.http://dx.doi.org/10.1080/10618600.2015.1044092网站。[12]Stan开发团队Cmdstan：stan的命令行界面。2018年，版本2.18。0的情况。[13]Cowles MK，Carlin BP.马尔可夫链蒙特卡罗收敛诊断：比较评论。美国统计学会杂志1996;91（434）：883网址：//dx.doi.org/10.1080/01621459.1996.10476956网站。[14] Hastie DI，绿色PJ。使用可逆跳马尔可夫链蒙特卡罗模型选择。Stat Neerl2012;66 （ 3 ） ： 309-38. http://dx.doi.org/10.1111/j 的 网 站 。 1467-9574.2012.00516.x。[15]卡梅隆AC，特里维迪PK。泊松模型中过度离散的回归检验。J Econometrics1990;46 （ 3 ） ： 347-64. http://dx.doi.org/10 的 网 站 。 1016/0304-4076 （ 90 ）90014-K。[16]卡斯特尔·JM，齐默尔曼·DL。可逆跳MCMC采样器的收敛性评估。统计和精算科学系，爱荷华大学，技术报告313; 2002。[17]放大图片作者：James S. MCMC模型选择的非参数收敛评估《计算图统计学杂志》2003;12（1）：1网址：//dx.doi.org/10.1198/1061860031347网站。[18]Sisson SA，Fan Y.跨维马尔可夫链的一个基于距离的诊断。Stat Comput2007;17（4）：357-67. http://dx.doi.org/10.1007/s11222-007-9025-z.[19]Deonovic BE，Smith BJ.分类变量MCMC绘制的收敛性诊断。2017年，arXiv预印本arXiv：1706.04919。[20]听到NA，Turcotte MJ。来自竞争样本的蒙特卡罗分布估计的收敛性。StatComput 2016;26（6）：1147-61. http://dx.doi.org/10.1007/s11222-015-9595-0.
我的内容管理
展开
