提名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 的狀態空間)一個馬可夫鏈。對於離散狀態的馬可夫鏈,遍歷性、常返性、周期性之類的還算是比較容易定義的(雖然本科非數學專業的可能已經有一定接受難度了),對於一般狀態空間上的馬可夫鏈的遍歷性,至少離開測度都很難定義常返性了。
總之,挺難的。