开始使用粒子Metropolis Hastings进行推理

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

英文标题：
《Getting Started with Particle Metropolis-Hastings for Inference in
  Nonlinear Dynamical Models》
---
作者：
Johan Dahlin and Thomas B. Sch\\\"on
---
最新提交年份：
2019
---
英文摘要：
  This tutorial provides a gentle introduction to the particle Metropolis-Hastings (PMH) algorithm for parameter inference in nonlinear state-space models together with a software implementation in the statistical programming language R. We employ a step-by-step approach to develop an implementation of the PMH algorithm (and the particle filter within) together with the reader. This final implementation is also available as the package pmhtutorial in the CRAN repository. Throughout the tutorial, we provide some intuition as to how the algorithm operates and discuss some solutions to problems that might occur in practice. To illustrate the use of PMH, we consider parameter inference in a linear Gaussian state-space model with synthetic data and a nonlinear stochastic volatility model with real-world data.
---
中文摘要：
本教程温和地介绍了用于非线性状态空间模型中参数推断的粒子Metropolis Hastings（PMH）算法，以及统计编程语言R中的软件实现。我们采用逐步方法与读者一起开发PMH算法（以及其中的粒子过滤器）的实现。这个最终的实现也可以作为包pmhtutorial在CRAN存储库中获得。在整个教程中，我们提供了一些关于算法如何运行的直觉，并讨论了一些解决实践中可能出现的问题的方法。为了说明PMH的使用，我们考虑了具有合成数据的线性高斯状态空间模型和具有真实数据的非线性随机波动率模型中的参数推断。
---
分类信息：

一级分类：Statistics        统计学
二级分类：Computation        计算
分类描述：Algorithms, Simulation, Visualization
算法、模拟、可视化
--
一级分类：Quantitative Finance        数量金融学
二级分类：Statistical Finance        统计金融
分类描述：Statistical, econometric and econophysics analyses with applications to financial markets and economic data
统计、计量经济学和经济物理学分析及其在金融市场和经济数据中的应用
--
一级分类：Statistics        统计学
二级分类：Machine Learning        机器学习
分类描述：Covers machine learning papers (supervised, unsupervised, semi-supervised learning, graphical models, reinforcement learning, bandits, high dimensional inference, etc.) with a statistical or theoretical grounding
覆盖机器学习论文（监督，无监督，半监督学习，图形模型，强化学习，强盗，高维推理等）与统计或理论基础
--

---
PDF下载：
--> Getting_Started_with_Particle_Metropolis-Hastings_for_Inference_in_Nonlinear_Dyn.pdf (1.11 MB)

2022-5-9 11:09:19 上传

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏0 回帖

关键词：Metropolis Polis Metro ASTIN Etro

非线性动力学模型中的粒子Metropolis Hastings推断入门约翰·达林林林克"oping大学托马斯·B·舍努帕拉大学摘要本教程温和地介绍了用于非线性状态空间模型中参数推断的粒子Metropolis Hastings（PMH）算法，以及统计编程语言中的软件实现R.我们采用逐步的方法与读者一起开发PMH算法（以及粒子过滤器）的实现。该最终实现也可以从综合R归档网络（CRAN）存储库中的PackagePMHTurial中获得。通过本教程，我们提供了一些关于算法如何运行的直观信息，并讨论了一些解决实践中可能出现的问题的方法。为了说明PMH的使用，我们考虑了具有合成数据的线性高斯状态空间模型和具有真实数据的非线性随机波动率模型中的参数推断。关键词：贝叶斯推理、状态空间模型、粒子滤波、粒子马尔可夫链蒙特卡罗、随机波动率模型。1.引言我们关注非线性状态空间模型（SSMs）中的贝叶斯参数推断。这是一个重要的问题，因为SSM在自动控制（Ljung 1999）、计量经济学（Durbin and Koopman 2012）、金融（Tsay 2005）和许多其他领域无处不在。Langrock（2011）对状态空间建模的一些具体应用进行了概述。在SSMs中，贝叶斯参数和状态推理的主要问题是，它是一个难以分析的问题，无法以封闭形式解决。因此，需要近似值，这些值通常分为两类：分析近似值和统计模拟。

本教程的重点是后一个类别的成员，即粒子大都会黑斯廷斯（PMH；Andrieu等人，2010），该类别通过使用特定的马尔可夫链对其进行采样来近似后验分布。这要求我们能够逐点计算后验概率的无偏估计。在PMH中，这些由粒子过滤器提供（Gordon等人1993；Doucet and Johansen 2011），这是另一个重要且相当普遍的统计模拟算法。本教程的主要目的和贡献是对PMH算法进行温和的介绍（包括方法论和软件方面）。我们假设读者熟悉Brockwell和Davis（2002）的传统时间序列分析和卡尔曼滤波。在Ross（2012）的水平上熟悉蒙特卡罗近似、重要性抽样和马尔科夫链蒙特卡罗（MCMC）也是有益的。2非线性动力学模型推理的PMH入门重点是在编程语言R（R Core Team 2018）中逐步实现粒子过滤器和PMH算法。我们使用识字编程来使用自上而下的方法构建代码，参见Pharr和Humphreys（2010）中的第1.1节。最终实施以R包pmhtutorial（Dahlin 2018）的形式提供，通过综合R档案网络（CRAN）存储库根据GPL-2许可证分发，可在https://CRAN.R-project.org/package=pmhtutorial.R源代码还通过GitHub提供了相应版本的MATLAB（MathWorks Inc.2017）和Python（van Rossum et al.2011），以便读者能够从许多不同的背景中学习，详情请参见第8节。有许多现有的教程与此相关。Arulampalam等人对颗粒过滤进行了全面介绍。

（2002），其中讨论了许多不同的实现，并提供了伪代码以促进实现。对于对粒子过滤特别感兴趣的读者来说，这是一个极好的开始。然而，他们既不讨论也不实施PMH。软件LibBi和Murray（2015）的相关教程包括高性能的实现以及粒子滤波和PMH的示例。本教程的重点有点不同，因为它试图解释算法的所有步骤，并提供几种不同编程语言的基本实现。最后，Ala Luhtala等人（2016年）介绍了一篇关于扭曲粒子过滤的高级教程论文。它还包括PMH的伪代码，但其重点是粒子过滤，没有讨论，例如，如何在PMH中调整提案。在本教程的第3.3、4.3、6和7节中，我们将继续讨论相关的工作和软件。本教程的配置如下所示。在第2节中，我们介绍了SSMs中的贝叶斯参数和状态推理问题，并讨论了PMH算法的工作方式和原因。在第3章和第4章中，我们对一个玩具模型实施了粒子滤波和PMH算法，研究了它们的特性，并给出了一个适合进一步研究的小调查结果。在第5和第6节中，PMH的实施适用于从财务角度解决现实问题。此外，还讨论了PMH的一些更高级的修改和改进，并举例说明。最后，我们在第7节讨论了相关的软件实现，以此结束本教程。2.PMH算法概述在本节中，我们将介绍SSMs中与参数和潜在状态相关的贝叶斯推理问题。然后引入PMH算法来解决这些棘手的问题。

从两个互补的角度讨论了将粒子过滤器纳入PMH的问题：作为接受概率的近似值，以及作为增强后验分布中的辅助变量。代码的框架被概括为由<variable>表示的不同变量的集合。然后使用类似C++的语法用内容填充代码的这些部分。因此，运算符=将一段代码分配给相关变量（覆盖其中已存储的内容）。另一方面，+=表示在变量中存储的当前代码之后添加一些额外代码的操作。强烈建议读者打印GitHub的源代码，并与教程并行阅读，以便快速学习材料。约翰·达林，托马斯·B·舍恩3··xt-1xtxt+1··yt-1ytyt+1··图1：以图形模型表示的SSM结构。2.1. 贝叶斯推断在SSMSI图1中，我们给出了一个SSM的图形模型，其中包含潜在状态（顶部）、观测（底部），但没有外部输入。潜在状态xt仅取决于前一状态xt-这是一个（一阶）马尔可夫过程。也就是说，过去的所有信息都是x0:t-1，{xs}t-1s=0在最新状态xt中汇总-1.观察结果y1：皮重是条件独立的，因为它们只通过状态相互关联。我们假设状态和观测值是实值，即yt∈ Y Rnyand xt∈ 十、Rnx。初始状态xis根据由一些未知静态参数θ参数化的密度μθ（x）分布∈ Θ  Rp.潜态和观测过程由密度fθ（xt | xt）决定-1） 和gθ（yt | xt）。

描述状态过程的密度给出了下一个状态为xt的概率给定了前一个状态xt-1.类似的解释适用于观察过程。有了适当的符号，一般的非线性SSM可以用asx表示~ μθ（x），xt | xt-1.~ fθ（xt | xt）-1） ，yt | xt~ gθ（yt | xt），（1），其中随机变量与其实现之间没有符号上的区别。贝叶斯参数推断的主要目标是通过计算πθ（θ），p（θ| y1:T）=p（θ）p（y1:T |θ）ZΘp（θ）p（y1:T |θ）p（y1:T |θ）dθ，γθ（θ）Zθ，（2）式中p（θ）和p（y1:T |θ），p（y1:T）分别表示参数的先验分布和数据分布。由γθ（θ）=p（θ）p（y1:T |θ）表示的非归一化后验分布是PMH构造中的一个重要组成部分。分母Zθ通常被称为边际似然或模型证据。可以直接修改针对向量值量提出的算法，也可以包括aknown外部输入ut。在本教程中，我们假设xt | xt-1和yt | x可以建模为具有密度函数的连续随机变量。然而，所讨论的算法也适用于确定性状态和离散测试/观测。4.非线性动力学模型推断的PMH入门由于状态序列未知，对于大多数有趣的SSM，数据分布pθ（y1:T）是难以处理的。通过研究分解pθ（y1:T）=pθ（y）TYt=2pθ（yt | y1:T），可以看出这个问题-1） ，（3）其中所谓的预测似然由pθ（yt | y1:t）给出-1） =Zgθ（yt | xt）fθ（xt | xt）-1） πt-1（xt）-1） dxtdxt-1.

（4） 边际滤波分布πt-1（xt）-1） ，pθ（xt）-1 | y1:t-1） 可以从贝叶斯滤波递归（Anderson和Moore 2005）中获得，对于固定θ，它由πt（xt）=gθ（yt | xt）pθ（yt | y1:t）给出-1） ZXfθ（xt | xt）-1） πt-1（xt）-1） dxt-1，对于0<t≤ T、 （5）以π（x）=μθ（x）作为初始化。理论上，通过迭代应用（5），可以构造一系列过滤分布。然而，在实践中，由于xt的边缘化，这一策略不能用于许多感兴趣的模型-1（用积分表示）不能以封闭形式进行。2.2. 构建马尔可夫链PMH算法通过利用统计模拟，为这两个棘手问题提供了一个优雅的解决方案。首先，使用粒子过滤器来近似πt-1（xt）-1） 并获得对可能性的逐点无偏估计。其次，基于粒子滤波器提供的似然估计，采用了一种称为Metropolis Hastings（MH）的MCMC算法（Robert and Casella 2004）来近似πθ（θ）。PMH是这两种算法的结合。MCMC方法通过构造特定的马尔可夫链从所谓的目标分布生成样本，马尔可夫链可用于形成经验近似。在本教程中，目标分布将是参数后验πθ（θ）或参数和状态的后验πθ，T（θ，x0:T）=p（θ，x0:T | y1:T）。我们在这里集中讨论前一种情况，后一种情况也类似。执行PMH算法会产生K个相关样本θ（1:K），{θ（K）}Kk=1，可用于构造πθ（θ）的蒙特卡罗近似。参数后验分布的经验近似由bπKθ（dθ）=KKXk=1Δθ（K）（dθ），（6）给出，它对应于一组位于θ=θ且权重相等的狄拉克δ函数Δθ（dθ）。

在实践中，使用直方图或核密度估计器来可视化从（6）中获得的参数后验估计。PMH中马尔可夫链的构造相当于执行两个步骤，从πθ（θ）中获得一个样本。第一步是从提议分布q（θ|θ（k））中提出所谓的候选参数θ-1） ）给定由θ（k）表示的马尔可夫链的先前状态-1). 用户可以自由选择此方案，但其支持范围应包括Johan Dahlin，Thomas B.Sch"on 5对目标发行版的支持。第二步是确定Markovchain的状态是否应该更改为候选参数θ，或者它是否应该保持在先前的状态θ（k）-1). 这个决定是随机的，候选参数被指定为马尔可夫链的下一个状态，即{θ（k）← θ} 具有一定的接受概率α（θ，θ（k-1） ）由α给出θ、 θ（k-1)= min（1，πθ）θπθθ（k）-1)Qθ（k）-1)θQθθ（k）-1))= min（1，γθ）θγθθ（k）-1)Qθ（k）-1)θQθθ（k）-1)), （7） 其中Zθ取消，因为它与马尔可夫链的当前状态无关。（7）中引入的随机行为有助于探索（理论上）整个后验概率，也是马尔可夫链实际将目标作为其平稳分布的必要条件。后者被称为详细平衡，它确保马尔可夫链是可逆的，并且具有正确的平稳分布。也就是说，确保PMH生成的样本实际上来自πθ（θ）。因此，随机接受决策可以看作是对提议分布生成的马尔可夫链的修正。修正马尔科夫链的平稳分布就是期望的目标分布。

在某种意义上，这类似于重要性抽样，其中来自提案分布的样本通过加权进行校正，以根据目标分布进行近似分布。接受概率（7）的直觉（忽略提议q的影响）是，如果πθ，我们总是接受候选参数θθ> πθθ（k）-1). 也就是说，与之前的状态θ（k）相比，θ增加了目标的值-1). 这导致了amode搜索行为，类似于估计后验分布最大值的优化算法。然而，从（7）中，我们也注意到，可以接受后验值的小幅降低，以便于整个后验的探查。这是PMH和优化算法之间的主要区别，前者侧重于映射整个后部，而后者只希望找到其模式的位置。2.3. 近似测试函数在贝叶斯参数推理中，人们的兴趣通常在于计算称为测试函数的aso的期望值，这是一个性能良好的（可积）函数→ Rn~n将参数映射到实空间上的值。关于参数后验概率的期望值由πθ[~n]、Eh~n（θ）给出y1:Ti=ZΘ~n（θ）πθ（θ）dθ，（8）其中，例如，选择Θ（θ）=θ对应于计算后验参数的平均值。（8）中的期望值可以通过（6）中的经验近似估计，根据bπKθ[~n]、ZΘΘ（θ）bπKθ（dθ）=KKXk=1ZΘΘΘ（θ）Δθ（K）（dθ）=KKXk=1Θθ（k）, （9） 最后一个等式来自狄拉克三角函数的性质。该估计器性能良好，可以建立一个大数定律（LLN）和一个中心极限6。从PMH开始，在非线性动力学模型（CLT）中进行推理，有关更多信息，请参见Meyn和Tweedie（2009）。

从LLN中，我们知道估计量是一致的（且渐近无偏为K）→ ∞). 此外，从CLT中，我们知道误差近似为高斯分布，方差减少1/K，这是通常的蒙特卡罗速率。请注意，LLN通常假设独立样本，但称为遍历定理的类似结果给出了类似的结果（在某些假设下），即使θ（1:K）是相关的。然而，渐近方差通常比样本不相关时大。2.4. 在PMH算法得以实现之前，还有一个主要障碍有待解决。接受概率（7）难以计算，因为可能性未知。综上所述，我们知道粒子滤波器可以为似然度提供无偏的逐点估计，因此也可以提供后验πθ（θ）。事实证明，无偏性对于PMH能够提供πθ（θ）的有效经验近似是至关重要的。这种方法被称为精确近似，因为似然被近似所取代，但算法仍然保持其有效性，更多信息请参见Andrieu和Roberts（2009）。粒子过滤器从πt（xt）为每个t生成一组样本，可用于创建经验近似值。这些样本是通过依次应用重要性抽样来近似（5）的解而生成的。然后，过滤分布的近似值由bπNt（dxt）、bpNθ（dxt | y1:t）=NXi=1v（i）tPNj=1v（j）t |{z}，w（i）tδx（i）t（dxt），（10）给出，其中粒子x（i）及其归一化权重w（i）构成在粒子过滤器运行期间生成的所谓粒子系统。未规范化的权重用v（i）t表示。