提名MCMC(Markov Chain Monte Carlo),用于从一个目标分布 \pi (x) 中抽样。通常抽样的目的是为了计算一个积分,在统计中最为常用的例子就是贝叶斯了。
可以毫不夸张的说,MCMC算法拯救了贝叶斯统计学派。
实现有多简单呢?MCMC中最为通用的一个算法是Metropolis-Hastings算法,非常简单:
其中 q(y|x) 为一个工具分布,通常设置为对称的即 q(y|x)=q(x|y) ,最简单的方法是设定一个随机游走: y=x+\epsilon 。
比如,一个最简单的 随机游走 的M-H算法的函数:
##随机游走的MCMC算法,输入:
## N_samples : 抽样次数
## pai(x) : 目标密度函数
## q_sampler(x): 给定x,从q中抽样的函数
## x0 : 初始值
def
MH_RW
(
N_samples
,
pai
,
q_sampler
,
x0
):
X
=
[]
x
=
x0
for
i
in
range
(
N_samples
):
y
=
q_sampler
(
x
)
rho
=
min
(
1
,
pai
(
y
)
/
pai
(
x
))
if
nprd
.
uniform
()
<=
rho
:
X
.
append
(
y
)
x
=
y
else
:
X
.
append
(
x
)
return
X
就这么简单
理论有多复杂呢?
MCMC构造了一条马尔可夫链,可以证明在适当的条件下,这个马尔可夫链的平稳分布就是\pi (x) ,但是这个马尔可夫链不是最简单的离散状态空间的马尔可夫链,而是一个一般状态空间上的(至少是 \Re^k 的状态空间)一个马尔可夫链。对于离散状态的马尔可夫链,遍历性、常返性、周期性之类的还算是比较容易定义的(虽然本科非数学专业的可能已经有一定接受难度了),对于一般状态空间上的马尔可夫链的遍历性,至少离开测度都很难定义常返性了。
总之,挺难的。