高效计算的可能性推理模型近似方法

Computationally efficient variational-like approximations ofpossibilistic inferential models

https://www.sciencedirect.com/science/article/pii/S0888613X25001471?via%3Dihub

打开网易新闻 查看精彩图片
打开网易新闻 查看精彩图片

摘 要

推断模型(Inferential Models, IMs)提供了一种可证明可靠的、数据驱动的、可能性(possibilistic)统计推断方法。然而,尽管IM框架在理论和基础层面具有优势,其高效计算仍是一个挑战。本文提出了一种简单而强大的数值策略,用于近似IM的可能性轮廓(possibility contour),或至少近似其在给定显著性水平 α ∈ ( 0 , 1 )下的 α α-截集( α α-cut)。我们的方法首先指定一个参数化族,该族在某种意义上近似覆盖了与IM可能性测度相关联的信度集(credal set)。类似于变分推断,我们随后调整该参数化族的参数,使其 100 ( 1 − α ) % 可信集(credible set)大致匹配IM轮廓的 α α-截集。这种参数化的 α α-截集匹配策略,能够以远低于以往方法的计算成本,实现对IM可能性轮廓的完整近似。

关键词:贝叶斯、置信区域、信度集、Fiducial(信念推断)、蒙特卡洛、随机近似

  1. 引言

长期以来,尽管贝叶斯方法在理论基础上具有优势,但实际使用贝叶斯方法的统计学家却寥寥无几——其计算负担实在过于沉重。这一局面在蒙特卡洛方法使贝叶斯解变得触手可及时发生了显著改变。而近年来,随着各类近似贝叶斯计算方法(尤其是文献[8]、[73]及其所引文献中的变分近似方法)的发展,情况再次发生转变。曾经清晰划分贝叶斯与非贝叶斯方法在计算可行性上的界限如今已变得模糊,这重新激发了贝叶斯方法在现代应用中的活力。丹尼斯·林德利(Dennis Lindley)曾预言:“到2020年,(统计学家)都将变成贝叶斯主义者”[65]——他的预言并未实现,但可以说,贝叶斯学派如今比以往任何时候都更强大。

目前,贝叶斯与频率学派是统计推断中的两大主流思想流派,但它们并非仅有的视角。例如,Dempster–Shafer 理论最初便是作为对贝叶斯推断和费希尔(Fisher)信念推断(fiducial argument)的改进与推广而提出的。我们在此特别关注的是推断模型(Inferential Models, IMs;参见[44, 47, 52, 54])方面的最新进展。IM 框架提供了一种类似贝叶斯的、数据依赖的、基于可能性(possibilistic)的未知量不确定性量化方式,同时内嵌了类似频率学派的可靠性保证。IM 及其他新兴/非传统框架目前正面临当年贝叶斯学派曾遭遇的相同计算挑战:我们知道想要计算什么以及为何要计算它,但目前缺乏高效完成该任务的工具。虽然传统的蒙特卡洛方法仍然有用(见第2节),但 IM 可靠性保证所依赖的“不精确性”(imprecision)意味着仅靠蒙特卡洛方法是不够的。要实现埃夫龙(Efron)关于信念推断类方法的预言——“也许费希尔最大的失误将在21世纪大获成功!”[24]——就必须在蒙特卡洛计算中取得能够容纳不精确性的新进展。本文的贡献正是朝向这一总体方向,提出了一种新的 IM 近似方法及一种高效计算该近似的算法。

我们的出发点是一个相对简单的想法,却能导向一个通用的工具,用于实现计算高效且统计可靠的可能性推断(possibilistic inference)。如第2节所述,我们聚焦于基于可能性理论的 IM(possibility-theoretic IMs),这类 IM 完全由其对应的轮廓函数(contour function)决定,或者等价地,由该轮廓函数的所谓 α-截集(α-cuts)决定。我们利用一个广为人知的刻画:一个可能性测度的信度集(credal set)可被表示为所有满足“对上述 α-截集赋予至少 1−α 概率”的概率测度的集合。在 IM 的语境中,这些 α-截集即为 100(1−α)% 的置信区域,因此 IM 信度集中的元素可合理地解释为“置信分布”(confidence distributions;例如[62, 77])。具体而言,信度集中“最弥散”的那个元素(即“内概率近似”,inner probabilistic approximation)会尽可能接近 1−α 地为每个 α-截集分配概率。如果我们能通过蒙特卡洛或其他方式近似这一特殊的信度集成员,那么我们就离完成大部分(甚至全部)相关的 IM 计算不远了。挑战在于,除了极少数特殊问题类别[48]外,这种“内概率近似”通常相当复杂。然而,如果我们仅需对某个单一的 α-截集(例如 α = 0.1)进行精确近似,那么我们可以获得一个相对简单的近似方案。

为此,我们提议在参数空间上引入一个简单的参数化概率分布族(例如高斯分布),其部分参数依赖于数据,然后调整该分布中未指定的参数,使得该分布(近似地)对指定的 α-截集赋予概率 1−α。这类似于变分贝叶斯方法:我们试图用一个适当选择的、相对简单的概率分布族中的成员,去近似一个复杂的概率分布——在我们的情形中,这个复杂分布是 IM 可能性测度的“内概率近似”,而非贝叶斯后验分布。本文所提方法的具体技术细节受到近期文献[35]以及看似无关的文献[69]中发展的启发。

本文其余部分安排如下:第2节简要回顾可能性 IM 及其性质;第3节提出我们基本但极具普适性的类变分 IM 近似方法,该方法结合蒙特卡洛与随机近似来调整所设参数化近似族的索引参数,特别适用于低维参数的统计推断问题,并给出若干示例。第4节则提出该类变分近似方法的一个更精细版本,更适合高维问题,但主要适用于高斯近似族——这在实践中并无限制。借助更具结构化的近似形式,我们可通过减少 IM 轮廓函数的蒙特卡洛评估次数来降低计算成本。我们在多个例子中展示了这一点,包括一个带 Lasso 惩罚的相对高维问题,以及分别涉及参数、非参数和半参数模型中冗余参数(nuisance parameters)的问题。需要强调的是,本文并非提出一种新的 IM 构造方法,而是提出一种对已在文献中发展和研究过的 IM 解的新颖且计算高效的近似方法。因此,本文的示例并不将我们的 IM 近似与其他方法(如贝叶斯方法)进行比较,而是聚焦于新近似方法的质量,展示其与通常计算成本高得多的精确 IM 解高度吻合。如果所提近似足够准确,那么无论使用旧的昂贵计算策略还是本文提出的新高效策略,IM 解与其他方法之间的比较结果(参见例如[10–14, 41, 47, 49])将保持一致。不过,我们仍包含了一个 IM 与贝叶斯方法的简要对比(见例5)。第6节对全文进行简明总结,并讨论若干具有实际意义的扩展方向。

  1. 可能性推断模型(IMs)的背景

IM 框架最初的表述(例如,[52, 54])严重依赖于(嵌套的)随机集及其相应的信念函数。最近在 [47] 中提出的 IM 公式,建立在 [41, 42] 的发展基础之上,通过将概率到可能性的转换应用于相对似然函数来定义 IM 的可能性轮廓。这一看似微小但重要的转变具有理论动机,但我们在此仅作简要提及。本文综述的重点在于可能性 IM 的公式、其关键性质以及现有的计算策略。

考虑一个由参数空间 T ⊆ ℝᵈ 索引的参数统计模型 {Pθ : θ ∈ T}。例子包括 Pθ = Ber(θ)、Pθ = N(θ, 1)、Pθ = Gamma(α, β)(其中 θ = (α, β)),以及其他许多模型;参见第 3.3 节。非参数问题——参见 [12]、[13] 和 [49, 第 5 节]——同样可以处理,但我们将其讨论推迟至第 5 节。假设可观测数据 Xⁿ = (X₁, ..., Xₙ) 是来自分布 Pθ 的独立同分布样本,其中 Θ ∈ T 是未知/不确定的“真实值”。该模型与观测数据 Xⁿ = xⁿ 共同确定了一个似然函数 θ ↦ Lₓⁿ(θ) 以及相应的相对似然函数

打开网易新闻 查看精彩图片
打开网易新闻 查看精彩图片

其中 H 表示关于 Θ 的一个假设。这种纯粹由似然驱动的可能性测度具有若干理想的性质:例如,基于对可能性得分进行阈值化(使用与模型无关的阈值)的推断满足似然原理(例如,[4,6]),并且在标准正则条件下,它在 PΘ-概率意义下渐近一致(即当 n → ∞ 时,收敛于集中在 Θ 上的可能性测度)。然而,这种纯粹基于似然的可能性测度所缺乏的是一个校准性质(相对于所设定的模型),该性质赋予其对未知 Θ 的假设所分配的“可能性”以意义或信念形成方面的推断权重。更具体地说,如果我们以可能性测度作为对统计不确定性的量化,则其对应的信度集应包含在统计意义上具有意义的概率分布。我们假设先验信息为空白,因此不存在有意义或特殊的贝叶斯后验分布。对信度集元素唯一其他自然的约束是它们应为“置信分布”。但是,正如在 (8) 中所述,这种解释要求相对似然的 α-截集(即 {θ : R(xⁿ, θ) > α})是 Θ 的 100(1−α)% 置信集,而通常情况并非如此。因此,仅凭相对似然本身是不够的。

幸运的是,至少从概念上讲,通过应用 Martin [46] 所称的“有效化”(validification)——一种概率到可能性的转换版本(例如,[22,31])——可以实现这种校准。具体而言,对于观测数据 Xⁿ = xⁿ,可能性 IM 的轮廓定义为

打开网易新闻 查看精彩图片

换言之,如果 IM 以不超过 α 的可能性赋予真实假设(作为数据 Xⁿ 的函数)的比例不超过 α,则该 IM 是有效的(或已校准的)。这正是 IM 上述“推断权重”的来源:公式 (6) 意味着,当 H 为真时,我们不期望 Π̄ₓₙ(H) 很小,因此我们倾向于怀疑那些 Π̄ₓₙ(H) 较小的假设 H 的真实性。此外,上述性质确保了可能性 IM 不会陷入虚假置信(false confidence)[2,43,50],而所有默认先验的贝叶斯方法和信念推断(fiducial)解则无法避免这一点。甚至更强的一个版本——在所有假设上一致成立的 (6) 版本——也成立,如 [15] 中所示/讨论:

打开网易新闻 查看精彩图片

“对于某些包含 Θ 的假设 H”这一事件,可被视为所有包含 Θ 的假设 H 的并集,这使得它比公式 (6) 中与任何单个固定假设 H 相关联的事件要宽泛得多。因此,无论评估的假设数量多少或它们被选择的方式如何——即使它们依赖于数据——从 IM 得出的任何一个建议具有误导性的概率仍被控制在指定水平内。关于可能性 IM 性质及其与贝叶斯/信念推断联系等的更多细节和讨论,请参见 [47–49]。

在贝叶斯分析中,推断基于数据依赖后验分布的摘要,例如科学相关假设的后验概率、损失/效用函数的期望等。而所有这些摘要最终都归结为涉及决定后验的概率密度函数的积分。对于可能性 IM 而言,几乎完全相同的情况也成立:科学相关假设的上下概率对、损失/效用函数的上下期望等,最终都归结为涉及可能性轮廓 πₓₙ 的优化问题。例如,若关注点是 Θ 的某个特征 Φ = g(Θ),则贝叶斯方法可通过积分 Θ 的后验密度来获得 Φ 的边缘后验分布。类似地,相应的机会主义 IM 具有一个通过优化 πₓₙ 获得的轮廓:

打开网易新闻 查看精彩图片

重要的是,与贝叶斯积分不同,IM 的优化操作确保了 πₓₙ 中固有的有效性性质被传递到 πₓₙᵐᵃʳᵍ,这意味着 IM 关于 Φ 的边缘推断不会陷入虚假置信。

上述应用于 IM 轮廓函数 (2) 的操作——用于获得上概率或消除冗余参数——均为 Choquet 积分的特例(例如 [74, 附录 C])。这些更一般的 Choquet 积分在正式决策情境中具有统计相关性 [45] 等。也就是说,若 ℓₐ(Θ) 表示当世界状态为 Θ 时采取行动 a 所关联的损失,则可能性 IM 通过 Choquet 积分,在给定 xⁿ 的前提下,提供了一种对行动 a 所关联风险的评估,该评估考虑了 Θ 的不确定性:

打开网易新闻 查看精彩图片

那么,例如,人们可能会选择行动 â(xⁿ),该行动使上述上期望损失最小化。关键在于,即使在 IM 的轮廓函数已经获得之后,仍需执行非平凡的操作——这为我们寻找尽可能简单且高效的轮廓近似方法提供了真实而实用的动力。

尽管 IM 的构造在概念上很简单,且其性质很强,但计算可能是一个挑战。问题在于,我们很少能在 Pθ 下获得相对似然 R(Xⁿ, θ) 的抽样分布的闭式表达,以方便精确计算 πₓₙ。因此,通常采用的策略是:在足够精细的 θ 网格上,使用蒙特卡洛方法在每个 θ 值处近似该抽样分布(例如,[32, 47])。也就是说,可能性轮廓被近似为:

打开网易新闻 查看精彩图片

其中,Xₘ,θⁿ 由 Pθ 生成的 n 个独立同分布样本构成,m = 1, ..., M。上述计算在单个或少数几个 θ 值上是可行的,但通常需要在覆盖(相关区域)的、往往为多维的参数空间 T 上进行足够精细的网格计算。例如,公式 (5) 中的置信集要求我们能求解方程 πₓₙ(θ) = α,或至少找出满足 πₓₙ(θ) ≥ α 的那些 θ 值;一种朴素的方法是在一个庞大的网格上计算轮廓,然后保留那些(近似地)满足前述方程的点。与此相关的计算复杂度为 O(Mgᵈ),其中 M 是蒙特卡洛样本量,d 是参数空间的维度,g 是参数空间 T 每一维上的网格点数;这导致大量计算被浪费。更一般而言,IM 输出的相关摘要涉及对轮廓函数的优化,而数值实现该优化需要多次评估轮廓函数。虽然在某些情况下可以对这种最朴素的方法进行简单调整(例如重要性抽样),但这些调整需要针对具体问题进行考虑,且不能期望在计算效率上带来实质性改进。这是一个严重的瓶颈,因此迫切需要新的、非朴素的计算策略。

  1. 基础类变分 IMs

3.1 设置

上述回顾的基于蒙特卡洛的策略并非数值近似可能性 IM 的唯一方法。另一种选择是基于现有大样本理论 [57] 的解析“高斯”近似(见下文)。此处的目标是在(或多或少精确但昂贵的)蒙特卡洛近似与(粗糙但廉价的)大样本近似之间取得平衡。为实现这一平衡,我们选择聚焦于可能性 IM 输出的一个特定特征,即公式 (5) 中的置信集 Cα(xⁿ),并选择一种近似方法,使其至少能精确匹配给定的置信集。我们的具体提案类似于目前在贝叶斯分析中广泛使用的变分近似:首先指定一个相对简单的候选概率分布族,然后通过寻找使该候选分布与精确后验分布之间的距离/散度(的上界)最小化的成员来获得近似解。我们的方法的不同之处在于,我们旨在通过(对适当选择的概率分布应用)概率到可能性的转换,来近似一个可能性测度。

根据 Destercke 和 Dubois [21]、Couso 等人 [17] 及其他学者的研究,可能性 IM 的信度集(credal set)(Π̄ₓₙ),即所有被 Π̄ₓₙ 支配的精确概率分布的集合,具有一个非常简单且直观的刻画:

为方便起见,我们将在下文中用下标“n”替代下标“xⁿ”——例如,用 Qₙ 和 Π̄ₙ 代替 Qₓₙ 和 Π̄ₓₙ——以简化符号表示。)也就是说,一个依赖于数据的概率分布 Qₙ 与 Π̄ₙ 一致,当且仅当对于每个 α ∈ [0, 1],它赋予 IM 的置信集 Cα(xⁿ)(见公式 (5))至少 1−α 的概率质量。此外,如果存在“最佳”的内概率近似,则该近似对应于一个 Qₙ*,使得对所有 α ∈ [0, 1],Qₙ*{Cα(xⁿ)} = 1−α。对于某一类特殊的统计模型,Martin [48] 证明了这种最佳内近似对应于费希尔的信念分布(fiducial distribution)和默认先验贝叶斯后验分布。但在这一特殊模型类别之外,尚不清楚如何找到最佳内近似。一个不那么雄心勃勃的目标是:对于固定的 α,寻找一个概率分布 Qₙ,α*,使得

打开网易新闻 查看精彩图片

3.2 提出的近似方法

我们建议从一个由通用参数空间 Ξ 索引的数据依赖型概率分布族 ᵛᵃʳ = {Qₙᵡ : ξ ∈ Ξ} 开始。一个重要的例子是 Qₙᵡ 为高斯分布的情形,其均值向量和/或协方差矩阵以某种特定方式依赖于(数据和)ξ。具体而言,由于可能性 IM 轮廓的峰值位于最大似然估计量 θ̂ₙ = θ̂ₓₙ 处,因此将高斯分布 Qₙᵡ 的均值向量固定在 θ̂ₙ 是合理的;而对于协方差矩阵,一个自然的选择是 ξ² Jₙ⁻¹,其中 Jₙ = Jₓₙ 是依赖于数据及所设定统计模型的观测费舍尔信息矩阵。

鉴于当 ξ = 1 时,Qₙᵡ 是对 Π̄ₙ 渐近最优的内概率近似(参见 [57]),该高斯族是一个非常合理且默认的选择。因此,我们的方案是在高斯近似中引入一些额外的灵活性,允许其离散程度根据 ξ > 1 或 ξ < 1 而扩展或收缩。尽管基于高斯的近似是自然的,但选择 ᵛᵃʳ 为高斯族并非唯一选项——参见下面的例 4。事实上,如果参数空间具有在通常欧几里得空间中不存在的结构(例如,若 T 是概率单纯形),则选择 ᵛᵃʳ 以尊重该结构是完全合理的。

对 ᵛᵃʳ 所施加的一个高层条件是:它必须足够灵活,即随着 ξ 的变化,Qₙ,αᵡ 对可能性 IM 的 α-截集的概率值能够小于或大于目标水平 1−α。上述高斯近似显然满足此条件,因为 ξ 控制着 Qₙ,αᵡ 的离散程度,通过取足够小或足够大的 ξ 值,可以使前述概率任意地变小或变大。这一温和条件对于几乎所有其他(合理的)近似族 ᵛᵃʳ 同样易于验证。

给定这样一个由参数 ξ ∈ Ξ 索引的合适的近似族 ᵛᵃʳ,我们提出的程序如下。定义一个目标函数:

打开网易新闻 查看精彩图片

其中 Θ₁ᵡ, ..., Θₖᵡ 独立同分布于 Qₓₙᵡ。可以想见,上述样本对每个 ξ 而言都是廉价的,因为分布族 已由用户指定;但我们仍需要 M 个蒙特卡洛样本来对每个 k 计算 πₓₙ(Θₖᵡ)。这导致计算复杂度为 O(MK),除了在低维情形下(例如 d ∈ {1,2}),该修改相比 (7) 中的朴素策略已带来显著的计算节省。接下来将描述一个额外的计算步骤,但其复杂度是有界的,因此整体复杂度仍保持在 O(MK)。

仅拥有目标函数的无偏估计量还需要对数值算法进行一些调整。特别是,我们不能使用假设函数值无噪声的牛顿-拉弗森(Newton–Raphson)算法,而必须使用一种适应于函数值存在噪声的随机近似算法(例如,[40,51,61,69,70])。基本的 Robbins–Monro 算法通过以下迭代更新来求解 (10) 的根:

打开网易新闻 查看精彩图片
打开网易新闻 查看精彩图片

我们所提近似方法的伪代码见算法 1。总结而言,我们提出了一种针对参数 Θ 的简单数据依赖型概率分布,其概率到可能性的轮廓在指定阈值 α 处至少能紧密匹配 IM 的轮廓。更具体地说,一个合理的概率分布选择是高斯分布,而算法 1 中所呈现的方法通过缩放协方差矩阵,使其对应的轮廓函数在阈值 α 处能精确逼近 IM 的轮廓。

在某些温和条件下,上述定义的序列 {ξ⁽ᵗ⁾ : t ≥ 0} 在概率意义上收敛于 fα 在 (10) 中的根。若 ξ̂ 是算法达到实际收敛(例如,当 |f̂α(ξ⁽ᵗ⁾)| 或变化量 |ξ⁽ᵗ⁺¹⁾ − ξ⁽ᵗ⁾| 小于某个指定阈值)时返回的值,则我们设 Q̂ₙ,α = Qₙ,αᵡ̂。该分布应是对 IM 可能性测度内近似(即 (Π̄ₓₙ) 中“最弥散”的成员)的一个合理准确的近似。因此,将 (2) 中的概率到可能性转换应用于 Q̂ₓₙ,α,应能合理准确地近似精确的可能性 IM 轮廓 πₓₙ,至少在它们的上 α-截集方面如此。下文各节中的示例证实了这一点,以及第 4 节中提出的更复杂近似方法确实具有合理的准确性。

如上所述,合适的 ᵛᵃʳ 选择取决于具体情境。在此选择中的一个重要考量是:能否对由内近似 Q̂ₙ,α 所确定的近似轮廓进行精确计算。对于前述正态变分族 ᵛᵃʳ(均值为 θ̂ₙ,协方差为 ξ̂² Jₙ⁻¹)而言,情况正是如此,因为

打开网易新闻 查看精彩图片

3.3 数值示例

我们在此的第一个目标是为所提出的近似方法提供一个概念验证。为此,我们呈现几个低维示例,以便能够同时可视化精确的 IM 轮廓和近似轮廓,并直接评估近似的质量。除下面的例 4 外,所有示例均使用上述描述的正态变分族 ,其均值为 θ̂ₙ、协方差为 ξ² Jₙ⁻¹(其中 ξ 待确定)。所有示例均展示基于 α = 0.1、M = 200 个蒙特卡洛样本、步长 wₜ = 2(1 + t)⁻¹ 以及收敛阈值 ε = 0.005 的类变分 IM 近似 Q̂ₙ,α。

例 1. 回顾第 2 节中探讨的二项分布例子,其中 Xⁿ 由 n 个独立同分布的 Ber(Θ) 随机变量组成。基于观测数据 Xⁿ = xⁿ 的精确 IM 可能性轮廓表达式见公式 (3),并在图 1(a) 中针对 n = 15 且 ∑ᵢ₌₁ⁿ xᵢ = 6 的情形予以展示。图 1(a) 同时也显示了所提出的基于高斯的变分近似对应的轮廓。请注意,这两个轮廓非常吻合,尤其是在专门针对的水平 α = 0.1 处。

例 2. 假设 Xⁿ 由独立同分布的二元正态配对组成,其均值为零、方差为 1,具有共同密度函数:

打开网易新闻 查看精彩图片
打开网易新闻 查看精彩图片

对未知相关系数 Θ 进行推断是一个出人意料地具有挑战性的问题(例如,[3,53,59])。事实上,尽管这是一个标量指数族模型,但它并没有一个一维的充分统计量;此外,还存在多种不同的辅助统计量可用于条件化,从而导致不同的解。图 1(b) 展示了基于 (2) 的朴素蒙特卡洛实现所得到的精确 IM 轮廓,数据为模拟样本,样本量 n = 50,真实相关系数为 0.5。精确轮廓具有一些不对称性,而正态近似无法完美容纳这种不对称性,但通过略微扩大上 0.1 水平集来弥补了这一不足。所提近似的计算效率与准确性的定量评估将在下面的例 6 中给出。

例 3. [28] 第 8.4 表中呈现的数据涉及氯乙酸暴露与小鼠死亡率之间的关系。拟合一个简单的逻辑回归模型,以将二元死亡指标 (y) 与氯乙酸暴露水平 (u) 关联起来,该数据集包含 n = 120 只小鼠。也就是说,Xⁿ 由独立配对 Xᵢ = (Uᵢ, Yᵢ) 组成,且在给定 Uᵢ 的条件下,Yᵢ 服从伯努利模型,其质量函数为:

其中 F(z) = (1 + e⁻ᶻ)⁻¹ 是逻辑分布函数。相应的似然函数无法以闭式最大化,但数值求解是常规操作。最大似然估计量及相应的观测信息矩阵可导出渐近有效的推断,这是标准统计软件包所报告的结果。然而,对于精确推断而言,计算负担更重:在足够精细的 θ 值网格上评估经验证的相对似然函数成本相当高。图 1(c) 展示了基于 (2) 的朴素蒙特卡洛实现所得到的回归系数的精确 IM 可能性轮廓的 0.1 水平集,并与所提出的变分近似并列展示。变分方法比朴素方法快近 2.5 倍,但两个轮廓几乎完全吻合。

打开网易新闻 查看精彩图片

其中,自然地,我们将均值取为最大似然估计量 θ̂ = n⁻¹X,精度取为 nξ(其中 ξ > 0 待确定)。当然,此处也可以使用高斯变分近似,但我们采用狄利克雷近似旨在突出我们方案的灵活性。图 1(d) 展示了基于 K = 3 和计数 X = (8, 10, 7) 的近似 IM 轮廓。精确的 IM 轮廓几乎不可能计算,因为朴素蒙特卡洛方法速度缓慢,当蒙特卡洛样本量过小时轮廓会呈现噪声,且数据的离散性质使其形状类似于图 1(a) 中的二项分布图。然而,在此处,我们仅需几秒钟即可获得一个平滑的轮廓近似。

本节的第二个目标是提供更深入的示例,以说明利用所提出的类变分近似 IM 可以开展何种分析。我们在针对计数数据的回归建模背景下进行这一说明。

例 5. 泊松对数线性模型被广泛用于分析基于计数的离散响应变量与一组固定解释变量之间的关系;即使解释变量并非由设计固定,也几乎总是假设其分布不依赖于任何相关参数;这使得解释变量成为辅助统计量,因此通常需以其观测值为条件。令 Xᵢ 表示第 i 个观测的响应变量,zᵢ₁, ..., zᵢₚ 表示该观测的 p 个解释变量,i = 1, ..., n。泊松对数线性模型假设在给定 zᵢ = (zᵢ₁, ..., zᵢₚ) 的条件下,Xᵢ 服从泊松分布,且跨 i = 1, ..., n 相互独立,其边际质量函数为:

打开网易新闻 查看精彩图片

考虑 [1] 第 3.2 表中提供的数据,这些数据来源于一项关于马蹄蟹筑巢习性的研究。在该研究中,每只 n = 173 只雌性马蹄蟹均有一只雄性附着在其巢穴旁,研究目标是探索影响附近是否存在额外雄性(称为“卫星雄性”)的因素。响应变量 X 是每只雌性蟹观测到的卫星雄性数量。此处,我们重点关注评估两个与雌性蟹体型相关的解释变量对这一响应的影响:z₁(体重,单位 kg)和 z₂(甲壳宽度,单位 cm)。

IM 轮廓的变分近似由公式 (12) 得到,其中 ξ 使用 M = 200 个蒙特卡洛样本估计,步长 wₜ = 2(1 + t)⁻¹,收敛阈值 ε = 0.005。或许首先应考虑的问题是:z₁ 或 z₂ 中至少有一个是否对 X 有影响?为回答此问题,构建了 (Θ₁, Θ₂) 的边际 IM,如图 2(a) 所示。值得注意的是,假设 “H: Θ₁ = Θ₂ = 0” 的上概率接近零,这提供了强有力的证据表明 Θ₁ 或 Θ₂ 至少有一个不为零。为分别评估 z₁ 和 z₂ 的影响,图 2(b,c) 展示了 Θ₁ 和 Θ₂ 各自的边际 IM。尽管有强证据支持 “Θ₁ > 0”,但假设 “Θ₂ = 0” 却非常合理。最后,鉴于存在明显证据支持 “Θ₁ > 0”,人们可能会进一步询问形如 “Hᵧ: Θ₁ > γ”(γ > 0)的哪些假设得到良好支持。这一问题可通过使用 Hᵧ 的边际必要性测度来解决,如图 2(d) 所示。我们可以看到,“Θ₁ > 0.1” 得到良好支持,表明每增加一公斤体重,雌性蟹平均拥有的卫星雄性数量大约增加 10%。重要的是,IM 的一致有效性保证了上述任何一条推断具有误导性的概率都可被控制在极小范围内。

打开网易新闻 查看精彩图片

为了突出 IM 的校准保障并提供与另一种框架的比较,我们按照上述设置进行了一项模拟研究。我们生成了 500 个数据集,每个数据集包含 n = 25 个观测值。所有数据集中,25 对解释变量保持固定,并通过从上述原始数据集的 173 对中随机抽取获得。重要的是,所选解释变量经过缩放,使得 ∑ⱼ₌₁ᵖ zᵢⱼ = 0 且 p⁻¹∑ⱼ₌₁ᵖ zᵢⱼ² = 1——这种缩放不影响解释变量之间的依赖关系,确保了各 Θⱼ 具有可比性,从而使如下所述的假设 H₃ 具有意义。响应变量 Xᵢ 独立地从 Pois(Λᵢ) 中抽样,其中

打开网易新闻 查看精彩图片
打开网易新闻 查看精彩图片
打开网易新闻 查看精彩图片

  1. 超越基础类变分 IMs

4.1 设置与方案

为加速计算,可考虑对第 3 节所述的基本程序进行多种调整。如果在第 3 节提出的方案中仍存在计算瓶颈,则意味着在随机近似算法的每次迭代中,IM 的可能性轮廓必须在 M 个点上进行评估。虽然在某些应用中(包括上文所呈现的例子),这并不算昂贵,但在其他应用中可能构成显著负担。另一个相关挑战是:迄今为止我们专注的变分族上的标量索引 ξ 可能不够灵活。通过考虑更高维的 ξ 来增加灵活性,同样会增加计算负担,因此需要谨慎处理。在此,我们的目标是同时解决上述两个挑战。

我们在此考虑的特定修改最适合于以下情形:变分族 ᵛᵃʳ 满足如下性质——对于每个 ξ ∈ Ξ,对应 Qₙᵡ 的 100(1−α)% 可信集可以以闭式表达(或至少能够简洁概括)。此处将要呈现的想法更具一般性,但为使细节尽可能简单具体,我们将重点放在高斯变分族的情形上。此时,可信集是 d 维空间中的椭球体。

作为对第 3 节引入的标量 ξ 索引高斯族的推广,令 ξ ∈ Ξ = ℝᵈ₊₀ 为一个 d 维向量索引,并取 ᵛᵃʳ 为具有均值向量 θ̂ₙ(即最大似然估计量)和协方差矩阵 Jₙ(ξ)⁻¹ 的 d 维高斯族,其定义如下:对观测费舍尔信息矩阵 Jₙ 进行特征分解,记作 Jₙ = UΨUᵀ,然后设定

打开网易新闻 查看精彩图片
打开网易新闻 查看精彩图片

4.2 数值示例

我们将通过三个例子来说明这一新版高斯变分族及其近似算法。第一个例子重新审视前述的二元正态相关系数案例,但提供关于该近似方法计算时间和精度的更具体细节;第二个例子是经典的双参数伽马模型;第三个例子是一个涉及惩罚项的相对高维模型,旨在作为进入高维问题 IM 解法的入口。

例 6. 此处我们重新审视上文例 2 中的二元正态相关系数示例,但这次对所提近似方法(即算法 2 所述版本)与朴素蒙特卡洛轮廓评估(如公式 (7) 所述)在计算时间和精度方面进行定量比较。我们沿用例 2 的设置,数据集从真实相关系数 Θ = 0.5 的标准二元正态分布中生成,样本量 n 各不相同。在此情形下,我们生成 100 个不同大小(n = 50、100 和 200)的数据集,并针对每个数据集,分别基于公式 (7) 中的朴素策略和新提出的近似方法评估 IM 轮廓函数,记为 π̂ⁿᵃⁱᵛᵉ 和 π̂ᵃᵖᵖʳᵒˣ。两种方法在评估轮廓时均基于 M = 500 个蒙特卡洛样本。在表 1 中,我们比较了相对计算时间——定义为朴素策略耗时除以所提策略耗时——以及两者之间的 L₁ 距离 ∫ |π̂ᵃᵖᵖʳᵒˣ(θ) − π̂ⁿᵃⁱᵛᵉ(θ)| dθ;这些数值是在每种样本量 n 下的 100 个数据集上取平均得到的。这里的解释是:朴素策略是“黄金标准”,因为它在每个网格点上都能准确评估 IM 轮廓。因此,理想情况下,π̂ᵃᵖᵖʳᵒˣ 与 π̂ⁿᵃⁱᵛᵉ 之间的距离应很小,而前者解法应比后者的暴力穷举法更高效。然而,随着样本量 n 增加,我们知道两点:

  • 两种策略的计算复杂度均随 n 线性增长,但朴素策略的增长速率快于所提近似方法,因此计算时间比率应大致恒定且大于 1;
  • 根据文献 [57] 的结果,π̂ᵃᵖᵖʳᵒˣ 中采用的高斯近似精度会提高,因此我们预期两个轮廓将趋于一致。

表 1 的结果证实了这些预期,即相对计算时间稳定在数值 2 附近,表明朴素策略的计算时间大约是所提策略的两倍,且两个轮廓之间的 L₁ 距离随 n 增大而减小。在这种情况下,朴素策略在计算时间对比上占优,因为其仅涉及一个标量参数,但性能差异仍相当显著。

打开网易新闻 查看精彩图片

例 7. 假设 Xⁿ 是来自形状参数 Θ₁ 和尺度参数 Θ₂ 的伽马分布的一个大小为 n 的独立同分布样本。我们模拟了 n = 25 的数据,其中 Θ₁ = 7 且 Θ₂ = 3,并在图 4(a) 中绘制了 (Θ₁, Θ₂) 的近似 IM 轮廓。该轮廓的构建方法是:首先构建 (log Θ₁, log Θ₂) 的高斯近似轮廓,然后将其映射回 (Θ₁, Θ₂) 空间。当映射到对数尺度时,参数的非负约束被消除,从而提高了高斯近似的质量;若直接应用于 (Θ₁, Θ₂) 空间,则近似效果较差。作为对比,图中也展示了相对似然的轮廓,其与高斯轮廓的相似性表明后者是对精确 IM 轮廓的良好近似——尽管精确计算成本高昂。事实上,[42] 中的例 1 考虑了完全相同的模拟设定,他也得到了一个类似图 4(a) 中的香蕉形轮廓。

我们还重复上述例子 1000 次,图 4(b) 展示了随机变量 πₓₙ(Θ) 的分布函数,分别基于精确轮廓和高斯近似。结果显示,基于精确轮廓的分布函数为 Unif(0,1)(忽略蒙特卡洛抽样变异),而基于高斯近似的分布函数在整个范围内与 Unif(0,1) 在经验上无法区分。

例 8. 最简单且最典型的高维推断例子是“多正态均值问题”(many-normal-means problem),可追溯至经典论文如 [68]、[34]、[9] 等。该模型假设数据 Xⁿ 包含 n 个独立但非同分布的观测值,其中 Xᵢ ~ N(Θᵢ, σ²),σ² 被假定已知,而向量 Θ = (Θ₁, ..., Θₙ) 未知且待推断。粗略地说,上述文献的核心观点是:在平方误差损失下,最佳无偏估计量 θ̂ = Xⁿ(也是最大似然估计量)是不可容许的。这一结果激发了对惩罚估计的研究努力,包括如今著名的 lasso 方法(例如,[71,72])。遵循这一思路,并秉承第 5 节示例的精神,我们此处提出使用一种相对惩罚似然函数。

打开网易新闻 查看精彩图片
打开网易新闻 查看精彩图片

尽管相对惩罚似然函数可以闭式表达,但仍需蒙特卡洛方法来评估该轮廓。当维度 n 甚至只是中等大小时,在足够精细的 n 维网格上执行这些计算以获得 Θ 的置信集在实际操作中几乎是不可能的。但上述描述的变分近似提供了一种比朴素蒙特卡洛方法计算效率更高的替代方案,能够处理中等到较大的 n。

在此,我们遵循前述策略,即采用高斯近似,其均值 θ̂ 等于 lasso 或最大惩罚似然估计量,协方差矩阵为 n×n 阶的 Jₙ(ξ)⁻¹,由 n 维向量 ξ 索引;在此情形下,初始的 Jₙ 取为无惩罚信息矩阵 Jₙ = σ⁻²Iₙ,与单位矩阵成比例。此处的直觉是,坐标特定的调整因子 ξᵢ 将允许高斯近似在某种程度上适应真实信号 Θ 中的稀疏性。为说明这一点,我们考虑 n = 50,且真实参数 Θ 的前五个分量等于 5,其余 45 个分量等于 0,即 X 中仅 10% 的坐标包含信号,其余 90% 仅为噪声。我们还固定 α = 0.1 用于近似。对于单个数据集,图 5(a) 展示了所提随机近似更新收敛后得到的估计值 ξ̂。黑色点对应信号(非零真实均值),灰色点对应噪声。关键观察结果是:对应信号的 ξ 值往往大于对应噪声的 ξ 值;在信号情形下几乎没有变异性,但在噪声情形下存在显著变异性。在噪声情形下 ξ 值趋于更小是符合预期的,因为在那些明显为 0 的均值周围,IM 的可能性轮廓需要更小的离散程度。我们将上述模拟重复 1000 次,并绘制了在真实 Θ 处的精确(使用朴素蒙特卡洛)和近似(使用高斯变分族)IM 轮廓的分布函数。再次地,精确轮廓在 Θ 处服从 Unif(0,1) 分布,结果如图 5(b) 所示。高斯近似仅被设计为在水平 α = 0.1 处校准,这一点显然已实现;但在较低水平上略显激进,在较高水平上则偏保守。关于该方案在高维问题中的进一步研究将在其他地方报告。

  1. 冗余参数问题

5.1 参数情形

上述视角假设存在一个未知的模型参数 Θ,且主要目标是基于观测数据 Xⁿ = xⁿ 从模型中对 Θ 的整体不确定性进行量化。当然,对 Θ 的不确定性量化意味着对任意特征 Φ = g(Θ) 的不确定性量化(如第 2 节所述)。然而,如果唯一目标是对某个特定特征 Φ = g(Θ) 的不确定性进行量化,那么很自然会问:我们能否做得比先量化关于完整 Θ 的不确定性、再推导出 Φ 的相应结果更好?存在提高效率的机会,但这需要消除冗余参数——即在某种意义上与 Φ 互补或正交的 Θ 的那些方面。一种相当通用的消除冗余参数的策略是剖面化(profiling)(例如,[49,63,66]),如下所述。

或许并不令人意外的是,尽管在第 2 节所呈现的 IM 构造中使用相对似然函数是非常自然且在某种意义上“最优”的,但它并非唯一选项。对于涉及冗余参数的情形,一种策略是将公式 (2) 中的相对似然函数替换为一个替代物,即相对剖面似然函数

打开网易新闻 查看精彩图片

这种构造的优势在于,它通常比第 2 节中所介绍的朴素冗余参数消除方法更高效;参见例如 [57]。上述外层上确界(supremum)的出现是因为 Φ 并非所设定模型的完整参数;更多细节请参见 [49, 第 3.2 节]。通常情况下,上述右侧的概率在 θ 上近似为常数(当 g(θ) = φ 时),但不能依赖于此——为了确保 IM 的强有效性性质成立,不幸的是,必须计算该上确界。

为具体说明,我们将聚焦于一个看似简单却颇具挑战性的问题,即对伽马分布均值进行高效推断。粗略而言,伽马均值是形状和尺度参数的一个高度非线性函数,这使得经典的一阶渐近近似在有限样本下表现较差。因此,伽马均值问题受到了相当多的关注,重点在于推导具有更高阶精度的渐近近似;我们建议读者参考 [27] 以获取更多细节。Martin 和 Liu [55] 提出了伽马均值问题的一个精确 IM 解法,而最近,在 [49, 例 6] 中提出了一种基于剖面化的可能性 IM 解法,并被证明优于各种现有方法。此处我们的重点是在这一新背景下展示变分近似的质量。

例 9. 设伽马模型由参数 θ = (θ₁, θ₂) 索引,其中 θ₁ 和 θ₂ 分别代表(正的)形状参数和尺度参数。在伽马模型中,最大似然估计量 θ̂₁ 和 θ̂₂ 没有闭式表达式,但可以数值最大化似然函数来求得;也可以数值或解析地获得观测费舍尔信息矩阵 J。对于剖面似然函数,将模型重新参数化为以均值参数 Φ = θ₁θ₂ 和形状参数 θ₁ 表示可能有助于简化问题。将密度函数用这种新参数化形式表示为:

打开网易新闻 查看精彩图片

在这种形式下,对于任何固定的 φ,可以数值最大化基于数据 Xⁿ 的似然函数以求得 θ₁,从而得到(相对)剖面似然函数。

Fraser 等人 [27] 提供了一个例子,其中 n = 20 只小鼠暴露于 240 拉德的伽马辐射,并记录了它们的存活时间。图 6 中展示了精确的基于剖面化的边际可能性 IM 轮廓(黑线)。该计算相对昂贵,因为在网格上的每个 φ 点处,我们的蒙特卡洛近似都需要针对不同的 θ₁ 值进行优化。作为对比,我们考虑一个高斯可能性轮廓,其均值为 φ̂ = θ̂₁θ̂₂,方差为 ξ² ĝ(θ̂)ᵀJ⁻¹ĝ(θ̂),其中 g(θ) = θ₁θ₂,梯度 ĝ(θ) = (θ₂, θ₁)ᵀ。图 6 展示了 ξ = 1 的高斯近似(如 [57] 中所讨论)以及根据第 4 节中变分近似确定的 ξ = 1.28。该近似仅需不到一秒即可获得,我们发现,正如预期的那样,它在目标水平 α = 0.1 处(右侧粗线部分)与精确轮廓紧密匹配,而在左侧(细线部分)则略显保守。显然,基本的大样本高斯近似在右尾部过窄,这证实了上述观点:一阶渐近理论在小样本情况下提供的近似效果相对较差。另一方面,我们的变分近似能够适当调整,在某些地方匹配精确轮廓,而在其他地方则略显谨慎或保守。

打开网易新闻 查看精彩图片

5.2 非参数情形

非参数问题是指,底层分布 P 并未被假定为具有由有限维参数索引的特定形式。在某些应用中,感兴趣的量本身就是分布 P 本身(或例如,其密度函数);而在另一些情况下,感兴趣的则是 P 的某个(有限维的)特征或泛函 Θ。我们此处的关注点在于后一种情况,因此它也符合涉及冗余参数问题的一般框架,因为一旦考虑了 Θ,P 中剩余的部分将被视为“冗余”并需被消除。

至少在原则上,可以采用类似于上述参数情形的方法来处理该问题,即通过剖面化去除 P 中的冗余部分。回顾一下,剖面化的目标是降低维度,以便可以直接评估数据与感兴趣量候选值之间的相容性。由于通常情况下,感兴趣的量 Θ 具有某种现实世界的解释,因此有机会利用这种解释进行相容性评估,而无需进行剖面化。这正是 [13] 所采取的方法,该方法建立在 M-估计的经典工作(例如,[33,67])以及关于 Gibbs 后验的较不经典的工作(例如,[7,29,56,78])基础之上,我们在下文简要总结。

设数据 Xⁿ = (X₁, ..., Xₙ) 由独立同分布的分量组成,其中 Xᵢ ~ P,且对 P 本身一无所知或未作任何假设。在这种更一般的情形下,未知的兴趣量 Θ = Θ(P) 是底层分布的一个泛函。例子包括 P 的分位数和矩。假设 Θ 可以表示为某个风险或期望损失函数的最小化者。也就是说,假设存在一个损失函数 (x, θ) ↦ LOSSθ(x),使得

打开网易新闻 查看精彩图片

外层上确界(supremum),类似于公式 (15) 中的情形,是对所有那些使 P 的相关特征 Θ 取值为 θ 的概率分布 P 进行最大化。该上确界的出现是因为 Rᴱᴿ(Xⁿ, θ) 的分布显然依赖于底层的 P,但 P 是未知的。这使得基于朴素蒙特卡洛方法直接评估 IM 轮廓变得不可行。幸运的是,有效性仅要求在唯一真实的 P 下满足特定校准条件,这提供了一条捷径。Cella 和 Martin [13] 提出将“对所有与 θ 相容的 P 进行独立同分布抽样”替换为从经验分布中进行独立同分布抽样——后者是真实 P 的一个良好估计量。这相当于使用自助法(bootstrap)(例如,[18,23,25])来近似上述轮廓,且 Cella 和 Martin 证明了相应的 IM 是渐近有效的。在此,我们将展示所提出的类变分 IMs 可以为这种基于自助法的轮廓提供良好的近似。

例 10. 假设我们感兴趣的是分布 P 的 τ 阶分位数,即精确点 Θ = Θ⁽ᵗ⁾,使得 P(X ≤ Θ⁽ᵗ⁾) = τ,其中 τ ∈ (0,1)。上述非参数 IM 构造中的关键组成部分是选择一个合适的损失函数。对于分位数估计,众所周知,损失函数由下式给出:

打开网易新闻 查看精彩图片

其中,p̂ 表示基于观测数据对密度 p 的核密度估计。与第 3.3 节中相同的设置被应用,其中 ξ 使用 M = 200 个蒙特卡洛样本估计,步长 wₜ = 2(1 + t)⁻¹,α = 0.1,收敛阈值 ε = 0.005。请注意,变分近似在除左侧小 α 值外的所有地方都是完美的,在左侧它略显保守。

为验证变分方法在非参数设定下能提供近似的有效性,我们通过重复上述情景 250 次进行了一项模拟研究。对于每个数据集,近似轮廓在 Θ = 2.53 处进行评估,这大致对应于当 P 服从 Gamma(4,1) 分布时的第一四分位数。该轮廓的经验分布如图 7(b) 所示,表明近似的有效性确实已实现。

5.3 半参数情形

前两小节所描述的参数与非参数情形之间的中间地带被称为半参数问题,即同时包含参数和非参数部分的问题。或许最简单的例子是误差分布未指定的线性回归模型:线性均值函数是参数部分,而误差分布是非参数部分。下面我们将聚焦于带删失数据的半参数模型,但当然其他例子也是可能的;更多细节请参见例如 [5]、[75] 和 [38]。

上文第 5.1 节所述的相同剖面化策略也可应用于半参数模型;[58] 是一个重要参考文献。为具体说明,我们将考虑一个涉及删失数据的常见情形。也就是说,假设我们正在测量土壤中某种特定化学物质的浓度,但我们的测量仪器有一个较低的检测限——即低于该限值的浓度无法被检测到。在这种情况下,浓度数据(左)被删失。我们心中可能有一个针对测量浓度的参数模型,但删失会破坏数据并最终改变该模型。令 Yᵢ 表示站点 i 处的实际化学物质浓度,其值我们可能观察到也可能观察不到;Yᵢ 被赋予一个统计模型 {Pθ : θ ∈ T},而该模型参数的真实但未知值 Θ 需要被推断。令 Cᵢ 表示删失水平,我们假设——不失一般性——其受抽样变异影响,即 Cᵢ 是随机变量。那么观测数据 Xⁿ 由独立同分布的配对 Xᵢ = (Zᵢ, Tᵢ) 组成,其中

打开网易新闻 查看精彩图片

这取决于真实未知模型参数 Θ 的通用值 θ(用于浓度)以及真实未知删失水平分布 G 的通用值 G。在上述表达式中,g 和 pθ 分别是删失分布和浓度分布的密度函数,而 G 和 Pθ 则是相应的累积分布函数。现在应该清楚为何这是一个半参数模型:除了明显的参数模型外,还有一个针对缺失水平的非参数模型。

该半参数模型的一个显著特征是,似然函数是“可分离的”,即它是涉及 θ 的项与涉及 G 的项的乘积。因此,如果我们对 G 进行优化,然后构造相对剖面似然比,则涉及 G 优化的部分将被消去。这意味着我们可以简单地忽略涉及 G 的部分,而直接使用如下形式的相对剖面似然函数:

打开网易新闻 查看精彩图片

通过蒙特卡洛方法评估右侧表达式,归结为从 Ĝ 中抽样删失水平、从 Pθ 中抽样浓度水平,然后根据公式 (17) 构造新的数据集。虽然该过程在概念上相对简单,但在足够精细的 θ 值网格上进行朴素实现却相当昂贵。幸运的是,我们第 4 节提出的类变分近似方法可直接应用,能够快速生成一个闭式近似轮廓。

例 11. 为说明这一点,我们使用从内布拉斯加州一口井收集的阿特拉津(Atrazine)浓度数据。该数据包含 n = 24 个观测值,这些观测值如上所述受到随机左删失的影响。这是一个相当极端的情况,其中近一半(11 个)的 24 个观测值被删失,但先前的研究表明,对阿特拉津浓度采用对数正态模型是合适的 [30]。在环境科学应用中,对数正态分布常被用于建模左删失数据(例如,[39])。对数正态模型的密度函数为:

打开网易新闻 查看精彩图片

其中 θ = (θ₁, θ₂) 表示 log Y 的均值和方差参数。同样,对数正态模型仅用于描述观测到的浓度——对于删失观测值未作任何模型假设。图 8(a) 显示了通过应用 Kaplan–Meier 估计量(并将删失标签取反:tᵢ ↦ 1 − tᵢ)所获得的删失数据分布的非参数估计量 Ĝ。该 Ĝ 随后被用于定义我们此处所称的“精确”IM 轮廓(通过公式 (18)),然后相应的高斯变分近似——首先应用于 (θ₁, log θ₂),再映射回 (θ₁, θ₂)——如图 8(b) 所示。该图与 [11] 中图 10 所示的基于朴素蒙特卡洛方法的结果非常相似,但在计算上却远不那么昂贵。

打开网易新闻 查看精彩图片

  1. 结论

秉承当前贝叶斯统计中广泛使用的变分近似思想,并基于文献 [35] 中提出的最新思路,本文发展了一种策略,利用普通的蒙特卡洛抽样与随机近似方法,来近似可能性推断模型(possibilistic IM)的轮廓函数——或至少近似其在指定显著性水平 α 下的 α-截集(即水平集)。我们展示了一系列应用场景,从简单的教科书式问题,到涉及冗余参数的(参数、非参数和半参数)问题,甚至包括一个相对高维的问题,以突显所提方法的灵活性、准确性以及整体适用性。

当然,本文所提出的方法也存在若干局限性。这些局限自然引出了一些开放性问题和未来的研究方向。

第一,本文提出的更复杂且更高效的近似方法(即第 4 节所述方法)专门针对高斯变分族设计。这在实践中并非严重限制,因为在样本量 n 为中等或较大时,高斯分布通常能提供良好的近似 [57]。然而,肯定还存在其他变分族,其密度水平集在模型参数下具有简洁、近乎闭式的表达形式。识别适用于此类高效近似的其他合适模型,将为用户提供更多灵活性,并最终带来更好、更精确的近似结果。

第二,所提出的近似方法依赖于在初始阶段预先指定一个 α 值,这意味着我们实际上仅近似了可能性 IM 轮廓的某些特定特征。然而,IM 推断本质上具有整体性(holistic),这表明人们更倾向于一种广义而非高度特化的近似。因此,如何将这些针对不同 α 的 IM 近似“拼接”起来,构成一个统一的整体近似,是一个重要的开放问题。受文献 [35] 最新进展的启发,我们相信这一问题的答案是肯定的,相关细节将在其他地方报告。

第三,本文所呈现的具体方法聚焦于对未知参数 Θ 无先验信息(vacuous prior)的情形。近期研究(例如 [47])已开始探索如何将不完整或部分先验信息整合到可能性 IM 的构建中;随着我们向更高维问题扩展,这类方法几乎肯定是必要的。但引入部分先验信息的一个缺点是:IM 轮廓的评估通常比本文所考虑的无先验情形更加复杂。这种额外的复杂性意味着,在部分先验设定下,高效的数值近似方法变得更加重要。幸运的是,我们预期本文所提出的方案可近乎直接地推广至这些情形。

最后,除例 8 外,本文主要关注未知参数维度相对较低的问题。如常理所示,将所提近似策略扩展至高维场景必然面临挑战。目前我们尚不清楚如何克服这些挑战,但受本文所展示的最新进展的鼓舞,我们有信心这些问题终将被解决。

原文链接: https://www.sciencedirect.com/science/article/pii/S0888613X25001471?via%3Dihub