跳到主要内容
研究论文
开放式访问

hIPPYlib-MUQ:一种贝叶斯推理软件框架,用于不确定性下数据与复杂预测模型的集成

出版:2023年6月15日出版历史

跳过抽象节

摘要

贝叶斯推理为数据与数学模型的集成提供了一个系统框架,以量化反问题解决方案中的不确定性。然而,由以下描述的复杂正向模型控制的贝叶斯逆问题的解偏微分方程仍然禁止使用黑盒马尔可夫链蒙特卡罗(MCMC)方法。我们提出了hIPPYlib MUQ,这是一个可扩展和可扩展的软件框架,包含了最先进算法的实现,旨在克服高维、PDE约束的贝叶斯逆问题的挑战。这些算法通过导数信息和低秩近似,利用参数空间的几何特性和内在低维性,加速MCMC采样。该软件集成了两个互补的开源软件包,即hIPPYlib和MUQ。hIPPYlib使用自动生成的基于伴随的导数解决PDE约束的逆问题,但它缺乏完整的贝叶斯功能。MUQ提供了一系列强大的贝叶斯反演模型和算法,但预计正向模型将配备梯度和Hessian,以允许大规模求解。通过将这两个互补库结合起来,我们创建了一个健壮、可扩展和高效的软件框架,该框架实现了每一个库的优点,并允许我们处理跨越广泛科学和工程学科的复杂大规模贝叶斯逆问题。为了说明hIPPYlib-MUQ的功能,我们比较了集成软件中可用于几个高维贝叶斯反问题的一些MCMC方法。这些包括以线性和非线性偏微分方程、各种噪声模型和不同参数维为特征的问题。结果表明,与传统的黑盒和基于梯度的MCMC算法相比,通过利用Hessian信息(来自log-posterior)可以获得较大的(~50×)加速比,这突出了集成hIPPYlib-MUQ框架的强大功能。

跳过1简介部分

1简介

随着观测和实验数据的迅速爆炸,一个突出的挑战是如何从这些数据中获得知识和见解,以做出更好的预测和高后果的决策。这个问题出现在科学、工程、技术和医学的所有领域,在许多情况下,有可用的数学模型来表示观察或测量数据的潜在物理系统。这些模型通常受到未知或不确定输入模型参数(例如,系数场、本构关系、源项、几何、初始和/或边界条件)以及噪声和有限观测值引起的相当大的不确定性的影响。目标是通过相应的偏微分方程模型,并量化此类反问题求解中的不确定性。

贝叶斯反演为数据与基于PDE的复杂模型集成提供了一个系统框架,以量化模型参数推断中的不确定性[Kaipio和Somersalo2005; 塔兰托拉2005]。在贝叶斯框架中,将噪声数据和可能不确定的数学模型与先验信息集成在一起,得出模型参数的后验概率分布。这个马尔可夫链蒙特卡罗(MCMC)方法是利用抽样技术研究后验分布的常用方法。然而,基于传统MCMC方法的复杂正演模型的贝叶斯反演面临着一些计算挑战。首先,表征模型参数的后验分布或随后的预测通常需要对昂贵的大规模PDE模型进行重复评估。其次,由于模型参数到观测量的非线性映射,后验分布往往具有复杂的结构。第三,参数通常是字段,经过离散化后,会产生非常高的后验概率。这些困难使得具有复杂大规模PDE正演模型的贝叶斯逆问题的求解在计算上非常困难。

为了克服大规模偏微分方程控制的贝叶斯反问题的局限性,人们进行了大量的研究工作。随着高性能计算的快速发展,以及可扩展PDE求解器的进步,对不同输入参数的正向PDE模型进行了重复评估[Balay等人。2018; Trilinos项目团队2020]变得容易驾驭。此外,利用结构的MCMC方法有效地促进了复杂后验分布的探索[Cotter等人。2012; Beskos等人。2017; Petra等人。2014; Bui-Thanh等人。2012]。最后,事实证明,降维方法可以显著降低MCMC模拟的计算成本[Cui等人。2016年b; Zahm等人。2022]。应用和组合这些先进的技术可能是非常具有挑战性的。因此,一个能够帮助计算和科学界应用、扩展和调整这些方法的计算工具将是非常有益的。

在本文中,我们提出了一个软件框架,用基于PDE的正向模型来解决大规模贝叶斯逆问题,该框架在广泛的科学和工程领域具有应用。该软件集成了两个开源软件包反向问题Python库(hIPPYlib)[维拉等人。2021]以及麻省理工学院不确定性量化库(MUQ)[Parno等人。2014]尊重他们有吸引力的互补能力。

hIPPYlib是一个可扩展的软件框架,用于解决受复杂PDE模型约束的确定性和线性化贝叶斯逆问题。基于FEniCS[Logg等人。2012]对于偏微分方程的有限元近似和PETSc[Balay等人。2014]对于高性能线性代数运算和求解器,hIPPYlib允许用户以相对简单的方式描述(和求解)基于PDE的底层正向模型(逆问题求解器所需)。hIPPYlib还包含求解确定性和线性化贝叶斯反问题的有效数值方法的实现。这些包括全球化的不精确纽顿-共轭梯度[Akçelik et al。2006; 波茨和舒尔茨2012],基于伴随的梯度和Hessian作用计算[Tröltzsch2010],随机线性代数[Halko等人。2011]以及大规模高斯场的可缩放采样。hIPPYlib中实现的最新算法有效地解决了线性化贝叶斯反问题。然而,hIPPYlib主要用于确定性和线性化贝叶斯逆问题,缺乏完整的贝叶斯反演能力。

MUQ对hIPPYlib的功能进行了补充,为贝叶斯推理问题的公式化和解决提供了更多支持。MUQ是一个模块化软件框架,旨在解决涉及复杂模型的不确定性量化问题。该软件提供了一个抽象建模界面,用于组合物理(例如PDE)和统计组件(例如,加性误差模型、高斯过程先验等),以灵活和半侵入的方式定义贝叶斯后验分布。MUQ还包含一套强大的不确定性量化算法,包括马尔可夫链蒙特卡罗(MCMC)方法[Parno和Marzouk2018],交通地图【Marzouk等人。2016],类似于通知子空间,稀疏自适应广义多项式混沌扩张[康拉德和马尔佐克2013]Karhunen-Loéve展开,高斯过程建模[Rasmussen和Williams2005; Hartikainen和Särkkä2010],以及能够进行全局灵敏度分析和优化实验设计的预测方法。然而,为了有效地将这些工具应用于贝叶斯反问题,MUQ需要配备hIPPYlib可以提供的梯度类型和/或Hessian信息。

通过将这两个软件库相互连接,我们旨在创建一个健壮、可扩展、高效、灵活和易于使用的软件框架,以克服复杂大规模贝叶斯反问题固有的计算挑战。软件的代表性特点总结如下:

  • 该软件结合了hIPPYlib和MUQ这两个软件包的优点,以支持大规模PDE管理的贝叶斯反问题的可扩展解决方案。

  • 可以使用各种先进的MCMC方法来利用问题结构(例如,对数后验的导数/Hessian信息)。

  • 该软件利用对数后验函数的稀疏性、低维性和几何结构来实现可扩展和高效的MCMC方法。

  • 实施收敛诊断以评估MCMC样本的质量。

在以下几节中,我们首先简要回顾了无穷维和有限维空间中由偏微分方程控制的反问题的贝叶斯公式(第节2). 然后,我们描述了用于表征后部的MCMC方法(第)并总结软件中可用的收敛诊断(第节3.2). 接下来,我们介绍hIPPYlib-MUQ的设计(第4). 最后,我们提出了许多基准测试问题和一个分步实施指南,以说明当前软件的关键方面(第节5). 章节6提供了结束语。

跳过2贝叶斯推理框架部分

2贝叶斯推理框架

在本节中,我们将简要讨论用于求解PDE控制的反问题的贝叶斯推理方法。我们首先概述了Stuart之后的无限维贝叶斯反问题的框架[2010]; Bui Thanh等人。[2013]; Petra等人。[2014]。然后我们简要讨论了先验分布和后验分布的有限维近似;在Bui-Thanh等人。[2013]。最后,我们给出了后验分布的拉普拉斯近似,这需要解PDE约束优化问题来计算最大后验概率(地图).

2.1无限维贝叶斯反问题

反问题的目标是确定未知的输入参数字段这将产生给定的观测(或实验)噪声数据\(\boldsymbol{d}\)通过(基于物理的)数学模型。换句话说,给定\(\boldsymbol{d}\in\mathbb{R}^q\),我们试图推断\(m\in\mathcal{m}\)这样(1)\(开始{方程式}\粗体符号{d}\近似\mathcal{F}(m),结束{方程式{)哪里\(\mathcal{F}:\mathcal{M}\rightarrow\mathbb{R}^q\)是可观测参数根据给定参数预测观测值的地图通过正向数学模型,以及\(mathcal{M})是定义在域上的函数的无穷维希尔伯特空间\(\mathcal{D}\subset\mathbb{R}^D\)。请注意,此图的评估涉及求解给定的正向PDE模型然后从正问题的解中提取观测值。在下文中,我们假设前向方程残差是连续Fréchet-可微的,其Jacobian是具有连续逆的连续线性算子[Ghattas和Willcox2021].

在贝叶斯方法中,反问题被定义为一个统计推理问题。不确定参数和观测数据\(粗体符号{d})被视为随机变量,其解是条件概率分布\(mu{mbox{post}}(m|\boldsymbol{d})),表示给定数据的参数估计的置信水平。该方法结合了先验模型反映我们在获取数据之前对参数的先验知识或假设,以及似然模型表示给定参数集可能产生观测数据的概率。

使用Radon-Nikodym导数[Williams1991]后部测量\(mu{mbox{post}}(m))关于先验测度\(mu_{mbox{previous}}(m)),无限维贝叶斯定理表述为(2)\(开始{等式}\frac{d\mu{mbox{post}}}{d\mu{mbox}previous}}}\propto\pi{mbox类}}(粗体符号{d}|m),结束{方程式})哪里\(\pi_{\box{like}})表示似然函数。关于后验测量的详细条件,我们建议读者参考斯图亚特[2010].

在不损失一般性的情况下,似然的概率密度函数可以表示为(3)\(\begin{equation}\pi_{mbox{like}}(\boldsymbol{d}|m)\propto\exp\left\lbrace-\Phi(m;\boldsymbol{d})\right\rbrace)。\结束{方程式}\)负对数似然函数\(\Phi(m;\boldsymbol{d})\)有不同的形式,这取决于如何对测量不确定性和/或建模误差产生的噪声进行建模;例如,在加性高斯噪声模型的情况下\带有高斯噪声随机变量的(\boldsymbol{d}=\mathcal{F}(m)+\boldsymbol{eta}\)\均值和协方差矩阵为零的(黑体符号{\eta}\in\mathbb{R}^q\)\(\boldsymbol{\Gamma}_{\!\mbox{noise}}\ in \mathbb{R}^{q\ timesq}\),它的形式是\(\Phi(m;\boldsymbol{d})=\frac{1}{2}\Vert\mathcal{F}(m)-\boldsymbol{d\Vert^2_{\boldsimbol{\Gamma}_{\!\mbox{noise}}^{-1}}\),其中\(\Vert\cdot\Vert_{\boldsymbol{\Gamma}_{\!\mbox{noise}}^{-1}}\)表示\由逆噪声协方差加权的(L^2)范数\(\boldsymbol{\Gamma}_{\!\mbox{noise}}^{-1}\)。

我们将先验作为高斯测度,1即。,\(\mu_{\mbox{previous}}=\mathcal{N}\!({m_{\text{pr}}},{mathcal{C}_{\text{{previous}}}),并假设来自先验分布的样本是域中的平方积分函数\(\mathcal{D}\),即属于\(L^2(\mathcal{D})\)。先验协方差算子\(\mathcal{C}_{\text{{previous}}})被构造成一个追踪类算子,以保证样本的先验分布方差有界以及贝叶斯反问题的适定性[Stuart2010; Bui-Thanh等人。2013; Villa等人。2021]详细说明。具体来说,我们将先验协方差算子作为v(v)类拉普拉斯算子的次幂,即\(\mathcal{C}_{\text{{previor}}}:=\mathcal{A}^{-v}=(-\gamma\Delta+\deltaI)^{-v};v\gt\frac{d}{2})。在这里\(伽玛射线)和\(δ)控制先验算子的相关长度和逐点方差[Villa等人。2021; Lindgren等人。2011]。这些先前的选择确保了\(\mathcal{C}_{\text{{previous}})是一个追踪类算子,保证了有界逐点方差和适定的无限维贝叶斯反问题[Stuart2010; Bui-Thanh等人。2013].

2.2贝叶斯公式的离散化

在这里,我们简要介绍了上一节中描述的贝叶斯公式的有限维近似。我们考虑有限维子空间\(\mathcal{M} 小时(_h)\)第页,共页\(\mathcal{M}\),由一组基函数的跨度定义\(\lbrace\phi_j\rbrace\{j=1}^n\)。例如,在hIPPYlib-MUQ中,这些基函数是域网格离散化的每个元素上的全局连续分段多项式\(mathcal{D})[贝克等人。1981; 绞合和固定1988]。这些是椭圆算子离散化的自然选择\(\mathcal{A}\)用于定义先验协方差。参数字段然后近似为\(m\approxix m_h=\sum_{j=1}^nm_j\phi_j\),在下面,\(粗体符号{m}=left(m_1,dots,m_n\right)^T\in\mathbb{R}^n)表示\(m_h\)。

先验测度的有限维逼近\(\mu{\mbox{previous}}\)现在由密度指定(4)\(开始{方程式}\pi_{mbox{previous}}(黑体符号{m})\propto\exp\Big(-\frac{1}{2}\Vert\boldsymbol{米}-{\boldsymbol{m}}_{\mbox{pr}}\Vert^2_{\bodsymbol{\Gamma}_{\!\mbox{previous}}^{-1}}\Big),\end{方程式}\)哪里\({\boldsymbol{m}}_{\box{pr}}\in\mathbb{R}^n\)和\(\boldsymbol{\Gamma}_{\!\mbox{previous}}\in\mathbb{R}^{n\timesn}\)是离散化\(m_{text{pr}}\)和\(\mathcal{C}_{\text{{previous}}\)。我们建议读者参考Bui-Thanh等人。[2013]; Villa等人。[2021]先验协方差矩阵的显式表示\(\boldsymbol{\Gamma}_{\!\mbox{previous}}\)。

然后是后验测度的有限维近似密度的Bayes定理\(\mu_{\mbox{post}}\)由以下公式给出(5)\(\begin{equation}\pi{\mbox{post}}(\boldsymbol{m})=\pi{\nbox{post{}}。\结束{方程式}\)有限维后验概率密度函数可以明确表示为(6)\(\begin{equation}\pi_{mbox{post}}(\boldsymbol{m})\propto\exp\Big(-\Phi(\boltsymbol}m};\boldsymbol{d})-\frac{1}{2}\Vert\boldsympol{米}-{\boldsymbol{m}}_{\mbox{pr}}\Vert^2_{\bolsymbol{\Gamma}_{\!\mbox{previous}}^{-1}}\Big)。\结束{方程式}\)这里,正在评估\(\Phi(\boldsymbol{m};\boldsymbol{d})\)需要构造参数到可观测映射的有限维离散化,\({\mathbf{{F}}}(\boldsymbol{m})\)。

2.3后验分布的拉普拉斯近似

一般来说,后验概率分布(6)由于参数-观测值映射的非线性,因此不是高斯的。在本节中,我们将讨论所谓的线性化的利用拉普拉斯近似的贝叶斯反问题。拉普拉斯近似相当于构造一个以最大后验概率(MAP)点。MAP点代表后验参数向量的最可能值,即。,(7)\(\开始{方程式}\粗体符号{{米}_{\mbox{MAP}}:=\text{参数}_{\boldsymbol{m}}(-\log\pi_{\mbox{post}})=\text{参数}_{\boldsymbol{m}}\左[\Phi(\boldsymbol{m};\boldsimbol{d})+\frac{1}{2}\Vert\boldsympol{米}-{\boldsymbol{m}}_{\box{pr}}\Vert^2_{\boldsymbol{\Gamma}_{\!\mbox{previor}}^{-1}}}\right]。\结束{方程式}\)拉普拉斯近似的协方差矩阵是MAP点处评估的负对数后验的海森矩阵的倒数。然后在拉普拉斯近似下,线性化贝叶斯反问题的解由下式给出(8)\(\begin{equation}\hat{pi}_{text{post}}(\boldsymbol{m})\sim\mathcal{N}(\foldsymbol){{米}_{\box{MAP}}},\boldsymbol{\Gamma}_{\!\box{post}}),\ end{equation}\)具有(9)\(\begin{equation}\boldsymbol{\Gamma}_{\!\mbox{post}}=\mathbf{{H}}^{-1}(\boldsymbol{{米}_{\mbox{MAP}})=\左(\mathbf{{H}}_{\mbax{mismit}}(\boldsymbol{{米}_{mbox{MAP}}})+\boldsymbol{\Gamma}_{\!\mbox{previous}}^{-1}\right)^{-1{,\end{方程式}\)哪里\(\mathbf{{H}}(\boldsymbol{{米}_{\mbox{MAP}}})和\(\mathbf{H}}_{\box{misfit})(\boldsymbol{{米}_{\mbox{MAP}}})表示分别在MAP点处评估的负对数后验和负对数似然的Hessian矩阵;参见Villa等人。[2021]用于使用基于伴随的方法推导此Hessian。

后验的高斯近似值的质量取决于参数-观测值映射中的非线性程度、噪声协方差矩阵和观测值的数量[Gelman等人。2004; 塔兰托拉2005; Bui-Thanh等人。2013; Isaac等人。2015; 埃文斯和斯瓦茨2000; 按下2003; 斯蒂格勒1986; 蒂尔尼和卡丹1986; Wong(王)2001]。当参数-观测值映射为线性且加性噪声和先验模型均为高斯时,拉普拉斯近似是精确的。即使参数到观测值的映射是显著非线性的,拉普拉斯近似也是用MCMC方法实现可扩展、高效和准确后验采样的关键因素,我们将在下一节中讨论。

请注意,拉普拉斯近似涉及负对数似然的黑森函数(黑森函数的数据不匹配部分),当参数维数较大时,可能无法明确构造。然而,数据通常仅提供有关参数字段的有限信息,因此,Hessian矩阵的本征谱通常衰减得很快。我们利用Hessian的紧致性来克服其高昂的计算成本,并使用无矩阵方法(例如随机子空间迭代[Halko et al。2011]).

具体来说,我们构造了数据不匹配Hessian的低阶近似,即。,\(\mathbf{{H}}_{mbox{missit}}\近似\boldsymbol{\Gamma}_{\!\mbox{priority}}^{-1}\mathbf{{V}}_r\boldsympol{\Lambda}_r\mathbf}{V}{_r^T\boldsymbol{\ Gamma}{{\\(\boldsymbol{\Lambda}_r=\text{diag}(\Lambda_1,\ldots,\Lambda _r)\in\mathbb{r}^{r\times r}\)和\(\mathbf{{V}}_r=[\boldsymbol){v} _1个,\ldot,\boldsymbol{v_r}]\in\mathbb{r}^{n\times r}\)包含广义对称特征值问题的最大特征值和相应的特征向量(10)\(\开始{方程式}\mathbf{{H}}_{\mbox{mismit}}\boldsymbol{v} _ i=\lambda_i\boldsymbol{\Gamma}_{\!\mbox{previous}}^{-1}\boldsymbol{v} _ i; \四元i=1,\ldots,n.结束{方程式}注意特征向量\(\粗体符号{v} _ i\)相对于…是正交的\(\boldsymbol{\Gamma}_{\!\mbox{previous}}^{-1}\),即\(\粗体符号{v} _ i^T\boldsymbol{\Gamma}_{\!\mbox{previous}}^{-1}\boldsymbol{v} _j(_j)=\delta_{ij}\),其中\(delta{ij})是克罗内克三角洲。在这种行库近似下,使用Sherman-Morrison-Woodbury公式[Golub和Van Loan1996],我们得到了中Hessian的逆(9),(11)\(开始{方程式}\mathbf{{H}}^{-1}=\左(\mathbf{H}{{mbox{missit}}}+\boldsymbol{\Gamma}{\!\mbox{priority}}^}{-1}\right)_r^T+\mathcal{O}\left(\sum_{i=r+1}^n\frac{\lambda_i}{1+\lambda _i}\right),\end{方程式}\)哪里\(\mathbf{{D}}_r=\operatorname{diag}(\lambda_1/(\lampda_1+1),\dots,\lambda _r/(\lambda_r+1))\in\mathbb{r}^{r\times r}\)。我们建议读者参考Villa等人。[2021]有关随机算法的详细描述[Halko等人。2011]继续构造Hessian的低阶近似,包括相关的计算复杂性。

我们可以从中的最后一个余项中看出(11)以获得精确的低阶近似\(mathbf{{H}}^{-1}),我们必须保持特征向量对应于大于1的特征值。该近似用于有效地执行与黑森函数相关的各种操作,例如,将黑森函数的平方逆应用于向量,这是从本节讨论的高斯近似中提取样本所必需的。我们注意到,我们的方法的效率和可扩展性是基于数据不匹配Hessian的低范围。我们向读者推荐加塔斯和威尔科克斯[2021]其中,这一论点是通过模型问题提出的,在模型问题中可以分析显示低秩,在引用更复杂的问题时可以实证显示。

跳过3MMC采样部分

MCMC取样

如上所述,当参数-观测值映射为非线性时,拉普拉斯近似可能是后验的较差近似。在这种情况下,需要应用基于抽样的方法来探索完整的后路。在本节中,我们将重点介绍当前软件中可用的几种高级马尔可夫链蒙特卡罗(MCMC)方法。我们概述了MCMC方法的一般结构,并简要讨论了它们的主要特征。然后,我们提出各种诊断方法来评估MCMC仿真的收敛性。

3.1马尔可夫链蒙特卡罗

MCMC为探索后向分布提供了一个灵活的框架。它从后验分布中生成样本,可用于后验期望的蒙特卡罗近似。例如,对兴趣量的后验期望\(mathcal{G}(m))可以近似为(12)\(开始{方程式}\int\mathcal{G}(m),结束{方程式其中每个\(mi\sim\mu{mbox{post}})是后验分布的样本。

MCMC技术构建了遍历马尔可夫链,其中后验分布是链的唯一平稳分布[Robert和Casella2005]。因此,马尔可夫链的状态是后验分布的精确样本,可用于(12). 马尔可夫链是根据转移核定义的,它是一个位置相关的概率分布\(K(\cdot|\boldsymbol){m} _ i)\)过度状态\(\粗体符号{米}_{i+1}\)在给定先前状态的链中\(\粗体符号{m} _ i\)即。,\(\粗体符号{米}_{i+1}\sim K(\cdot|\boldsymbol{m} _ i)\). 注意,实际中必须使用有限长的链,因此蒙特卡罗估计量的统计精度高度依赖于过渡核有效探索参数空间的能力。

有几个适合于MCMC的构建转换内核的框架,包括众所周知的大都会黑斯廷斯(MH)规则【Metropolis等人。1953; 黑斯廷斯1970]吉布斯采样器(例如,[卡塞拉和乔治1992])、和延迟拒绝(DR)[Mira等人。2001]。MUQ提供了这些框架的实现,以及广义大都市黑斯廷斯(gMH)核[卡尔德黑德]2014]以及【Dodwell等人。2019]。大多数框架都是从一个或多个提案分发中抽取样本开始的\(q_1(\cdot|\boldsymbol{m} _ i),\,\ldots,\,q_K(\cdot|\boldsymbol{m} _ i)\)很容易从(例如,高斯)中采样,然后“校正”建议的样本,以获得准确但相关的后验样本。在MH和DR内核中,校正采取接受或拒绝所提出的点的形式。在gMH核中,校正涉及对多个建议点上的有限状态马尔可夫链进行分析采样。直观地说,提案分布捕捉到了后部的形状,或者是周围的局部形状\(mi)或全局参数空间,往往需要较少的“校正”,并产生更有效的算法。

建议书分发.让\(q(\cdot|\boldsymbol{m} _ i)\)表示由链的当前状态“参数化”的提案分布\(\粗体符号{m} _ i\). 我们要求建议分布很容易抽样,并且可以有效地评估其密度。MH规则【Metropolis等人。1953; 黑斯廷斯1970]定义转换内核\(K_{MH}(\cdot|\boldsymbol{m} _ i)\)通过两步过程:首先抽取随机样本\(\boldsymbol{m}^\prime\sim q(\cdot|\boldsymbol{m} _ i)\)从提议分发中,然后接受提议的样本\(\boldsymbol{m}^\prime\)作为链中的下一步\(\粗体符号{米}_{i+1}\)具有概率\(\alpha=\min\lbrace 1,\gamma\rbrace)其中\(\gamma=\frac{\pi{\mbox{post}}(\boldsymbol{m}^\prime)}{\pi_{\mbax{post{}}{m} _ i)}\frac{q(粗体符号{m} _ i|\boldsymbol{m}^\prime)}{q(\boldsymbol{m}^\prime |\bolsymbol{m} _ i)}\). 如果拒绝,则设置\(\粗体符号{米}_{i+1}=\粗体符号{m} _ i\). 在温和的技术条件下,提案分发(参见Roberts等人[2004]),MH规则定义了一个遍历的马尔可夫链\(mu_{mbox{post}})作为平稳分布,从而使链中的状态能够用于蒙特卡罗估值器。注意,详细的平衡条件(参见例如Owen[2013])通常用于验证马尔可夫链\(mu{mbox{post}})作为平稳分布,但仅此条件不足以保证链收敛到平稳分布。参见Roberts等人。[2004]详细讨论MH收敛和收敛速度。

虽然MH规则将为一大类提案分布生成有效的MCMC核,但提案对前一状态的依赖性,加上提案状态可能被拒绝,导致马尔可夫链中的样本间相关性。由于这些相关性(12)使用MCMC时,将比使用独立样本的经典蒙特卡罗设置中的更大。相关性较大的马尔可夫链将导致较大的估计方差。为了减少马尔可夫链中的相关性,我们寻求能够以较高的接受概率采取较大步骤的提案分布。从MH规则中的接受概率来看,当提案密度\(q(粗体符号{m}|\boldsymbol{m} _ i)\)是一个很好的近似值\(\pi_{\mbox{post}}(\boldsymbol{m})),因此\(\gamma\)接近1。

我们现在来描述hIPPYlib-MUQ中使用的特定提案分发。首先,我们从描述利用梯度和曲率信息加速有限维空间采样的常见提议机制开始。这些算法构成了图中立方体的左侧1然后,我们展示了如何将这些思想扩展到构建性能独立于网格细化(即独立于维度)的提案,从而将派生加速提案“提升”到无限维设置。这个“提升”操作改变了图中左侧的提案1与提案立方体右侧的维度相关的类比。

图1。

图1。各种MCMC方案分布相对于网格细化独立性(蓝色箭头)、梯度感知(绿色箭头)和曲率感知(红色箭头)的关系。缩写代表以下MCMC提案:RW代表随机行走,pCN代表预处理Crank-Nicolson,MALA代表大都会调整Langevin算法,H-pCN表示曲率信息pCN,H-MALA代表曲率信息MALA,\(infty)-MALA代表无限维MALA,以及H-\(inffy)-MALA表示曲率无限维MALA。

利用梯度和曲率信息。也许最简单、最常见但通常不高效的提案分配形式是以链中当前状态为中心的高斯分布,(13)\(开始{方程式}q_{text{RW}}(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\left(粗体符号{m} _ i,\boldsymbol{\Gamma}_{\text{prop}}\right),\end{equation}\)哪里\(\boldsymbol{\Gamma}_{\text{prop}}\in\mathbb{R}^{n\timesn}\)是用户定义的协方差矩阵。当与MH规则一起使用时随机游走(RW)该提案产生了一种MCMC算法,通常称为随机行走大都会算法。这个自适应大都市(AM)算法采用了该建议的一种变体,其中协方差\(\boldsymbol{\Gamma}_{\text{prop}}\)根据之前的示例进行了改编[Haario等人。2001]。提案协方差\匹配后验协方差的(\boldsymbol{\Gamma}{\text{prop}})可以提高效率,但随机游走建议仍然是后验密度的较差近似值。

通过对朗之万随机微分方程进行一步Euler-Maruyama离散化,可以获得一个更有效的建议[Roberts和Stramer2003]。最终的朗之万提案采用以下形式(14)\(开始{方程式}q_{text{MALA}}(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\left(粗体符号{m} _ i+\tau\boldsymbol{\Gamma}_{\text{prop}}\nabla\log\pi_{mbox{post}}(\boldsymbol{m} _ i),\,2\tau\boldsymbol{\Gamma}_{\text{prop}}\right),\end{方程式}\)哪里\(\tau)是步长参数。有此建议的MH采样器称为都市调整的朗之万算法(马来西亚)与AM算法一样,调整MALA方案的协方差也可以提高性能[Attchadé2006; 马歇尔和罗伯茨2012]。通常也会使用后验协方差的近似值,例如对数后验Hessian的倒数,以帮助MALA提案捕获后验相关性。例如,在这项工作中,我们在MAP点采用对数后验Hessian的低秩近似(c.f.方程(11))(15)\(开始{方程式}q_{text{H-MALA}}(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\left(粗体符号{m} _ i+\tau\mathbf{{H}}^{-1}\nabla\log\pi_{mbox{post}}(粗体符号{m} _ i),2\tau\mathbf{{H}}^{-1}\right)。\结束{方程式}\)该指标与Martin等人使用的指标类似。[2012]并等同于中的预处理MALA提案(14)使用拉普拉斯近似的协方差(9).

两者(13)和(14)在参数空间中使用常量协方差。允许此协方差适应后验密度的局部相关结构,可以获得更高阶的近似值,从而产生更高效的MCMC算法。在Girolma和Calderhead[2011]利用微分几何的观点定义了黎曼流形上的一系列建议机制。调整MALA提案(14)这种流形设置并忽略流形的曲率,会导致(16)\(开始{方程式}q_{text{sMMALA}}(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\left(粗体符号{m} _ i+\tau\mathbf{G}^{-1}(粗体符号{m} _ i)\nabla\log\pi_{\mbox{post}}(\boldsymbol{m} _ i),2\tau\mathbf{G}^{-1}(粗体符号{m} _ i)\右),\结束{方程式}\)哪里\(\mathbf{G}(\boldsymbol{m})\)是一个位置相关的度量张量。这被称为简化歧管MALA(sMMALA)建议。吉洛米和卡尔德黑德[2011]定义了度量张量\(\mathbf{G}(\boldsymbol{m})\)使用预期的Fisher信息度量,它提供了点处后验协方差的正定近似值\(\粗体符号{m}\)。

哈密顿蒙特卡罗方法,包括无U形旋转取样器(螺母)[霍夫曼和盖尔曼2014]是MCMC的另一类重要提案。这些技术近似地解决了哈密顿系统在参数空间中发生大跳跃的问题。虽然在许多情况下都很有效(参见Neal[2010])特别是对于纯统计模型,我们发现求解哈密顿系统涉及对我们感兴趣的基于PDE的问题进行难以处理的后验梯度评估。Parno和Marzouk的运输图MCMC算法[2018]由于构建高维变换的挑战,这里也没有考虑。

维度独立提案分发对于有限维参数,上面定义的随机游走和MALA建议可以与MCMC的MH规则一起使用。然而,它们的性能并不是离散不变的。作为函数的离散化对有限维后验数据的采样器性能进行了改进\(\pi{\mbox{post}}(\boldsymbol{m}))将恶化。随着维数的增加,MCMC跃迁核的最大两个特征值之间的差异(即谱间隙)趋于零,链的混合时间无限增长;参见Hairer等人。[2014]; Cotter等人。[2012]了解详细信息。为了获得“尺寸相关”性能,有必要对提案进行一些修改。科特等人的作品。[2012],Beskos等人。[2017]和Bardsley等人。[2020]例如,修改现有的有限维建议,以确保算法性能独立于网格细化。

RW提案的尺寸相关模拟为预处理曲柄尼科尔森(pCN)Cotter等人提出的建议。[2012]。它采取的形式(17)\(\开始{方程式}q_{text{pCN}}(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\左({\boldsymbol{m}}_{\mbox{pr}}+\sqrt{1-\beta^2}(\boldsymbol{m} _ i-{\boldsymbol{m}}{\mbox{pr}}),\,\beta^2\boldsympol{\Gamma}{\!\mbox{previous}}\right)。\结束{方程式}\)请注意,当\(β=1),pCN提案等于先前的分布。Cotter等人也对MALA提案进行了修改。[2012]获得无限维MALA(\(infty)-MALA)提案(18)\(开始{方程式}q_{text{MALA}}^{infty}(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\left(\sqrt{1-\beta^2}\boldsymbol{m} _ i+\beta\frac{\sqrt{h}}{2}\left({\boldsymbol{m}}_{\bmbox{pr}}-\boldsymbol{\Gamma}_{\!\bmbox{previor}}\nabla\Phi(\boldsymbol{m} _ i)\右),\,\β^2\粗体符号{\Gamma}_{\!\mbox{previous}}\右),\end{equation}\)哪里\(β=4\sqrt{h}/(4+h))和小时是一个可以调整的参数。而pCN和\(infty)-MALA提案导致离散化非变Metropolis-Hastings算法,它们与有限维RW和MALA类似物存在相同的缺陷,即它们没有捕捉到后验几何。

已经做出了一些努力来最大限度地减少这种不足,例如参见Beskos等人。[2017]; 鲁道夫和斯普伦克[2018]; Pinski等人。[2015]; Petra等人。[2014]。我们考虑了Pinski等人。[2015]。它将MAP点和该点的后曲率信息合并到pCN提案中,该提案由H-pCN表示,形式如下(19)\(开始{方程式}q_{text{H-pCN}}(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\left(粗体符号{{米}_{\box{MAP}}+\sqrt{1-\beta^2}(\boldsymbol{m} _ i-\粗体符号{{米}_{mbox{MAP}}}),\,\beta^2\mathbf{{H}}^{-1}\right)。\结束{方程式}\)另一种可以利用后向几何的方法是\(infty)-Beskos等人讨论的MALA提案。[2017]:(20)\(开始{方程式}q_{text{sMMALA}}^{infty}(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\left(\mu^\prime(\boldsymbol{m} _ i),\,\伽马^\素数(\粗体符号{m} _ i)\右),\结束{方程式}\)哪里(21)\(开始{eqnarray}\mu^\prime(\boldsymbol{m} _ i)&=&\sqrt{1-\beta^2}\boldsymbol{m} _ i+\beta\frac{\sqrt{h}}{2}\左(\boldsymbol{m} _ i-\mathbf{G}^{-1}\boldsymbol{\Gamma}_{\!\mbox{previous}}^{-1}(\boldsymbol{m} _ i-{\boldsymbol{m}}_{\mbox{pr}})-\mathbf{G}^{-1}\nabla\Phi(\boldsymbol{m} _ i)\右)\结束{eqnarray}\) (22)\(\开始{eqnarray}\Gamma^\prime(\boldsymbol{m} _ i)&=&\beta^2\mathbf{G}^{-1}(\boldsymbol{m} _ i). \结束{eqnarray}\)这个\(\infty)-sMMALA提案简化为\(infty\)-当\(\mathbf{G}^{-1}(\boldsymbol{m} _ i)=\boldsymbol{\Gamma}_{\!\mbox{previous}}\)。什么时候?\(\mathbf{G}(\boldsymbol{m} _ i)\)是拉普拉斯近似值Hessian来自(9),的\(\infty)-sMMALA提案简化为(23)\(begin{equation}q_{\text{H-MALA}}^\infty(\boldsymbol{m}|\boldsymbol{m} _ i)=\mathcal{N}\left(\sqrt{1-\beta^2}\boldsymbol{m} _ i+\beta\frac{\sqrt{h}}{2}\左(\boldsymbol{m} _ i-\mathbf{{H}}^{-1}\boldsymbol{\Gamma}_{\!\mbox{previous}}^}{-1}(\boldsymbol{m} _ i-{\boldsymbol{m}}_{\mbox{pr}})-\mathbf{{H}}^{-1}\nabla\Phi(\boldsymbol{m} _ i)\右),\β^2\mathbf{{H}}^{-1}\右)、\结束{方程式}\)我们用H表示-\(\infty\)-马来亚。

替代过渡内核。上述建议分发是在Metropolis-Hasting内核的上下文中进行经典考虑的。然而,有一些替代的转换核也会导致遍历马尔可夫链。这里我们考虑由Mira等人的延迟拒绝方法构造的过渡核。[2001]以及Metropolis-within-Gibbs核函数,它在不同条件下对后验分布切片重复使用Metropolin-Hastings规则来构造马尔可夫链。特别是,我们认为尺寸相关信息(DILI)方法[Cui等人。2016年a; 崔和扎姆2021],它定义了Metropolis-within-Gibbs采样器,该采样器从适当的维度依赖建议继承维度依赖属性。这里的维数依赖性意味着当问题的维数增加时,接受率和混合特性不会恶化。

延迟拒绝核允许在马尔可夫链的每个步骤中尝试多个建议。当使用具有互补属性的多个提案时,这可能是有利的。例如,可以从一个提案开始,该提案试图在参数空间中进行大幅度的跳跃,但接受概率可能较低,而回退到一个更保守的提案,该提案采取较小的步骤,接受概率较高。同样,可以从计算效率更高(例如,不需要梯度信息)但不太可能被接受的提案开始,同时在第二阶段使用更昂贵的提案机制,以确保链探索空间。在这两种情况下,如果第一个提议的行动被大都会黑斯廷斯规则拒绝,则可以在调整后的接受概率下尝试另一个更可能被接受的更昂贵的提议。也可以采用两个以上的阶段。

DILI将参数空间划分为一个有限维子空间,该子空间可以用标准提议机制进行探索,而一个互补的无限维空间可以用维度相关方法进行探索,如上述方法。生成的转换核比Metropolis-Hastings规则更复杂,但继承了互补空间方案的维度相关属性。利用广义特征值问题计算似然信息子空间(10). 如果特征值大于1,则表明似然函数在该方向上支配先验密度。因此,用于近似后验海森函数的相同低秩结构可用于将参数空间分解为似然信息子空间(LIS)由的列跨越\(\mathbf{{V}}_r)和一个正交互补空间在每个子空间中,使用标准的Metropolis-Hastings核。只要CS中的内核使用维度依赖建议(通常是pCN),那么DILI采样器将保持维度依赖。与Cui等人所述的原始实现不同。[2016年a],MUQ实现不使用白化变换,从而避免计算先验协方差的任何对称分解。一般来说,黑森语用于(10)可以适应捕获更多的相关结构。然而,我们在下面的数值实验中没有发现这一必要。

汇编MCMC算法。几乎可以将上述任何建议和内核结合起来,从而产生无数可能的MCMC算法。如图所示2MCMC算法有三个基本构建块。该链跟踪以前的点,并允许计算蒙特卡洛估计值。内核定义了一种对下一个状态进行采样的机制\(\粗体符号{米}_{i+1}\)给定当前状态的值\(\粗体符号{m} _ i\)以及一个或多个提案分发。该提案定义了一个位置特定的概率分布,该概率分布可以很容易地采样,并且具有可以有效评估的密度。我们在软件设计中模拟这些抽象接口,以定义和测试大量内核-应用程序组合。

图2。

图2。hIPPYlib-MUQ的灵活框架允许使用多种不同的转换内核和建议分发组合。请注意,每个内核都可以与任何提案分发进行交互,这使得许多不同的MCMC算法可以从相同的基本组件构建。

3.2MCMC诊断

分析长度时自然会出现两个问题马尔可夫链\([\mathbf{m} _1个,\ldot,\mathbf{米}_{一} ]\)由MCMC生产。首先,链是否收敛到平稳分布?第二,链的统计效率是什么,也就是说,链中有多少独立样本实际上有助于蒙特卡罗估值器的准确性?大多数理论保证是渐近的,当使用有限长MCMC链时,定量地回答这些问题很重要。基于这些考虑,本节描述了hIPPYlib-MUQ中实现的诊断,以检查高维MCMC链的收敛性和统计效率。

3.2.1评估融合。

为了评估收敛性,我们计算了两种不同的后验协方差的渐近无偏估计量:一种是对有限协方差的高估这是对有限的低估当这两个估计值的比率接近1时,我们可以确信MCMC链已经收敛(例如,见Brooks和Gelman[1998]; Gelman等人。[2004]; Vehtari等人。[2020]).

估算基于运行J型从随机选择的点开始的独立链,这些点比后验分布更分散\(\mu{\mbox{post}}),其中我们将“分散”分布定义为协方差大于\(\mu_{\mbox{post}}\)。每条链条的长度相同.

出租\(\mathbf{米}_{ij}\)成为链中第个MCMC样本j个,我们定义了序列内协方差矩阵\(mathbf{W})和序列间协方差矩阵\(\mathbf{B}\)作为(24)\(开始{align}\mathbf{W}&=\frac{1}{J(I-1)}\sum_{J=1}^J\sum_{I=1}^I(\mathbf{米}_{ij}-\bar{\mathbf{m}}{.j})(\mathbf{米}_{ij}-\bar{\mathbf{m}}_{.j})^T;\四边形&&\bar{\mathbf{m}}_{.j}&&=\frac{1}{I}\sum_{I=1}^I\mathbf{米}_{ij},\结束{align}\) (25)\(开始{align}\mathbf{B}&=\frac{I}{J-1}\sum_{J=1}^J(\bar{\mathbf{m}}_{.J}-\bar{\ mathbf}m}}{..})&&\bar{\mathbf{m}}{..}&&=\frac{1}{J}\sum{J=1}^J\bar{\mathbf{m}}{.J}。\结束{align}\)正如布鲁克斯和盖尔曼所指出的那样[1998],\(\mathbf{W}\)和\(\mathbf{B}\)可以组合生成一个估计值\采用以下形式的后验协方差的(\widehat{\mathbf{V}})(26)\(开始{方程式}\widehat{\mathbf{V}}=\frac{I-1}{I}\mathbf{W}+\frac}J+1}{JI}\mathbf{B}.\end{方程式{)每条链中初始点的过度分散导致\(\widehat{\mathbf{V}})高估有限的后验协方差另一方面,链内协方差的平均值\(mathbf{W})会低估协方差,因为链没有探索整个参数空间。比较\(\mathbf{W}\)和\因此,(widehat{mathbf{V}})提供了一种评估收敛性的方法。

这个\Gelman等人的统计。[2004]Vehtari等人。[2020]是一种常用的比较方式\(\mathbf{W}\)和\(\widehat{\mathbf{V}})。它使用的是\(\widehat{\mathbf{V}}\)和\(\mathbf{W}\)构造组件收敛诊断。对于高维问题,考虑多元收敛诊断更为自然。因此,我们将使用多元潜在标度缩减因子布鲁克斯和盖尔曼[1998],这是组件的自然扩展\({R}\)统计。MPSRF定义如下(27)\(\begin{方程式}\begin{aligned}\text{MPSRF}&=\sqrt{\max_a\frac{a^T\widehat{\mathbf{V}}a}{a^T \mathbf{W}a}}&=\sqrt{\ frac{I-1}{I}+\frac{J+1}{JI}\lambda_{\text{max}}}},\end{aligned}\end{方程式{)哪里\(lambda{text{max}})是满足广义特征值问题的最大特征值\(\mathbf{B}\boldsymbol{v}=\lambda\mathbf{W}\bodsymbol{v}\)。

请注意\(\text{MPSRF}\ge 1\)当\(lambda{text{max}}/gt1),当链的起始点过于分散时会发生这种情况,从而导致链间方差\(\mathbf{B}\)大于链内方差\(\mathbf{W}\)。当MPSRF接近1时,每个序列内的方差接近序列间的方差,从而表明每个链已收敛到目标分布。文献中包含了几个关于MPSRF值的建议,这些值表明收敛;例如,盖尔曼和鲁宾[1992]建议常用值\(文本{MPSRF}),而Vehtari等人。[2020]主张更保守的门槛\(\text{MPSRF})。

3.2.2统计效率。

MCMC链中的样本通常是相关的,这增加了用MCMC样本构造的蒙特卡罗估计量的方差。对于一定数量的利息\(\mathcal{G}(\boldsymbol{m})\)有效样本量(ESS)马尔可夫链的概率定义为估计所需的后验独立样本数\(\mathbb{E}[\mathcal{G}]\)具有与马尔可夫链估计相同的统计精度。因此,ESS是MCMC链中包含多少信息的度量。在这项工作中,通常假设ESS是为后验平均值的估计量推导的,即。,\(\mathbb{E}[\mathcal{G}]=\mathbb{E}[\boldsymbol{m}]\),这里我们在这个常见假设下导出ESS。其他ESS变体,如Vehtari等人所述。[2020],更适合于涉及尾部概率的问题,但这些方法在hIPPYlib-MUQ中的实现留待以后的工作。

有几种估算ESS的方法。例如,光谱方法使用MCMC链的积分自相关来估计有效样本量(例如,见Gelman等人[2004]; 沃尔夫等人。[2004]). 其他常用方法使用小样本批次的统计数据(如Flegal和Jones[2010]; Vats等人。[2019]). MUQ提供光谱和批处理方法的实现。然而,这里我们将重点放在ESS的光谱公式上,因为它为MCMC链的结构提供了额外的见解。组件的ESSk个属于\(\boldsymbol{m}\)由定义(28)\(\开始{方程式}\text{ESS}_k(_k)=\frac{JI}{1+2\sum{t=1}^\infty\rho{kt}},\end{方程式}哪里\(\rho{kt})是分量的自相关函数k个在MCMC链中处于滞后状态t吨这里是自相关函数\(rho{kt})通过以下公式进行估算【Gelman等人。2004]:(29)\(开始{方程式}\rho{kt}\approx\hat{rho}{kt}=1-\frac{v{kt}}{2\hat{垂直}_{kk}},结束{方程式})哪里\(\帽子{垂直}_{kk}\)是k个中定义的后验协方差估计的第个对角分量(26)和\(v{kt}\)是由(30)\(开始{方程式}v{kt}=\frac{1}{J(I-t)}\sum{J=1}^J\sum_{I=t+1}^I(m_{ij,k}-m_{(I-t)J,k{)^2。\结束{方程式}\)实际上,\({\rho}{kt})对于大值t吨我们截断了总和(28)有点滞后\(t^{\prime}\)。按照惯例,我们选择\(t^{prime}\ge0)为和连续自相关估计的滞后\(hat{\rho}_{2t^{\prime}}+hat{\rho}_{2t^{\prime}+1})为负[Gelman等人。2004].

跳过4软件框架部分

4软件框架

hIPPYlib-MUQ是一个Python接口,它将两个开源软件库hIPPYlib和MUQ集成到一个独特的软件框架中,使用户可以无缝地实现最先进的贝叶斯反演算法。我们在图中概述hIPPYlib和MUQ的主要功能及其互补组件的集成。

图3。

图3。hIPPYlib和MUQ的功能及其接口的描述。橙色和红色方框分别表示hIPPYlib和MUQ功能。绿色方框表示提供有限元离散化和求解器并行实现的外部软件库FEniCS和PETSc。箭头表示单向或双向交互。

hIPPYlib基于FEniCS和PETSc构建,用于PDE的离散化和求解,提供了可扩展的基于伴随的算法的Python实现,用于解决由PDE控制的大规模确定性和线性化贝叶斯逆问题。hIPPYlib模型组件提供了一组库,用户可以通过这些库描述FEniCS形式语言中的正向PDE、先前模型和似然模型[Logg等人。2012]。hIPPYlib算法组件包含优化算法、随机线性代数和高斯场的可扩展采样,这些都是有效解决确定性和线性化贝叶斯反问题所需的。hiPPYlib用户手册可在Villa等人。[2020],其中包括有关软件安装、文档和教程的详细信息。

MUQ是一个易于使用的软件框架,用于定义和解决不确定性量化问题。MUQ建模工具允许用户以一种半侵入式的方式轻松灵活地构建复杂模型,包括贝叶斯层次模型,从而实现高效的梯度和Hessian评估。MUQ还实施了多种先进的不确定性量化技术,包括MCMC采样方法、替代模型(例如高斯过程)和预测工具(例如全局敏感性分析)。我们建议读者参考muq.mit.edu有关软件的详细说明以及安装说明和分步教程。

我们注意到,通过利用这两个软件库的互补性,可以获得显著的协同效应:hIPPYlib输出,如梯度和Hessian评估和Laplace近似,以及MUQ的高级MCMC采样模块和灵活的建模能力。在集成框架中,hIPPYlib用于定义正向模型、先验和似然,以计算最大后验概率点,并基于后验协方差的近似值构建后验分布的拉普拉斯近似值,作为先验的低阶更新[Bui-Thanh等人。2013]。MUQ用于利用先进的MCMC方法来充分表征非高斯/非线性环境中的后验分布。hIPPYlib-MUQ提供了一组包装器,以MUQ可以访问hIPPYlib的各种功能的方式封装hIPPY1lib的功能。hIPPYlib-MUQ的一个关键方面是它能够使用曲率信息的MCMC方法,这对于大规模贝叶斯反问题的后验分布的有效和可扩展的探索至关重要。

在hIPPYlib-MUQ环境中,hIPPYlib提供了以下工具:(i)自动实现伴随梯度和Hessian动作,(ii)有效采样高斯马尔可夫随机场(GMRF),和(iii)用低阶黑森函数构造拉普拉斯近似。另一方面,MUQ提供了两个重要功能:(i)一个抽象的图形建模框架,它为实现模型组件(例如,先前的分布或正向模型)提供了接口,并使模型或逆向问题的多个组件易于组合和区分,以及(ii)一套MCMC算法,包括曲率变换和离散化变换方法。hIPPYlib的伴随技术使MUQ能够有效地计算具有基于PDE的组件的图形模型的梯度和Hessian动作。高效的GMRF采样和低阶拉普拉斯近似加速了MUQ的离散化变MCMC技术,该技术在MCMC提案中使用了这些hIPPYlib工具。下面提供了我们用于无缝混合这些MUQ和hIPPYlib工具的面向对象方法的详细信息。

4概述了hIPPYlib-MUQ接口实现的Python类。接口类继承自MUQ类,封装了实现曲率信息MCMC采样方法所需的hIPPYlib功能。其中包括:

图4。

图4。hIPPYlib MUQ框架的类层次结构。hIPPYlib、MUQ和接口的类分别用橙色、红色和蓝色表示。虚线箭头表示两个类之间的继承关系:箭头连接到超类,另一个连接到子类。

(1)

先验高斯界面(拉普拉斯高斯双拉普拉斯高斯)启用hIPPYlib早期模型(拉普拉斯先验BiLaplacian优先)在MUQ概率分布模型中(高斯基数). 该接口允许MUQ使用高斯先验和协方差作为v(v)-类拉普拉斯算子的次幂(\(v=1\)对于拉普拉斯高斯\(v=2\)对于BiLaplacian优先)高斯随机场的可缩放采样技术。

(2)

可能性界面(参数2对数似然)合并hIPPYlib似然函数(模型)MUQ贝叶斯模型(ModPiece公司). 该接口使MUQ模型能够有效地执行模型评估(参数-观测值映射),包括正向PDE解和基于伴随的梯度和Hessian作用计算。

(3)

拉普拉斯近似界面(LA后验高斯)让MUQ MCMC模块获得hIPPYlib生成的后验分布的拉普拉斯近似(高斯LR后). 该接口为MUQ MCMC提案提供MAP点和/或MAP点处Hessian的低阶近似值,从而实现有效且可扩展的后验抽样。

然后可以使用这些接口类,通过MUQ的图形建模接口,形成由PDE控制的贝叶斯后验模型(工作图)如图所示5以及制定MCMC提案。

图5。

图5。使用hIPPYlib-MUQ软件框架(左)和示例代码片段(右)对贝叶斯后验建模进行图形化描述。在左图中,MUQ和接口的类名分别用红色和蓝色表示。MUQ的WorkGraph类通过其成员函数AddNode和AddEdge提供了一种组合所有贝叶斯后验模型组件的方法。MUQ的IdentityOperator类标识输入参数,输入参数dim表示参数维度。MUQ的DensityProduct类定义了先验密度和似然密度的乘积,输入参数2表示输入密度的数量。函数AddEdge的第二个和第四个参数分别表示第一个参数的输出索引和第三个参数的输入索引。例如,图形。AddEdge('Parameter',0,'Prior',0)表示索引为0的“Parameter”的输出连接到索引为0“Prior”的输入。

hIPPYlib-MUQ还实现了第节所述的MCMC收敛诊断3.2这些包括潜在的尺度缩减因子及其对多元参数情况的扩展[布鲁克斯和盖尔曼1998]自相关函数和有效样本量。hIPPYlib-MUQ的所有类和功能的详细描述也可以在https://hippylib2muq.readthedocs.io/en/latest/modules.html。

跳过5数字图解部分

5数字插图

本节的目标是通过一个分步实施过程展示前面章节中讨论的集成软件框架的应用程序。我们重点比较了软件框架中可用的几种MCMC方法的性能。为了进行说明,我们首先回顾了Villa等人考虑的模型问题。[2021]二维椭圆型偏微分方程中对数扩散系数场的重构逆问题。然后我们考虑一个非线性第页-三维泊松问题,其中在Robin边界条件中推导了参数场。在本节中,我们总结了示例问题的贝叶斯公式,并给出了使用建议的软件框架hIPPYlib-MUQ 0.2.0版获得的数值结果,该软件框架建立在hIPPYlib 3.1.0版和FEniCS 2019版以及MUQ 0.3.0版的基础上;看见https://hippylib2muq.readthedocs.io/en/latest/installation.html软件安装说明及其依赖项。附带的Jupyter笔记本提供了hIPPYlib-MUQ实现的详细描述;看见https://hippylib2muq.readthedocs.io/en/latest/tutorial.html。其他示例,包括一些具有层次模型的示例,也可以在muq.mit.edu.

5.1二维泊松线性偏微分方程中的系数场反演

我们首先考虑给定逐点噪声状态测量的泊松偏微分方程中的系数场反演。我们首先描述了正向模型的设置和感兴趣的数量,然后定义了先验分布和似然分布。接下来,我们给出了后验的拉普拉斯近似,并应用几种MCMC方法来表征后验分布以及感兴趣标量的预测后验分布。最后,在网格细化研究中评估了所提方法在参数维方面的可扩展性。

5.1.1远期模型。

\(\Omega\in\mathbb{R}^d(d=2,3)\)是一个有边界的开有界域\(\partial\Omega=\partial \ Omega _D \ cup\ partial \Omega_N,\partial/Omega-D \ cap \ partial\ Omega-N=\ emptyset\)。给出了不确定参数域的一种实现,状态变量u个由…管辖(31)\(开始{align}-\nabla\cdot\left(e^m\nabla u\right)&=f&\quad&\text{in}\Omega\nonumber\nonumber,\\u&=g&&\text{on}\partial\Omega _D,\\e^m\nabla u\cdot\ mathbf{n}&=h&\text哪里(f)是卷源项,小时分别是规定的Dirichlet和Neumann边界数据,以及\(mathbf{n})是向外的单位法向量。

的弱形式(31)内容如下:查找\(u\in\mathcal{五} g(_g)\)这样的话(32)\(开始{方程式}\langle e^m\nabla u,nabla p\rangle=\langle f,p\range+\langle h,p\langgle_{\partial\Omega_N}\quad\forall p\in\mathcal{五} _0(0),\结束{方程式}\)哪里(33)\(\开始{align}\mathcal{五} g(_g)&=\left\lbrace v\ in H^1(\Omega)|v=g\\text{on}\\partial\Omega_D\right\rbrace,\nonnumber\nonumber\\\mathcal{五} _0(0)&=\left\lbrace v\在H^1(\Omega)|v=0\\text{on}\\partial\Omega_D\right\rbrace中。\结束{align}\)在上面,我们表示\(L^2 \)-内积超过\(\Omega\)作者\(\langle\cdot,\cdot\rangle)和以上\(部分\Omega _N)\(\langle\cdot,\cdot\rangle_{\partial\Omega_N}\)。

作为一个有趣的量,通过底部边界的法向通量的对数\考虑(\partial\Omega _b\subset\partial \Omega_D\)。具体来说,我们定义了利息数量\(\mathcal{G}(m)\)作为(34)\(开始{方程式}\mathcal{G}(m)=\ln\left\lbrace-\int_{\partial\Omega_b}e^m\nabla-u\cdot\mathbf{n}\,ds\right\rbrace。\结束{方程式}\)

在这个例子中,我们考虑\(\mathbb{R}^2\)没有源项(\(f=0\),无正常通量(\(h=0\)),以及施加在顶边界上的Dirichlet条件(\(g=1)和底部边界(\(g=0))。

对于空间离散化,我们使用二次有限元作为状态变量(也用于伴随变量),使用线性有限元作为参数变量。使用2048个三角形元素的规则网格离散计算域。这导致状态变量和参数变量分别有4225和1089个自由度。在第节所示的可扩展性结果中5.1.6,然后使用多达四个级别的均匀细化对网格进行细化,从而在最精细的级别上分别为状态变量和参数变量提供263169和66049个自由度。

5.1.2先前型号。

如第节所述2,我们选择先验为高斯分布\(\mathcal{N}\!({m_{text{pr}}},{mathcal{C}_{\text{{previous}}}})\(\mathcal{C}_{\text{{previous}}=\mathcal{A}^{-2}\)其中\(\mathcal{A}\)是一个类拉普拉斯算子,如下所示(35)\(\begin{equation}\mathcal){A} 米={\left\lbrace\begin{array}{ll}-\gamma\nabla\cdot(\Theta\nabla-m)+\delta m&\text{in}\Omega,\\Theta\nab la-m\cdot\mathbf{n}+\beta m&text{on}\partial\Omega。\结束{数组}\right。}\结束{方程式}\)

在这里,\(\beta\propto\sqrt{\gamma\delta})是引入的最佳罗宾系数,用于缓解不良边界效应[Daon和Stadler2018]和各向异性张量\(\Theta\)的形式为(36)\(begin{方程式}\Theta=\begin{bmatrix}\Theta_1\sin^2\alpha+\Theta_2\cos^2\ alpha&(\Theta_1-\Theta_2)\sin\alpha\cos\alpha\\(\Theta _1-\ttheta _2)\sin \alpha\ cos\alha&\Theta 1\cos^2 \alpha+Theta_2\sin^2 \阿尔法\end{bmatricx}。\结束{方程式}\)

对于这个示例,我们采用\(伽马=0.1,δ=0.5,β=sqrt{伽马δ}/1.42,θ_1=2.0,θ_2=0.5)和\(α=pi/4)。6显示了之前的平均值\(m{text{pr}})和先验分布中的三个样本。

图6。

图6。泊松问题的先验分布的先验平均值(最左侧)和三个样本场。

5.1.3噪声和可能性观测。

我们在以下位置生成状态观测均匀分布的随机位置\([0.05,0.95]^2)通过求解具有真参数场的最细网格上的正问题\(m_{\scriptscriptstyle\text{true}}\)(这里使用前面的样本),然后向结果状态值添加随机高斯噪声;参见图7.观察次数在本例中设置为300。合成观测值的矢量由下式给出(37)\(\begin{equation}\boldsymbol{d}=\mathcal{B}u+\boldsymbol{eta},\end{equation}\)哪里\(\mathcal{B}\)是一个线性观测算子,将状态解限制为观测点。加性噪声矢量\(\boldsymbol{\eta}\)具有相互独立的分量,这些分量以零平均值和标准偏差正态分布\(σ=0.005)。然后通过以下公式给出似然函数(38)\(\begin{equation}\pi-{\rm like}(\boldsymbol{d}\,|\,m)\proto\exp\left(-\frac{1}{2}\Vert\mathcal{B}\,u(m)-\boldssymbol{d}\Vert^2_{\boldsymbol{\Gamma}_{\!\mbox{noise}}^{-1}}\right),\end{equation}\)哪里\(\boldsymbol{\Gamma}_{\!\mbox{noise}}=\sigma^2\mathbf{I}\)。

图7。

图7。泊松问题的真参数字段(左)和相应的状态字段(右)。观测点的位置在右图中标记为白色方块。

5.1.4拉普拉斯后验近似。

接下来,我们构造后验高斯分布的拉普拉斯近似\({\hat{\mu}_{\mbox{post}}}\sim\mathcal{N}\!({m_{\scriptscriptstyle\text{MAP}}}},{\mathcal{H}(m_{\scriptscriptstyle\text{MAP}})^{-1}})}),平均值等于MAP点,Hessian给出的在MAP点评估的负对数后验的协方差。MAP点是通过最小化负对数后验来获得的,即。,(39)\(\开始{方程式}\下集{m\in\mathcal{m}}{\min}\;\mathcal{J}(m)\;:=\;\裂缝{1}{2}\Vert\mathcal{B}\,u(m)-\boldsymbol{d}\Vert^2_{{bf\Gamma}_{text{nois}}^{-1}}+\frac{1}}{2{\Vertm-m{text{pr}}\Vert_22_{mathcal{C}_{\text{{previous}}}^{-1}}。\结束{方程式}\)

我们使用了实现于hIPPY库解决上述优化问题。我们建议读者参考Villa等人。[2021]对于负对数后验的梯度和Hessian作用的算法和表达式的详细描述\(\mathcal{J}(m)\)。

如第节所述2,显式计算Hessian对于大规模问题是禁止的,因为这需要求解两个与参数数量相同的类前向PDE。为了使使用Hessian的操作在参数维上具有可伸缩性,我们调用Hessian数据不匹配部分的低阶近似,仅保留特征向量是数据中信息最丰富的方向【Villa等人。2021].

8显示了先决条件数据失配Hessian的特征谱。由提供的双通过随机算法\使用过采样因子为20的(\texttt{hIPPYlib})来精确计算主特征对。我们看到,特征值在\(第60位)特征值,表明保持60个特征对足以进行低阶近似。8还显示了四个特征向量,正如预期的那样,这说明对应于较小特征值的特征向量显示出更多的波动。

图8。

图8。与泊松问题的第一、第四、第十六和第64个最大特征值相对应的前-条件数据失配Hessian和特征向量的主导特征值(r=100)的对数图。

在图中9,我们描述了MAP点和从后验拉普拉斯近似中提取的三个样本。

图9。

图9。从泊松问题后验分布的拉普拉斯近似中得出的MAP点(最左侧)和三个样本场。

5.1.5使用MCMC方法探索后验。

在本节中,我们实现了本节中讨论的高级MCMC算法探索后效并比较其表现。

特别是,我们认为pCN、MALA、,\(\infty)-MALA、DR、DILI和他们的Hessian知情同行。对于每种方法,我们模拟了20个独立的MCMC链,每个链有25000个样本(在丢弃2500个样本作为老化样本后),因此从后面总共提取了500000个样本。从拉普拉斯近似值中选取一个样本作为每条链的起点。

为了检验MCMC链的收敛性和统计效率,我们考虑由中广义特征系统的主导特征向量(10),而不是参数向量的所有组件\(\粗体符号{m}\)。具体来说,我们计算MPSRF、自相关时间和关于系数向量的ESS\(\mathbf{c}\in\mathbb{R}^R\)由定义(40)\(\begin{equation}\mathbf{c}=\mathbf{五} _r(r)^粗体符号{\Gamma}_{\!\mbox{previous}}^{-1}\mathbf{m}。\结束{方程式}\)

1显示了MCMC样本的收敛诊断和计算效率。MPSRF和ESS是根据参数样本沿MAP点处的前25个主特征向量的投影来计算的。1报告所有25个预测的最小、最大和平均ESS。

表1。
方法AR(%)MPSRF公司最小ESS(指数)最大ESS(指数)平均ESS核动力源/应急
中国邮政编码(5.0E-3)242.62925 (24)225 (8)845,952
马六甲(6.0E-6)482.64226 (22)874 (5)14810,135
\(infty)-马来半岛(1.0E-5)572.94325 (23)1,102 (5)1609,375
H-pCN(4.0E-1)271.19264 (1)3,598 (15)2,314216
H-MALA(6.0E-2)601.014545 (1)8,868 (19)6,459232
H(H)-\(\infty)-马来半岛(1.0E-1)711.016582(1)8,417 (18)5,905254
DR(H-pCN(1.0E0),H-MALA(6.0E-2))(4, 61)1.013641 (1)12,522 (17)9,222215
DR(H-pCN(1.0E0),高-\(信息)-MALA(2.0E-1))(4, 48)1.011613 (1)12,812 (17)9,141213
DILI-PRIOR(0.8,0.1)(60, 33)1.064314 (1)4,667 (13)3,216548
DILI-LA(0.8,0.1)(83, 36)1.017562 (1)10,882 (17)7,192245
DILI-MAP(0.8,0.1)(77, 22)1.0061,675 (1)10,271 (20)8692个202
  • 验收率(AR)多元潜在标度缩减因子、和有效样本量(ESS)报告以进行收敛诊断。MPSRF和ESS是根据参数样本沿MAP点处的前25个主特征向量的投影来计算的。DR和DILI-MAP中列出了AR的两个值,分别用于第一个和第二个提案步骤。我们还提供每个有效样本的正向和/或伴随PDE解算数(NPS/ES)采样效率。我们使用20个MCMC链,每个链有25000次迭代(总共500000个样本)。每个方法名称中括号中的数字表示使用的参数值(\pCN的(β),\(τ)代表MALA,\(h)用于\(\ infty \)-马来半岛,以及\(β)和\(\tau\)for和DILI)。最小ESS和最大ESS括号内的数字表示相应的特征向量指数。

表1。几种MCMC方法对泊松问题的性能比较:pCN、MALA、(infty)-MALA、DR、DILI及其Hessian通知版本

  • 验收率(AR)多元潜在标度缩减因子、和有效样本量(ESS)报告以进行收敛诊断。MPSRF和ESS是根据参数样本沿MAP点处的前25个主特征向量的投影来计算的。DR和DILI-MAP中列出了AR的两个值,分别用于第一个和第二个提案步骤。我们还提供每个有效样本的正向和/或伴随PDE解算数(NPS/ES)采样效率。我们使用20个MCMC链,每个链有25000次迭代(总共500000个样本)。每个方法名称中括号中的数字表示使用的参数值(\pCN的(β),\(τ)代表MALA,\(h)用于\(\ infty \)-马来半岛,以及\(β)和\(\tau\)for和DILI)。最小ESS和最大ESS括号内的数字表示相应的特征向量指数。

表中的最后一列1表示绘制单个独立样本所需的正向和/或伴随PDE解的数量(使用平均ESS)。该数量可用于测量采样效率,并根据计算效率对方法进行排序。在这个指标下,DILI-MAP是最有效的方法,只需要202个PDE解即可得到有效样本。DR(213 NPS/ES用于H-\(\infty)-MALA和215 NPS/ES用于H-MALA)和H-pCN(216 NPS/ES)接近秒数。

接下来,我们评估感兴趣数量的MCMC样本的收敛性(34)预测的后验分布\(\mathcal{G}(m)\):感兴趣量的自相关函数估计(34)如图所示10(这里,我们使用公式(29)图中描述了三条独立MCMC链的迹线图以及pCN和H-pCN的所有MCMC样品的直方图11.

图10。

图10。自相关函数估计(29)利息数量\(\mathcal{G}\)(34)用于多种MCMC方法。请注意,pCN、MALA和\(\infty\)-MALA的自相关函数图似乎没有变化并重叠。

图11。

图11。左:感兴趣量的跟踪图\(\mathcal{G}\)(34)pCN和H-pCN模拟的三条MCMC链(20条独立链中的);不同的颜色(这里是蓝色、绿色和红色)代表每条链的轨迹。右:感兴趣量的概率密度函数估计(34)根据pCN和H-pCN样本计算;所有50万个样本,20个链,每个链有25000个样本,在直方图中被拉在一起;计数的数量是标准化的,因此图表示概率密度函数;黑色虚线表示从真参数字段(m\text{true})计算出的关注量。

最后,我们比较了不同采样策略对感兴趣量矩的估计。对于每个MCMC链k个第个(\(k=1,2,3\))从参数样本计算的利息量矩\(\mathbf{m} _ i\) (\(i=1,2,\ldots,N;N=25000))计算如下(41)\(开始{方程式}\hat{mathcal{G}}_k=\frac{1}{N}\sum_{i=1}^{N}\mathcal}G}^k(\mathbf{m} _ i). \结束{方程式}\)结果如图所示12如方框和胡须图所示。

图12。

图12。感兴趣量的一阶、二阶和三阶矩估计(k=1,2,3)的Whisker图(41)使用多种MCMC方法进行计算。中心标记为中间值;下四分位数和上四分位数分别代表第25和75个百分位数。晶须延伸到从下四分位或上四分位到四分位间距1.5倍(上四分位数和下四分位数之间的距离)范围内的极端数据点;所有其他数据点都绘制为离群值。每种方法的数据点数量为20,即独立MCMC链的数量。

根据本节中的结果,我们得出以下结论:

  • MAP点的Hessian信息在提高MCMC方法的采样性能方面起着重要作用。事实上,没有Hessian信息的MCMC链并没有在链的整个长度上收敛,而是围绕起点进行定位。只有当MCMC提案利用包含Hessian信息的后验拉普拉斯近似时,才能实现或接近实现收敛。

  • 就每个有效样本的正向和/或伴随PDE解的数量而言,DILI-MAP显示了最佳采样效率。注意MCMC方法中使用的参数值(例如。,\(β)和/或\(τ))不是最优的,不同的参数值可能会得到不同的结果。

我们进一步研究了MCMC方法在不同问题设置下的性能,以便更深入地了解hIPPYlib-MUQ框架的实际使用。

5.1.6Hessian告知pCN的可扩展性。

在这里,我们研究了网格分辨率对采样性能的影响。一种曲率感知的MCMC方法,使用\(β=0.4\)。从粗网格(网格1)到细网格(网格4)的参数和状态变量的维数分别为(10894225)、(422516641)、(1664166049)和(66049263169)。

我们遵循与之前相同的问题设置,并对所有级别使用相同的合成观测值(从最精细的网格生成的真实参数场中获得)。图的左上部分13显示了\(r=100\)先决条件数据的主导特征值与Hessian不匹配。有人观察到,本征谱实际上与网格细化无关。

图13。

图13。左:先-后条件数据失配Hessian(顶部)的主特征值(r=100)的对数图,以及自相关函数估计(29)利息数量\(\mathcal{G}\)(34)(底部)。右:感兴趣量的概率密度函数估计(34); 所有样本,20条链,每条链有25000个样本,总计500000个,在直方图中集中在一起;计数的数量是标准化的,因此图表示概率密度函数。我们考虑四种不同的网格,它们从最粗(网格1)逐渐细化到最细(网格4)。我们使用H-pCN方法(β=0.4)来提取样本。

评估MCMC方法的收敛性,见表2我们报告了后验样本的接受率、MPSRF和ESS。MPSRF和ESS是根据参数样本沿MAP点处的前25个主特征向量的投影来计算的。我们给出了自相关函数估计(29)在图的左下方13,并显示感兴趣数量的直方图\(\mathcal{G}\)(34)在图的右侧13结果表明,样本的收敛性与MPSRF和自相关函数几乎无关。

表2。
尺寸(状态、参数)AR(%)MPSRF公司最小ESS(指数)最大ESS(指数)平均ESS
(4,225, 1,089)271.19264 (1)3,598 (15)2,314
(16,641, 4,225)241.33363 (1)3,221 (18)1830年
(66,049, 16,641)231.075209 (1)3,073 (11)1,940
(263,169, 66,049)221.117102 (2)3,276 (15)1767年
  • 我们使用\对于H-pCN方法(β=0.4\),共抽取500000个样本(20个MCMC链,每个链有25000次迭代)。MPSRF和ESS是根据样本沿MAP点上的前25个主特征向量的投影计算的。最小ESS和最大ESS括号内的数字表示相应的特征向量指数。

表2。不同维度H-pCN方法生成的后验样本的接受率(AR)、多变量潜在标度缩减因子(MPSRF)和有效样本量(ESS)

  • 我们使用\对于H-pCN方法(β=0.4\),共抽取500000个样本(20个MCMC链,每个链有25000次迭代)。MPSRF和ESS是根据样本沿MAP点上的前25个主特征向量的投影计算的。最小ESS和最大ESS括号内的数字表示相应的特征向量指数。

5.2三维Robin边界条件下的系数场反演第页-泊松非线性偏微分方程

到目前为止,我们只关注加性高斯噪声模型。虽然这种加性噪声模型是最常用的模型,但在某些情况下,如在合成孔径雷达图像。在这个例子中,我们考虑了一个不同的噪声模型,其中噪声与观测值成比例。

示例正向模型是一个三维非线性偏微分方程,我们试图在Robin边界条件下推断未知系数场。具体来说,正向控制方程如下所示(42)\(\开始{align}-\nabla\cdot\left(|\nabla u|^{p-2}_\epsilon\nabla u\right)&=f&\quad&\text{in}\Omega,\nonumber\nonumble\\|\nabla u|^{p-2}_\epsilon\nablau\cdot\mathbf{n}+e^mu&=0&&\text{on}\partial\Omega_R,\\|\nablau|^{p-2}_\epsilon\nabla u\cdot\mathbf{n}&=0&&\text{on}\partial\Omega\setminus\partial \Omega _R,\nonnumber\nonnumber \end{align}\)具有\(1)。请注意第页-拉普拉斯学派,\(\nabla\cdot(|\nablau|^{p-2}\nablau)\),当\(p\lt 2)和退化\(p\gt 2)点\(\nabla u=0\)[Lindqvist2017; 棕色2010],所以正则化术语\(\ε\)(这里我们取\(epsilon=1.0\乘以10^{-8}\))被引入到上述等式中,如下所示\(|\nabla u|_\epsilon=\sqrt{|\nablau |^2+\epsilen}\)。这个第页-拉普拉斯算子是拉普拉斯算符的非线性对应物,出现在许多非线性扩散问题(例如,非牛顿流体)中,其中非线性扩散被建模为幂律类型。

我们假设\(p=3\),并考虑一个薄砖域\(\Omega=[0,1]^2\倍[0,0.05]\),带有体积源项(\(f=1))和一个混合边界条件,例如,我们在底部边界表面上施加了一个Robin边界条件,而在剩余的边界表面上没有法向通量。

问题域\(\Omega)使用规则四面体网格离散化,并对所有状态、伴随和参数变量使用线性有限元。离散化后,状态变量和伴随变量的维数为66564,参数变量为16641。

对于前一种情况,我们使用与前一示例中使用的高斯分布相同的高斯分布。我们在以下位置创建合成状态观测\首先用真参数场求解正问题,在顶边界面上均匀分布(l=300)个随机位置\(m{\textrm{true}})从先验样本中获得,然后乘以伽马分布噪声。具体来说,合成观测值的向量的形式如下(43)\(\begin{equation}\boldsymbol{d}=\boldsymbol{eta}\odot\mathcal{B}u,\end{equation}\)哪里\(\odot\)表示组件式乘法\(\boldsymbol{\eta}\)是独立且相同的伽马分布形状\(\kappa\)和刻度\(\ nu\),即。,\(eta_i=\text{Gamma}(\kappa,\nu),i=1,\ldots,q\)。我们接受\(\kappa=1/\nu\),在这种情况下,负对数似然函数的形式为(44)\(\begin{方程式}\Phi(m;\boldsymbol{d})=\kappa\sum_{i=1}^{q}\left(\log(\mathcal{B}u)_i+\frac{d_i}{(\mathcal{B}u).i}\right),\end{方程}\)其中下标方法对应向量的第个分量。在这个例子中,我们设置\(\kappa=10^{4}\)。

14示出了底部边界上的真实参数场、顶部表面上的观测位置以及通过解决最小化负对数后验的优化问题而获得的MAP点。然后,基于MAP点处数据失配Hessian的低秩分解,构造后验函数的Laplace近似。先验-后验数据失配Hessian的谱表明,主特征值(大于1)约为55。

图14。

图14。左:底部曲面上的True参数字段。中间:对应的状态字段和顶面上的\(l=300\)观察点(白色正方形标记)。右:底面上的MAP点。

5.2.1MCMC后部特征分析结果。

这里我们给出了不确定边界系数向量的MCMC模拟结果。在这个例子中,我们考虑H-pCN方法\(β=0.2)并运行40个独立的MCMC链,每个链在丢弃2500个样本作为老化后进行25000次迭代(总共生成1000000个样本)。对于每次MCMC运行,取拉普拉斯近似的后验值作为起点。

列出了MCMC仿真的收敛诊断。将参数样本投影到MAP点处先验-后验数据失配Hessian的前25个特征向量上,并基于该投影对MPSRF和ESS进行评估。

表3。
AR(%)MPSRF公司最小ESS(指数)最大ESS(指数)平均ESS
231.041243 (0)2,891 (22)1,586
  • 我们使用H-pCN方法(\(β=0.2\)),40条链,每个链有25000次迭代(总共1000000个样本)。最小ESS和最大ESS括号内的数字表示相应的特征向量指数。

表3。(p)-Poisson问题的收敛诊断:MAP点上参数样本沿前25个特征向量投影的接受率(AR)、多元势标度缩减因子(MPSRF)和有效样本量(ESS)

  • 我们使用H-pCN方法(\(β=0.2\)),40条链,每个链有25000次迭代(总共1000000个样本)。最小ESS和最大ESS括号内的数字表示相应的特征向量指数。

15显示了后验MCMC样本和拉普拉斯近似的边缘分布。如前所述,这些边缘是根据参数样本到特征向量的投影来计算的。我们观察到,MCMC样本的边缘分布与拉普拉斯近似的边缘分布之间存在明显差异,特别是对于对应于较大特征值的特征向量。

图15。

图15。MCMC样本的边缘分布(红色)和拉普拉斯近似(黑色),均用于将参数样本投影到特征向量上{v} _1个,\粗体符号{v} _2,\粗体符号{v}(v)_{10} ,\text{和},\boldsymbol{v}(v)_{25}\) . 在一维边缘图中,蓝色虚线表示真实参数向量到特征向量的投影。对于二维边缘图,我们只显示了分别代表5%、50%和95%分布的三条等高线。我们注意到比例上的差异,因此在第一列中,估计值和拉普拉斯近似值之间出现了更大的差异。

已经对勘探与发现多环境计算机(MERCED)UC Merced的集群。虽然根据作业分配给的集群CPU节点,性能会有一些变化,但第一个示例大约需要0.07秒的单个PDE解算,第二个三维示例大约需要4秒(两者都是经过的时间)。

跳过6结论部分

6结论

我们提出了一个健壮且可扩展的软件框架,用于解决由PDE控制的大规模贝叶斯反问题。该软件集成了两个互补的开放源码软件库hIPPYlib和MUQ,形成了一个独特的软件框架,解决了PDE控制的反问题的贝叶斯解的禁止性。拟议软件框架的主要目标是

(1)

为领域科学家提供一套复杂且计算效率高的MCMC方法,利用贝叶斯反问题结构;

(2)

允许研究人员轻松实施新方法,并与现有技术进行比较。

这两个库的集成使得先进的MCMC方法能够利用参数空间的几何形状和内在低维性,从而高效且可扩展地探索后验分布。特别是,后验的拉普拉斯近似被用来生成高质量的MCMC提案。这种近似基于对数后验函数的Hessian逆,通过对数似然函数的Hesian低阶近似变得容易处理。基于线性和非线性PDE的贝叶斯反问题的数值实验表明,基于Laplace的方案能够通过以下因素加速MCMC采样\(\sim 50\次)。

尽管这些先进的结构扩展MCMC方法具有快速且依赖于维数的收敛性,但许多由代价高昂的PDE控制的贝叶斯反问题仍然遥不可及。例如,第节的结果5.1.5对于泊松系数反问题,表明\即使使用最有效的MCMC方法,也可能仍然需要(O(10^6)\)PDE解。在这种情况下,hIPPYlib-MUQ可以用作原型环境,以研究进一步利用问题结构的新方法,例如通过使用各种简化模型(例如,Cui等人[2016年b])或者通过超越低阶Alger等人的高级Hessian近似。[2019]; Ambartsumyan等人。[2020].

我们还谈到了我们的软件框架的局限性。由于我们依赖FEniCS进行PDE的有限元近似,hIPPYlib-MUQ继承了该软件的所有局限性和挑战。然而,FEniCS、hIPPYlib和MUQ开发人员和用户社区提供了丰富的支持基础,hIPPYlib-MUQ的用户可以在此基础上进行构建。对替代有限元实现的支持是未来工作的主题。hIPPY lib-MUQ当前的另一个限制是缺少MCMC方法的并行实现。因此,hIPPYlib-MUQ未来版本的目标是合并多级并行。这将包括现在可用的并行PDE解算器和额外的并行MCMC链,并将允许使用高维参数空间解决更复杂的基于PDE的贝叶斯反问题。

软件可用性.hIPPYlib-MUQ在GNU通用公共许可证版本3(GPL3)下发布。hIPPYlib-MUQ项目托管在GitHub上(https://github.com/hippylib/hippylib2muq)并使用Travis CI进行持续集成。hIPPYlib-MUQ使用语义版本控制。本工作中的结果是通过hIPPYlib-MUQ 0.3.0版、hIPPYlib 3.1.0版和MUQ 0.3.5版获得的。Docker形象[默克尔2014]包含预安装的软件和示例,请访问https://hub.docker.com/r/ktkimyu/hippylib2muq.hIPPYlib-MUQ文档托管在ReadTheDocs上(https://hippylib2muq.readthedocs.io). 鼓励用户加入Slack上的hIPPYlib和MUQ工作区,与其他用户联系,获取帮助,并讨论新功能;看见https://hippylib.github.io/松弛-通道https://mituq.bitbucket.io有关如何加入的更多信息。

跳过ACKNOWLEDGMENTS部分

致谢

我们感谢匿名审稿人对我们手稿的仔细阅读,以及他们的宝贵意见和建议,这有助于我们提高手稿的质量。

脚注

  1. 1 在反演参数是空间变化场的无限维反问题中,使用高斯测度(通常使用Matérn协方差核)是一种常见的选择。在hIPPYlib-MUQ中,非高斯先验概率分布可以通过引入非线性可逆的改变变量,将期望的先验分布映射为高斯分布。

    脚注

参考文献

  1. 阿克利克·沃尔坎比罗斯·乔治加塔斯·奥马尔希尔·朱迪思凯斯·戴维、和Waanders Bart van Bloeman公司.2006.并行PDE约束优化.英寸科学计算的并行处理希罗克斯·M·。拉加文·P·。、和西蒙·H。(编辑)。暹罗.谷歌学者谷歌学者
  2. 阿尔杰·N。拉奥五世。迈耶斯A。Bui-Thanh T公司。、和加塔斯·O。.2019.局部平移不变算子的无矩阵自适应乘积卷积逼近.SIAM科学计算杂志 414(2019),A2296–A2328.https://arxiv.org/abs/1805.06.2018.谷歌学者谷歌学者数字图书馆数字图书馆
  3. Ambartsumyan Ilona公司布卡拉姆·瓦吉布汉坦加塔斯·奥马尔凯斯·戴维斯塔德勒·乔治土基亚·乔治、和扎皮尼·斯特凡诺.2020.偏微分方程反问题中Hessian的层次矩阵逼近.SIAM科学计算杂志 425(2020),A3397–A3426.谷歌学者谷歌学者数字图书馆数字图书馆
  4. AtchadéYves F。.2006.具有截断漂移的Metropolis调整Langevin算法的自适应版本.应用概率论中的方法论与计算 8(2006),235254.谷歌学者谷歌学者数字图书馆数字图书馆
  5. 巴莱·萨蒂什阿比扬卡尔·施里朗亚当斯·马克·F。布朗·杰德布鲁恩·彼得布歇尔曼·克里斯达尔辛·利桑德罗德纳阿尔卑斯山埃伊霍特·维克多格罗普·威廉·D·。考希克·迪内什内普利·马修·G。5月Dave A。麦克因斯·路易斯·柯夫曼Mills Richard Tran公司蒙森·托德鲁普·卡尔萨南·帕特里克史密斯·巴里·F·。扎皮尼·斯特凡诺、和张红.2018.PETSc网页.http://www.mcs.anl.gov/petsc.http://www.mcs.anl.gov/petsc.谷歌学者谷歌学者
  6. 巴莱·萨蒂什阿比扬卡尔·施里朗亚当斯·马克·F。布朗·杰德布鲁恩·彼得布歇尔曼·克里斯埃伊霍特·维克多格罗普·威廉·D·。考希克·迪内什内普利·马修·G。麦克因斯·路易斯·柯夫曼鲁普·卡尔史密斯·巴里·F·。、和张红.2014.PETSc网页.http://www.mcs.anl.gov/petsc.http://www.mcs.anl.gov/petsc.谷歌学者谷歌学者
  7. 巴兹利·乔纳坦(Bardsley Johnathan M.)。崔天刚马尔祖克·优素福M。、和王正.2020.函数空间上基于可扩展优化的采样.SIAM科学计算杂志 422(2020),A1317–A1347.谷歌学者谷歌学者数字图书馆数字图书馆
  8. 贝克尔E.B。凯里·G.F。、和奥登·J·T。.1981.有限元:导论,第一卷.普伦蒂斯·霍尔新泽西州恩格尔伍德悬崖.谷歌学者谷歌学者
  9. 贝斯科斯·亚历山德罗斯吉罗米·马克蓝世伟法雷尔·帕特里克·E。、和斯图亚特·安德鲁·M。.2017.无穷维反问题的几何MCMC.J.计算。物理学。 335(2017),327351.谷歌学者谷歌学者交叉引用交叉引用
  10. 波兹·阿尔菲奥舒尔兹·沃尔克.2012.偏微分方程控制系统的计算优化.暹罗.谷歌学者谷歌学者数字图书馆数字图书馆
  11. 布鲁克斯·斯蒂芬·P。盖尔曼·安德鲁.1998.监测迭代模拟收敛性的一般方法.计算与图形统计杂志 74(12月。1998),434455.内政部:谷歌学者谷歌学者交叉引用交叉引用
  12. 布朗·杰德.2010.三维节点高阶有限元的高效非线性求解器.科学计算杂志 451(2010),4863.内政部:谷歌学者谷歌学者数字图书馆数字图书馆
  13. 布汉坦卡斯滕火山爆发加塔斯·奥马尔马丁·詹姆斯斯塔德勒·乔治、和威尔科克斯·卢卡斯(Wilcox Lucas C.)。.2012.偏微分方程控制的贝叶斯反问题的极值尺度UQ.英寸SC12:高性能计算、网络、存储和分析国际会议记录.戈登·贝尔奖决赛选手。谷歌学者谷歌学者数字图书馆数字图书馆
  14. Bui-Thanh T公司。加塔斯·O。马丁·J·。、和斯塔德勒·G。.2013.无限维贝叶斯反问题的计算框架第一部分:线性化情况及其在全球地震反演中的应用.SIAM科学计算杂志 356(2013),A2494–A2523.谷歌学者谷歌学者数字图书馆数字图书馆
  15. 卡尔德黑德·本.2014.Metropolis-Hastings算法并行化的一般构造.美国国家科学院院刊 11149(2014),1740817413.谷歌学者谷歌学者交叉引用交叉引用
  16. 卡塞拉·乔治乔治·爱德华一世。.1992.解释吉布斯采样器.美国统计学家 46(1992),167174.谷歌学者谷歌学者交叉引用交叉引用
  17. 康拉德·帕特里克·R。马尔祖克·优素福M。.2013.自适应Smolyak伪谱近似.SIAM科学计算杂志 356(2013),A2643–A2670.内政部:谷歌学者谷歌学者数字图书馆数字图书馆
  18. 科特S.L。罗伯茨G.O。斯图亚特·A.M。、和白色D。.2012.函数的MCMC方法:修改旧算法使其更快. (2012).提交。谷歌学者谷歌学者
  19. 崔T。法律K.J.H。、和马尔佐克Y.M。.2016年a.独立于维度的类似知情MCMC.J.计算。物理学。 304(2016),109137.谷歌学者谷歌学者数字图书馆数字图书馆
  20. 崔天刚马尔佐克·优素福、和威尔科克斯·卡伦.2016年b.基于似然参数和状态约简的大规模贝叶斯反问题的可缩放后验近似.J.计算。物理学。 315(2016),363387.谷歌学者谷歌学者数字图书馆数字图书馆
  21. 崔天刚扎姆·奥利维尔.2021.贝叶斯反问题的无数据似然信息降维.反问题 374(2021),045009.谷歌学者谷歌学者交叉引用交叉引用
  22. 道尼尔斯塔德勒·乔治.2018.减轻边界条件对椭圆偏微分方程协方差算子的影响.反问题和成像 125(2018),10831102.arXiv:1610.05280谷歌学者谷歌学者交叉引用交叉引用
  23. 多德维尔·蒂姆·J。凯特尔森克里斯蒂安谢利克·罗伯特、和Teckentrup Aretha L.公司。.2019.多级马尔可夫链蒙特卡罗.SIAM版本。 61(2019),509545.谷歌学者谷歌学者交叉引用交叉引用
  24. 埃文斯·M。斯瓦茨T。.2000.蒙特卡罗和确定性方法逼近积分第卷。20.牛津大学出版社.谷歌学者谷歌学者
  25. Flegal James M。琼斯·加林(Jones Galin L.)。.2010.马尔可夫链蒙特卡罗中的批均值和谱方差估计.统计年鉴 382(2010),10341070.谷歌学者谷歌学者交叉引用交叉引用
  26. 盖尔曼·安德鲁卡林·约翰·B。斯特恩·哈尔·S。、和鲁宾·唐纳德B。.2004.贝叶斯数据分析.谷歌学者谷歌学者
  27. 盖尔曼·安德鲁鲁宾·唐纳德B。.1992.使用多序列的迭代模拟推断.统计科学(1992),457472.谷歌学者谷歌学者交叉引用交叉引用
  28. 加塔斯·O。威尔科克斯K。.2021.从数据中学习基于物理的模型:从反问题和模型简化的角度.数字学报 30(2021),445554.内政部:谷歌学者谷歌学者交叉引用交叉引用
  29. 吉罗米·马克卡尔德黑德·本.2011.黎曼流形Langevin和Hamilton蒙特卡罗方法.英国皇家统计学会杂志:B辑(统计方法) 732(2011),123214.谷歌学者谷歌学者交叉引用交叉引用
  30. Golub基因H。贷款Charles F.Van.1996.矩阵计算(第三版).约翰霍普金斯大学出版社马里兰州巴尔的摩.谷歌学者谷歌学者数字图书馆数字图书馆
  31. 哈里奥·海基萨克斯曼·埃罗、和塔米尼·约翰纳.2001.一种自适应metropolis算法.伯努利 72(九月。2001),223242.内政部:谷歌学者谷歌学者交叉引用交叉引用
  32. 海尔·马丁斯图亚特·安德鲁·M。、和沃尔默·塞巴斯蒂安J。.2014.无限维Metropolis–Hastings算法的光谱间隙.应用概率年鉴 246(2014),24552490.谷歌学者谷歌学者交叉引用交叉引用
  33. 哈尔科·内森Martinsson Per Gunnar公司、和特罗普·乔尔A。.2011.寻找随机结构:构造近似矩阵分解的概率算法.SIAM版本。 532(2011),217288.谷歌学者谷歌学者数字图书馆数字图书馆
  34. Hartikainen Jouni公司萨科卡·西莫.2010.时间高斯过程回归模型的卡尔曼滤波和平滑解.英寸2010年IEEE信号处理机器学习国际研讨会IEEE标准,379384.内政部:谷歌学者谷歌学者交叉引用交叉引用
  35. 黑斯廷斯·W·基思.1970.马尔可夫链蒙特卡罗抽样方法及其应用.生物特征 571(1970),97109.谷歌学者谷歌学者交叉引用交叉引用
  36. 霍夫曼·马修·D·。盖尔曼·安德鲁.2014.无转取样器:在哈密顿蒙特卡罗中自适应设置路径长度。 机器学习研究杂志 151(2014),15931623.谷歌学者谷歌学者
  37. 艾萨克·托宾佩特拉·诺埃米斯塔德勒·乔治、和加塔斯·奥马尔.2015.通过推理将数据中的不确定性传播到大规模问题预测的可扩展高效算法,并应用于南极冰盖的流动.J.计算。物理学。 296(九月2015),348368.内政部:谷歌学者谷歌学者数字图书馆数字图书馆
  38. 凯皮奥·贾里萨默萨洛·埃尔基.2005.统计和计算反问题.应用数学科学,卷。160.纽约斯普林格-弗拉格.内政部:谷歌学者谷歌学者交叉引用交叉引用
  39. 林格伦·芬恩赫瓦尔德街、和林德斯特伦·约翰.2011.高斯场和高斯-马尔可夫随机场之间的显式联系:随机偏微分方程方法.英国皇家统计学会杂志:B辑(统计方法) 734(2011),423498.内政部:谷歌学者谷歌学者交叉引用交叉引用
  40. 林克维斯特·彼得.2017.关于p-Laplace方程的注记161号。Jyväskylä大学.谷歌学者谷歌学者
  41. 洛格·安德斯马尔代尔·肯特·安德雷、和Wells Garth N.公司。(编辑)。2012.微分方程的有限元自动求解.计算科学与工程课堂讲稿,卷。84.施普林格.内政部:谷歌学者谷歌学者交叉引用交叉引用
  42. 马歇尔·特里斯坦罗伯特·加雷斯.2012.Langevin MCMC的自适应方法.统计与计算 225(9月。2012),10411057.内政部:谷歌学者谷歌学者数字图书馆数字图书馆
  43. 马丁·詹姆斯威尔科克斯·卢卡斯(Wilcox Lucas C.)。卡斯滕火山爆发、和加塔斯·奥马尔.2012.大规模统计反演问题的随机Newton-MCMC方法及其在地震反演中的应用.SIAM科学计算杂志 34(2012),A1460–A1487.谷歌学者谷歌学者数字图书馆数字图书馆
  44. 马尔佐克·优素福莫塞利·塔雷克帕诺·马修、和斯潘蒂尼·阿莱西奥.2016.通过测量传输采样:简介.施普林格国际出版公司141.内政部:谷歌学者谷歌学者交叉引用交叉引用
  45. 默克尔·德克.2014.Docker:用于一致开发和部署的轻量级Linux容器.Linux J。 2014239,文章2(2014).http://dl.acm.org/citation.cfm?id=2600239.2600241.谷歌学者谷歌学者数字图书馆数字图书馆
  46. 尼古拉斯大都会罗森布鲁斯·阿里安娜·W。罗森布鲁斯·马歇尔N。特勒·奥古斯塔H。、和出纳员爱德华.1953.快速计算机器的状态方程计算.化学物理杂志 216(1953),10871092.内政部:谷歌学者谷歌学者交叉引用交叉引用
  47. 米拉·安东尼塔.2001.具有延迟拒绝的Metropolis-Hastings算法.Metron公司 593-4(2001),231241.谷歌学者谷歌学者
  48. 尼尔·R·M。.2010.马尔可夫链蒙特卡罗手册.查普曼和霍尔/CRC出版社,第MCMC章使用哈密顿动力学。谷歌学者谷歌学者
  49. 欧文艺术B。.2013.蒙特卡罗理论、方法和示例. (2013).谷歌学者谷歌学者
  50. 帕诺·马修戴维斯·安德鲁康拉德·帕特里克、和马尔佐克Y.M。.2014.麻省理工学院不确定性量化(MUQ)库.https://muq.mit.edu.谷歌学者谷歌学者
  51. 帕诺·马修·D·。马尔祖克·优素福M。.2018.运输图加速马尔可夫链蒙特卡罗.SIAM/ASA不确定性量化杂志 62(2018),645682.内政部:谷歌学者谷歌学者交叉引用交叉引用
  52. 佩特拉·诺埃米马丁·詹姆斯斯塔德勒·乔治、和加塔斯·奥马尔.2014.无限维贝叶斯反问题的计算框架:第二部分。随机牛顿MCMC及其在冰盖反问题中的应用.SIAM科学计算杂志 364(2014),A1525–A1555.谷歌学者谷歌学者数字图书馆数字图书馆
  53. 平斯基·弗兰克·J·。辛普森-吉迪恩斯图亚特·安德鲁·M。、和韦伯·亨德里克.2015.无限维概率测度的Kullback–Leibler近似算法.SIAM科学计算杂志 376(2015),A2733–A2757.谷歌学者谷歌学者数字图书馆数字图书馆
  54. 按S.J。.2003.主客观贝叶斯统计:原理、方法和应用.纽约威利.谷歌学者谷歌学者
  55. 拉斯穆森·卡尔·爱德华威廉斯·克里斯托弗·K·I。.2005.机器学习的高斯过程.麻省理工学院出版社.谷歌学者谷歌学者交叉引用交叉引用
  56. 罗伯特·克里斯蒂安·P。卡塞拉·乔治.2005.蒙特卡洛统计方法(统计学中的斯普林格文本).Springer-Verlag纽约公司。美国新泽西州塞考克斯.谷歌学者谷歌学者数字图书馆数字图书馆
  57. 罗伯茨·加雷思·O。罗森塔尔·杰弗里。,等人.2004.一般状态空间马尔可夫链和MCMC算法.概率调查 1(2004),2071.谷歌学者谷歌学者交叉引用交叉引用
  58. 罗伯茨·加雷思·O。Stramer Osnat公司.2003.Langevin扩散和Metropolis-Hastings算法.应用概率论中的方法论与计算 4(2003),337357.谷歌学者谷歌学者数字图书馆数字图书馆
  59. 鲁道夫·D·。斯普伦克B。.2018.预处理Crank-Nicolson Metropolis算法的推广.计算数学基础 18(2018),309343问题2。谷歌学者谷歌学者数字图书馆数字图书馆
  60. 斯蒂格勒。S.M.公司。1986.拉普拉斯1774年关于逆概率的回忆录.统计师。科学。 1(081986),359363.内政部:谷歌学者谷歌学者交叉引用交叉引用
  61. 绞合线G。固定G.J。.1988.有限元法分析.韦尔斯利-剑桥出版社马萨诸塞州韦尔斯利.谷歌学者谷歌学者
  62. 斯图亚特·安德鲁·M。.2010.反问题:贝叶斯视角.数字学报 19(2010),451559.内政部:谷歌学者谷歌学者交叉引用交叉引用
  63. 塔兰托拉·阿尔伯特.2005.模型参数估计的反问题理论与方法.暹罗宾夕法尼亚州费城.xii+342页。谷歌学者谷歌学者交叉引用交叉引用
  64. 蒂尔尼L。卡丹·J·B。.1986.后验矩和边缘密度的精确近似.J.Amer。统计师。协会。 81393(1986),8286.内政部:谷歌学者谷歌学者交叉引用交叉引用
  65. Trilinos项目团队.2020(2020年5月22日访问)。Trilinos项目网站.https://trilinos.github.io.谷歌学者谷歌学者
  66. 弗雷迪(Tröltzsch Fredi).2010.偏微分方程的最优控制:理论、方法和应用.数学研究生课程,卷。112.美国数学学会.谷歌学者谷歌学者交叉引用交叉引用
  67. Vats斗提卡Flegal James M。、和琼斯·加林(Jones Galin L.)。.2019.马尔可夫链蒙特卡罗的多元输出分析.生物特征 1062(2019),321337.谷歌学者谷歌学者交叉引用交叉引用
  68. 维塔里·阿基盖尔曼·安德鲁辛普森-丹尼尔木匠鲍勃、和布尔克纳·保罗-克里斯蒂安.2020.秩规范化、折叠和局部化:一种改进的\(\widehat{R})用于评估MCMC的收敛性(与讨论).贝叶斯分析 162(2020),126.内政部:arxiv:arxiv:1903.08008v5谷歌学者谷歌学者交叉引用交叉引用
  69. 翁贝托别墅佩特拉·诺埃米、和加塔斯·奥马尔.2020.hIPPYlib用户手册:版本2. (62020).内政部:谷歌学者谷歌学者交叉引用交叉引用
  70. 翁贝托别墅佩特拉·诺埃米、和加塔斯·奥马尔.2021.hIPPYlib:由PDE管理的大规模反问题的可扩展软件框架:第一部分:确定性反演和线性贝叶斯推理.ACM事务处理。数学。柔和。 472,文章16(四月2021),共34页。内政部:谷歌学者谷歌学者数字图书馆数字图书馆
  71. 威廉姆斯·大卫.1991.鞅概率.剑桥大学出版社.谷歌学者谷歌学者交叉引用交叉引用
  72. Wolff Ulli公司协作Alpha,等人.2004.蒙特卡罗误差较小.计算机物理通信 1562(2004),143153.谷歌学者谷歌学者交叉引用交叉引用
  73. Wong R。.2001.积分的渐近逼近.工业和应用数学学会.内政部:arXiv:https://epubs.siam.org/doi/pdf/101137/1.9780898719260谷歌学者谷歌学者交叉引用交叉引用
  74. 扎姆·奥利维尔崔天刚法律科迪斯潘蒂尼·阿莱西奥、和马尔佐克·优素福.2022.非线性贝叶斯反问题的证明降维.数学。公司。 91336(2022),17891835.谷歌学者谷歌学者交叉引用交叉引用

索引术语

  1. hIPPYlib-MUQ:一种贝叶斯推理软件框架,用于不确定性下数据与复杂预测模型的集成

                  建议

                  评论

                  登录选项

                  检查您是否可以通过登录凭据或您的机构访问本文。

                  登录

                  完全访问权限

                  • 发布于

                    封面图片ACM数学软件汇刊
                    ACM数学软件汇刊 第49卷第2期
                    2023年6月
                    275页
                    国际标准编号:2009年8月35日
                    EISSN公司:1557-7295
                    内政部:10.1145/3604595
                    期刊目录

                    如果复制品不是为了盈利或商业利益而制作或分发的,并且复制品的第一页载有本通知和完整引文,则允许免费制作本作品的全部或部分数字或硬拷贝以供个人或课堂使用。必须尊重作者以外的其他人对本作品组成部分的版权。允许用信用证进行摘要。要以其他方式复制或重新发布,在服务器上发布或重新发布到列表,需要事先获得特定许可和/或付费。从请求权限[电子邮件保护].

                    出版商

                    计算机协会

                    美国纽约州纽约市

                    出版历史

                    • 出版:2023年6月15日
                    • 在线AM:2023年2月9日
                    • 认可的:2022年12月15日
                    • 修订日期:2022年10月30日
                    • 收到:2021年11月29日
                    发布于汤姆斯第49卷第2期

                    权限

                    请求有关此文章的权限。

                    请求权限

                    检查更新

                    限定符

                    • 研究论文

                  PDF格式

                  以PDF文件查看或下载。

                  PDF格式

                  电子阅读器

                  使用eReader联机查看。

                  电子阅读器