拉普拉斯近似

拉普拉斯近似(Laplace Approximation)

对于后验分布,往往是难以解析表达的概率分布形式,进而其性质难以描述,例如期望的计算需要积分,而后验分布常常是难以积分的。

对此,一个简单的想法是利用某种性质良好的概率分布去近似后验分布。利用高斯分布近似后验分布的方法就是拉普拉斯近似。

公式

后验分布来源于贝叶斯公式。
p ( θ ∣ y ) = p ( y ∣ θ ) p ( θ ) p ( y ) p(\theta|y)=\frac{p(y|\theta)p(\theta)}{p(y)} p(θy)=p(y)p(yθ)p(θ)
在此公式中,似然函数 p ( y ∣ θ ) p(y|\theta) p(yθ)通常是已知的,先验分布 p ( θ ) p(\theta) p(θ)也是已知的,而分母 p ( y ) = ∫ θ p ( y ∣ θ ) p ( θ ) d θ p(y)=\int_{\theta}p(y|\theta)p(\theta)d\theta p(y)=θp(yθ)p(θ)dθ通常是难以计算的。其原因在于需要在 θ \theta θ上进行积分,而两个密度函数的乘积会非常难以积分,特别是 θ \theta θ常常是多维的。
另一方面,公式左侧是后验分布,是参数 θ \theta θ的函数,与变量 y y y无关。因此,不妨将公式右侧的分子部分记作:
p ^ ( z ) = p ( y ∣ θ ) p ( θ ) \hat{p}(z)=p(y|\theta)p(\theta) p^(z)=p(yθ)p(θ)
省略 y y y,进而将归一化的后验分布记作:
p ( z ) = p ( z ) ^ Z p(z)=\frac{\hat{p(z)}}{Z} p(z)=Zp(z)^其中, Z Z Z是归一化因子。注意这里之所以将 y y y省略是因为 y y y是某一固定取值,而非变量。这里将后验概率写为 p ( z ) p(z) p(z)是为了和先验分布相区分。
在这样的设定下,拉普拉斯近似首先需要通过某种优化或其他方法,进行最大后验估计,从而得到MAP估计值。
z M A P = arg max ⁡ z ∈ Θ p ( z ) z_{MAP}=\argmax_{z \in {\Theta}}p(z) zMAP=zΘargmaxp(z)
然后,我们需要在MAP值处对后验概率进行二阶的泰勒展开。进行二阶泰勒展开的原因是高斯分布是二次型,因此二阶泰勒展开会自然而然的导出一个高斯分布。
l n p ( z ) ^ ≈ l n p ( z ) ^ − 1 2 ( z − z M A P ) T H ( z − z M A P ) ln\hat{p(z)}\approx ln\hat{p(z)}-\frac{1}{2}(z-z_{MAP})^TH(z-z_{MAP}) lnp(z)^lnp(z)^21(zzMAP)TH(zzMAP)其中, H H H是Hessian矩阵。那么一次项没有的原因是因为 z M A P z_{MAP} zMAP是最大后验处,雅可比是零蛋。
两边同时取指数,仍然成立,有:
p ( z ) ^ = p ( z M A P ) ^ e x p ( − 1 2 ( z − z M A P ) T H ( z − z M A P ) ) \hat{p(z)}=\hat{p(z_{MAP)}}exp(-\frac{1}{2}(z-z_{MAP})^TH(z-z_{MAP})) p(z)^=p(zMAP)^exp(21(zzMAP)TH(zzMAP))
这是对未归一化概率分布的一个二次项,则高斯分布前面的常数就很容易的通过凑高斯的方法获得,协方差矩阵的逆由Hessian矩阵给出,维度已知,则常数可以写出。
p ( z ) = 1 ( 2 π ) k 2 ∣ H ∣ − 1 2 e x p ( − 1 2 ( z − z M A P ) T H ( z − z M A P ) ) p(z)=\frac{1}{(2\pi)^{\frac{k}{2}}|H|^{-\frac{1}{2}}}exp(-\frac{1}{2}(z-z_{MAP})^TH(z-z_{MAP})) p(z)=(2π)2kH211exp(21(zzMAP)TH(zzMAP))
上式就是完整的后验分布的拉普拉斯近似。

证据项的估计

基于对后验分布的拉普拉斯近似,对于贝叶斯公式上难以求得的分母也有了表达:
∫ θ p ( θ ∣ y ) p ( θ ) d θ ≈ ∫ z p ( z M A P ) ^ e x p ( − 1 2 ( z − z M A P ) T H ( z − z M A P ) ) d z \int_{\theta}p(\theta|y)p(\theta)d\theta\\\approx \int_{z}\hat{p(z_{MAP)}}exp(-\frac{1}{2}(z-z_{MAP})^TH(z-z_{MAP}))dz θp(θy)p(θ)dθzp(zMAP)^exp(21(zzMAP)TH(zzMAP))dz
注意, θ \theta θ z z z是一个变量,只不过上面为了区分后验和先验分布的表达,写作了两个变量。
把前面的 p ( z M A P ^ ) \hat{p(z_{MAP}}) p(zMAP^)提出来,然后积里面的高斯分布,实际上就会得到高斯分布前面归一化常数的倒数,也就是 ( 2 π ) k 2 ∣ H ∣ − 1 2 (2\pi)^{\frac{k}{2}}|H|^{-\frac{1}{2}} (2π)2kH21这个东西,因此,将积分结果再乘以提出去的 p ( z M A P ^ ) \hat{p(z_{MAP}}) p(zMAP^),就得到: ( 2 π ) k 2 ∣ H ∣ − 1 2 p ( z M A P ^ ) (2\pi)^{\frac{k}{2}}|H|^{-\frac{1}{2}}\hat{p(z_{MAP}}) (2π)2kH21p(zMAP^)
p ( z M A P ^ ) \hat{p(z_{MAP}}) p(zMAP^)完整写出,在整理,就是:
p ( y ) = ( 2 π ) k 2 ∣ H ∣ − 1 2 p ( y ∣ θ M A P ) p ( θ M A P ) p(y)=(2\pi)^{\frac{k}{2}}|H|^{-\frac{1}{2}}p(y|\theta_{MAP})p(\theta_{MAP}) p(y)=(2π)2kH21p(yθMAP)p(θMAP)
这在贝叶斯实验设计中可用作采样方法的一个近似。(Ryan, 2003)

Kenneth J Ryan (2003) Estimating Expected Information Gains for Experimental Designs With Application to the Random Fatigue-Limit Model, Journal of Computational and Graphical Statistics, 12:3, 585-603.

评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值