目录
Metropolis-Hasting抽样算法
1. 随机模拟的基本思想
假设我们有一个矩形区域 R R R,面积为 S 0 S_0 S0。在此区域中,有一个不规则区域 M M M,其面积 S S S待求。
方法1:把不规则区域 M M M划分为多个小的规则区域,由这些规则区域的面积总和 S ′ S' S′近似。
方法2:抓一把黄豆,均匀铺在 R R R中,再统计落在不规则区域 M M M中的黄豆比例。该比例可近似 S S 0 \frac {S} {S_0} S0S。
方法2采用的是抽样(采样)解决问题的思想,即随机模拟。
因此,随机模拟的关键,在于如何将实际复杂问题,转化为抽样可以解决的问题。当然,这非常灵活,也考验我们的创造力。
我们接下来的核心问题,是利用计算机可直接实现的均匀抽样,借助随机模拟的方法,实现满足特定概率分布 p ( x ) p(x) p(x)的抽样。
2. 拒绝抽样
当我们的概率分布很简单时,如 [ 0 , 1 ] [0,1] [0,1]上的均匀分布,抽样是很好实现的。现成的随机算法非常多。
但是,对于其他更复杂的分布,如果我们还想借助简单的均匀随机抽样,就需要换换思路了。
例如上图中的 p ( z ) p(z) p(z),是较为复杂的概率密度函数。
我们先构造比较简单的、容易抽样的高斯分布 q ( z ) q(z) q(z),乘一个常数 k k k,使之满足:
k q ( z ) ⩾ p ( z ) , ∀ z kq(z) \geqslant p(z), \forall z kq(z)⩾p(z),∀z
拒绝抽样的方法是:
- 对高斯分布的样本进行采样,得到一个 z 0 z_0 z0,如图;
- 在 [ 0 , k q ( z 0 ) ] [0,kq(z_0)] [0,kq(z0)]上均匀采样,得到 u 0 u_0 u0;
- 如果 u 0 > p ( z 0 ) u_0>p(z_0) u0>p(z0),就拒绝该采样结果;反之接受;
- 重复直到接受足够多的样本。
我们理解一下为什么拒绝抽样是可行的,这对我们之后运用随机模拟思想非常关键!
- 该简单分布与高斯分布的趋势大致相同,中间大两头小,因此取到中间的 z z z的概率比较大;
- 通过 u 0 u_0 u0进一步筛选 z z z, p ( z ) p(z) p(z)较低时拒绝率较大,反之接受率较高,大致符合 p ( z ) p(z) p(z)分布要求。
一句话,通过简单抽样和设置接受率,“强迫”结果趋于复杂分布。
这样,我们就通过简单的高斯分布采样,以及计算机可直接实现的均匀采样,间接实现了复杂的 p ( z ) p(z) p(z)采样。
然而,在高维情况下,该方法的劣势很明显:
- 合适的简单分布 q ( z ) q(z) q(z)很难找到。
- 合适的 k k k值很难找到。
这两个问题会导致拒绝率极高,计算效率很低。其根本缺陷在于:样本之间是独立无关的。
3. Metropolis-Hastings抽样
3.1 引入思想
结合马氏链知识(参考信息论相关书籍),我们知道:
对于遍历的马尔可夫链,当其转移次数足够多时,将会进入稳态分布 π ( x ) \pi (x) π(x),即各状态的出现概率趋于稳定。
进一步,假设变量 X X X所有可能的取值,构成某遍历马氏链的状态空间。
则无论从什么状态出发,只要转移步数足够大,该马氏链将趋于稳定,即逐渐开始依次输出符合 p ( x ) p(x) p(x)分布的状态 X i X_i Xi。
这样,通过收集马氏链收敛后产生的 X i X_i Xi,我们就得到了符合分布 p ( x ) p(x) p(x)的样本。
此时稳态分布即