变分贝叶斯方法:修订间差异

来自testwiki
跳转到导航 跳转到搜索
imported>InternetArchiveBot
补救8个来源,并将0个来源标记为失效。) #IABot (v2.0.9.5
 
(没有差异)

2024年7月11日 (四) 10:16的最新版本

Template:贝叶斯统计 Template:Proofreader needed 变分贝叶斯方法是近似贝叶斯推理机器学习中较难积分的一系列技术,通常用于由观测变量及未知参数潜变量组成的复杂统计模型,之间存在的各种关系可用图模式描述。贝叶斯推理中,参数和潜变量一般归为“未观察变量”。变分贝叶斯方法主要用于两个目的:

  1. 提供未观测变量后验概率的分析近似,以便对变量统计推断
  2. 得出观测数据的边缘似然下界(即给定模型的数据的边缘概率,边缘化未观测变量)。通常用于模型选择,一般認為给定模型的边缘似然越高,拟合效果越好,产生数据的可能性越大(另见贝叶斯因子)。

就目的1(近似后验概率),变分贝叶斯是蒙特卡洛法(特别是马尔可夫链蒙特卡洛,如吉布斯采样)的替代,可采用完全贝叶斯的方法,对难以取值或抽样的复杂分布进行统计推断。蒙特卡洛技术利用一组样本对精确后验值进行数值近似,而变分贝叶斯可求近似后验的局部最优的精确解析解。

变分贝叶斯可视作最大期望算法(EM)的推广,从估计参数的单一最似然值的最大后验概率(MAP估计),推广到计算参数及潜变量的整个近似后验分布的完全贝叶斯估计。贝叶斯估计同EM一样,也能找到一组最优参数值,且与EM有相同的交替结构——基于不能求解析解的相依方程組。

很多应用中,变分贝叶斯能更快得到与吉布斯采样精度相当的解,但推导迭代方程组需要更多工作。即便是很多概念上非常简单的模型也如此,下面以只有2个参数、没有潜变量的基本不分层模型为例说明。

数学推导

问题

变分推断中,数据𝐗的一组未观测变量𝐙={Z1Zn}的后验分布近似于所谓变分分布Q(𝐙)

P(𝐙𝐗)Q(𝐙).

分布Q(𝐙)属于在形式上比P(𝐙𝐗)简单的分布族(如高斯分布族),这是为使Q(𝐙)与真实后验P(𝐙𝐗)更相似。

相似性用不相似函数d(Q;P)衡量,推断(inference)通过求分布Q(𝐙)使d(Q;P)最小化进行。

KL散度

最常见的变分贝叶斯法以QPKL散度作为不相似函数,很适合最小化。KL散度的定义是

DKL(QP)𝐙Q(𝐙)logQ(𝐙)P(𝐙𝐗).

注意QP是相反的。反向KL散度在概念上类似于期望最大化算法(KL散度则产生期望传播算法)。

难解性

变分技术通常用于近似以下参数:

P(𝐙𝐗)=P(𝐗𝐙)P(𝐙)P(𝐗)=P(𝐗𝐙)P(𝐙)𝐙P(𝐗,𝐙)d𝐙

𝐙边缘化(以计算分母中的P(𝐗))通常是难解的,因为𝐙的搜索空间在组合上可能很大。因此,我们用Q(𝐙)P(𝐙𝐗)求近似。

证据下界

Template:Main 考虑到P(𝐙𝐗)=P(𝐗,𝐙)P(𝐗),KL散度也可写作

DKL(QP)=𝐙Q(𝐙)[logQ(𝐙)P(𝐙,𝐗)+logP(𝐗)]=𝐙Q(𝐙)[logQ(𝐙)logP(𝐙,𝐗)]+𝐙Q(𝐙)[logP(𝐗)]

因为Q(𝐙)是分布,所以𝐙Q(𝐙)=1;又由于P(𝐗)对于𝐙是常数,所以

DKL(QP)=𝐙Q(𝐙)[logQ(𝐙)logP(𝐙,𝐗)]+logP(𝐗)

根据(离散随机变量期望的定义,上式可以写作

DKL(QP)=𝔼𝐐[logQ(𝐙)logP(𝐙,𝐗)]+logP(𝐗)

重排为

logP(𝐗)=DKL(QP)𝔼𝐐[logQ(𝐙)logP(𝐙,𝐗)]=DKL(QP)+(Q)

由于对数证据logP(𝐗)对于Q是定值,最大化末项(Q)会使QP的KL散度最小化。适当选择Q可以让(Q)更容易计算,更容易最大化,因此既有后验P(𝐙𝐗)的解析近似Q,也有对数证据logP(𝐗)的下界(Q)(KL散度非负)。

下界(Q)热力学自由能类似,也称作(负)变分自由能,因为它也可以表为负能量EQ[logP(𝐙,𝐗)]Q(Q)项也称作证据下界(ELBO),是数据的对数证据的下界。

证明

根据布雷格曼散度(KL散度是其特例)的广义勾股定理,可以证明[1][2]

布雷格曼散度的广义勾股定理[2]
DKL(QP)DKL(QQ*)+DKL(Q*P),Q*𝒞

其中𝒞是凸集,且

Q=Q*argminQ𝒞DKL(QP).

时等号成立。这时,全局最小值Q*(𝐙)=q*(𝐙1𝐙2)q*(𝐙2)=q*(𝐙2𝐙1)q*(𝐙1),其中𝐙={𝐙𝟏,𝐙𝟐}可以如下求得:[1]

q*(𝐙2)=P(𝐗)ζ(𝐗)P(𝐙2𝐗)exp(DKL(q*(𝐙1𝐙2)P(𝐙1𝐙2,𝐗)))=1ζ(𝐗)exp𝔼q*(𝐙1𝐙2)(logP(𝐙,𝐗)q*(𝐙1𝐙2)),

其中归一化常数为

ζ(𝐗)=P(𝐗)𝐙2P(𝐙2𝐗)exp(DKL(q*(𝐙1𝐙2)P(𝐙1𝐙2,𝐗)))=𝐙2exp𝔼q*(𝐙1𝐙2)(logP(𝐙,𝐗)q*(𝐙1𝐙2)).

ζ(𝐗)项在实践中常称作证据下界(ELBO),因为P(𝐗)ζ(𝐗)=exp((Q*))[1]如上所示。

交换𝐙1, 𝐙2的角色,可以迭代计算真实模型边际P(𝐙1𝐗), P(𝐙2𝐗)的近似q*(𝐙1), q*(𝐙2)。这种迭代保证单调收敛,[1]但收敛到的Q*只是DKL(QP)的局部极小值。

若约束空间𝒞被限制在独立空间内,即q*(𝐙1𝐙2)=q*(𝐙𝟏),则上述迭代将成为所谓平均场近似(mean field approximation)Q*(𝐙)=q*(𝐙1)q*(𝐙2),

平均场近似

变分分布Q(𝐙)常被假定为潜变量的某划分的因式分解,即对潜变量𝐙的划分𝐙1𝐙M

Q(𝐙)=i=1Mqi(𝐙i𝐗)

变分法可以证明,对每个因子qj的最佳分布qj*(就KL散度最小的分布而言,如上所述)满足[3]

qj*(𝐙j𝐗)=eEqj*[lnp(𝐙,𝐗)]eEqj*[lnp(𝐙,𝐗)]d𝐙j

其中Eqj*[lnp(𝐙,𝐗)]是数据与潜变量的对数联合概率期望值,是相对于不属于划分的所有变量的q*而取。关于qj*(𝐙j𝐗)的推导见[4]的引理 4.1。

实践中,通常用对数表示:

lnqj*(𝐙j𝐗)=Eqj*[lnp(𝐙,𝐗)]+constant

上式中的常数与归一化常数qj*上述表达式中的分母)有关,通常通过检查来恢复,因为表达式的其余部分通常是已知的分布(如正态分布伽马分布等)。

利用期望的性质,式Eqj*[lnp(𝐙,𝐗)]通常可简为潜变量先验分布的固定超参数及不属于当前划分的潜变量(即不属于𝐙j的潜变量)的期望值的函数(有时还包括方差等高阶)。这就在某一划分的变量分布参数同其他划分的变量的期望值之间产生了循环依赖,自然就需要类似期望最大化算法迭代法,以某种方式(也许随机)初始化潜变量期望(可能还有高阶矩),再利用当前期望依次计算每个分布的参数、适当设置新算出的分布的期望。这种算法保证收敛。[5]

换一种说法,对每个变量划分,简化其中变量分布的表达式、研究分布与相关变量的函数依赖关系,通常就能确定分布的族(反过来确定了常数)。分布的参数公式可用先验分布的超参数(已知常数)表示,也可用其他划分中变量函数的期望表示。通常这些期望可简为变量本身的期望值函数(即平均值);有时也会出现变量平方(与变量方差有关)或更高次幂()高阶)的期望。其他变量的分布一般来自已知族,相关期望的公式可查到;但公式取决于分布的参数,参数又取决于其他变量的期望,所以每个变量的分布参数公式都可表为一群变量间非线性相互依赖的方程。通常不能直接求解,不过可以用简单的迭代算法,大多数时候都能保证收敛。

变分推理对偶公式

用对偶公式图解坐标上升变分推理算法[4]

下面的定理称作变分推理对偶公式,[4]解释了变分贝叶斯方法中变分分布的一些重要性质。

Theorem 考虑两概率空间(Θ,,P), (Θ,,Q),其中QP。假设存在共同的主导概率测度λ使Pλ, Qλ。令h(Θ,,P)上的任意实值随机变量,满足hL1(P),则下列等式成立

logEP[exph]=supQP{EQ[h]DKL(QP)}.

此外,当且仅当(supremum)

q(θ)p(θ)=exph(θ)EP[exph],

关于概率测度Q几乎是确定的,其中p(θ)=dP/dλ, q(θ)=dQ/dλ分别表示概率测度PQλ的Radon–Nikodym导数,右式才取上界。

基本例子

考虑简单的不分层贝叶斯模型,由来自正态分布均值方差未知)的独立同分布观测量组成。[6]下面我们将详细介绍这模型,以说明变分贝叶斯方法的工作原理。

为了数学上的方便,下例中我们用精度(即方差倒数;多元正态分布时是协方差矩阵的逆。理论上精度和方差等价,因为两者间是双射)。

数学模型

共轭先验分布置于未知均值μ、精度τ,即均值也遵循高斯分布,而精度遵循伽马分布

τGamma(a0,b0)μ|τ𝒩(μ0,(λ0τ)1){x1,,xN}𝒩(μ,τ1)N=number of data points

先验分布超参数μ0,λ0,a0b0是已知定值,可设为小正数,以得到较宽的先验分布,表明我们对μτ的先验分布一无所知。

已知N个数据点𝐗={x1,,xN},而我们的目标是推断参数μ, τ后验分布q(μ,τ)=p(μ,τx1,,xN)

联合概率

所有变量的联合概率可重写为

p(𝐗,μ,τ)=p(𝐗μ,τ)p(μτ)p(τ)

个体因子是

p(𝐗μ,τ)=n=1N𝒩(xnμ,τ1)p(μτ)=𝒩(μμ0,(λ0τ)1)p(τ)=Gamma(τa0,b0)

其中

𝒩(xμ,σ2)=12πσ2e(xμ)22σ2Gamma(τa,b)=1Γ(a)baτa1ebτ

因式分解近似

假设q(μ,τ)=q(μ)q(τ),即后验分布因式分解为μ, τ的独立因子。这种假设是变分贝叶斯方法的基础。事实上,真正的后验分布并非如此(这种简单情形我们知道它是正态-伽马分布),因此所得结果将是近似值。

lnqμ*(μ)=Eτ[lnp(𝐗μ,τ)+lnp(μτ)+lnp(τ)]+C=Eτ[lnp(𝐗μ,τ)]+Eτ[lnp(μτ)]+Eτ[lnp(τ)]+C=Eτ[lnn=1N𝒩(xnμ,τ1)]+Eτ[ln𝒩(μμ0,(λ0τ)1)]+C2=Eτ[lnn=1Nτ2πe(xnμ)2τ2]+Eτ[lnλ0τ2πe(μμ0)2λ0τ2]+C2=Eτ[n=1N(12(lnτln2π)(xnμ)2τ2)]+Eτ[12(lnλ0+lnτln2π)(μμ0)2λ0τ2]+C2=Eτ[n=1N(xnμ)2τ2]+Eτ[(μμ0)2λ0τ2]+Eτ[n=1N12(lnτln2π)]+Eτ[12(lnλ0+lnτln2π)]+C2=Eτ[n=1N(xnμ)2τ2]+Eτ[(μμ0)2λ0τ2]+C3=Eτ[τ]2{n=1N(xnμ)2+λ0(μμ0)2}+C3

上面的推导中,C, C2, C3指对μ为常数的值。注意Eτ[lnp(τ)]项不是μ的函数,无论μ的值是多少,它都不变。因此,第三行中可以把它吸收到末尾的常数项,第七行也如此。

最后一行是μ的二次多项式。由于这是qμ*(μ)的对数,我们可以看到qμ*(μ)本身是正态分布

通过一些繁琐的计算(展开括号里的平方、分离出涉及μμ2的项、对μ应用配方法),可得到正态分布的参数:

lnqμ*(μ)=Eτ[τ]2{n=1N(xnμ)2+λ0(μμ0)2}+C3=Eτ[τ]2{n=1N(xn22xnμ+μ2)+λ0(μ22μ0μ+μ02)}+C3=Eτ[τ]2{(n=1Nxn2)2(n=1Nxn)μ+(n=1Nμ2)+λ0μ22λ0μ0μ+λ0μ02}+C3=Eτ[τ]2{(λ0+N)μ22(λ0μ0+n=1Nxn)μ+(n=1Nxn2)+λ0μ02}+C3=Eτ[τ]2{(λ0+N)μ22(λ0μ0+n=1Nxn)μ}+C4=Eτ[τ]2{(λ0+N)μ22(λ0μ0+n=1Nxnλ0+N)(λ0+N)μ}+C4=Eτ[τ]2{(λ0+N)(μ22(λ0μ0+n=1Nxnλ0+N)μ)}+C4=Eτ[τ]2{(λ0+N)(μ22(λ0μ0+n=1Nxnλ0+N)μ+(λ0μ0+n=1Nxnλ0+N)2(λ0μ0+n=1Nxnλ0+N)2)}+C4=Eτ[τ]2{(λ0+N)(μ22(λ0μ0+n=1Nxnλ0+N)μ+(λ0μ0+n=1Nxnλ0+N)2)}+C5=Eτ[τ]2{(λ0+N)(μλ0μ0+n=1Nxnλ0+N)2}+C5=12(λ0+N)Eτ[τ](μλ0μ0+n=1Nxnλ0+N)2+C5

注意上面所有步骤用二次和公式都能化简。

qμ*(μ)𝒩(μμN,λN1)μN=λ0μ0+Nx¯λ0+NλN=(λ0+N)Eτ[τ]x¯=1Nn=1Nxn

qτ*(τ)的推导大致相同,简洁起见略去细节。

lnqτ*(τ)=Eμ[lnp(𝐗μ,τ)+lnp(μτ)]+lnp(τ)+constant=(a01)lnτb0τ+12lnτ+N2lnττ2Eμ[n=1N(xnμ)2+λ0(μμ0)2]+constant

两侧取指数,可见qτ*(τ)服从伽马分布

qτ*(τ)Gamma(τaN,bN)aN=a0+N+12bN=b0+12Eμ[n=1N(xnμ)2+λ0(μμ0)2]

计算参数

回忆前几节的结论:

qμ*(μ)𝒩(μμN,λN1)μN=λ0μ0+Nx¯λ0+NλN=(λ0+N)Eτ[τ]x¯=1Nn=1Nxn
qτ*(τ)Gamma(τaN,bN)aN=a0+N+12bN=b0+12Eμ[n=1N(xnμ)2+λ0(μμ0)2]

每种情形下,某变量的分布参数取决于对另一变量的期望。可用正态分布和伽马分布矩期望的标准公式推广期望:

E[τaN,bN]=aNbNE[μμN,λN1]=μNE[X2]=Var(X)+(E[X])2E[μ2μN,λN1]=λN1+μN2

大多数时候将这些公式应用到上述方程不困难,但bN需要更多工作:

bN=b0+12Eμ[n=1N(xnμ)2+λ0(μμ0)2]=b0+12Eμ[(λ0+N)μ22(λ0μ0+n=1Nxn)μ+(n=1Nxn2)+λ0μ02]=b0+12[(λ0+N)Eμ[μ2]2(λ0μ0+n=1Nxn)Eμ[μ]+(n=1Nxn2)+λ0μ02]=b0+12[(λ0+N)(λN1+μN2)2(λ0μ0+n=1Nxn)μN+(n=1Nxn2)+λ0μ02]

这样就可以写出参数方程如下,不需要期望:

μN=λ0μ0+Nx¯λ0+NλN=(λ0+N)aNbNx¯=1Nn=1NxnaN=a0+N+12bN=b0+12[(λ0+N)(λN1+μN2)2(λ0μ0+n=1Nxn)μN+(n=1Nxn2)+λ0μ02]

注意λNbN存在相互依赖关系,自然产生了类似最大期望算法的算法:

  1. 计算n=1Nxn, n=1Nxn2.用所得值计算μN, aN.
  2. 初始化λN为任意值。
  3. λN的现值和其它参数的已知值计算bN
  4. bN的现值和其它参数的已知值计算λN
  5. 重复最后两部,直到收敛(即直到两值的更新量不超过阈值)。

然后就有了后验参数近似分布的超参数值,可用于计算后验的任何属性,如均值、方差、95%最高密度区域(包含总概率95%的最小区间)等等。

可以证明这种算法保证收敛到局部极大值。

还要注意,后验分布与相应的先验分布有同样的形式。只需假设分布可以分解,分布形式便随之而来。事实证明,后验、先验形式相同并非巧合,而是先验分布为指数族成员时的一般结果,大多数标准分布都如此。

进一步讨论

分步法

上面例子展示了给定贝叶斯网络中推导后验概率密度的变分贝叶斯近似的方法:

  1. 图模式描述网络,确定观察变量(数据)𝐗和未观察变量(参数Θ潜变量𝐙)及其条件概率分布。之后,变分贝叶斯构建后验概率p(𝐙,Θ𝐗),这近似具有可分解分布的基本特性,即是多个独立分布在不相交的未观测变量子集上的积。
  2. 将未观测变量划分为多个子集,在其上推导出独立因子。这种方法没有通用的程序,子集太多会使近似结果不佳,子集过少会使整个变分贝叶斯方法变得难以实现。第一种分割方法是将参数和潜变量分开,一般足以产生可操作的结果。划分记作𝐙1,,𝐙M
  3. 对给定的划分𝐙j,用基本方程lnqj*(𝐙j𝐗)=Eij[lnp(𝐙,𝐗)]+constant写出最佳近似分布qj*(𝐙j𝐗)
  4. 用图模式填写联合概率分布公式。任何不涉及𝐙j中变量的分量条件分布都可忽略,它们将被折叠到常数项中。
  5. 按上面例子,简化公式并应用期望算子。理想情况下这应该简化为不属于𝐙j的变量的基本函数的期望(如第一或第二原始、对数期望等等)。为让变分贝叶斯方法顺利运行,这些期望值通常应可以解析地表为变量分布的参数和/或超参数的函数。这些期望项对当前划分的变量都是常数。
  6. 式对当前划分中变量的函数形式表明了分布的类型。特别地,取公式指数,便得到分布的概率密度函数(PDF)(或至少是与之成正比的函数,带未知归一化常数)。要使方法可操作,应能识别出函数形式所属的已知分布;要将公式转为匹配已知分布的PDF的形式,可能需要大量计算。若能做到,就可根据定义恢复归一化常数,并提取公式中的适当部分得到已知分布的参数方程。
  7. 期望都可用非当前划分变量的函数进行解析替换、并将PDF转为可与已知分布相对应的形式,此时结果应是一组方程,将最优参数值表为其他划分变量参数的函数。
  8. 若这程序适用于所有划分,就会产生相互依赖的方程组,指明所有参数的最优值。
  9. 最大期望算法型程序可为每个参数选取初值,再迭代。每一步中,我们都会在方程中循环,依次更新每个参数。可以确证收敛。

另见

参考文献

Template:Reflist

外部链接