作為一種隨機(jī)采樣方法,馬爾科夫鏈蒙特卡羅(Markov Chain Monte Carlo,以下簡(jiǎn)稱MCMC)在機(jī)器學(xué)習(xí),深度學(xué)習(xí)以及自然語(yǔ)言處理等領(lǐng)域都有廣泛的應(yīng)用,是很多復(fù)雜算法求解的基礎(chǔ)。
從名字我們可以看出,MCMC由兩個(gè)MC組成,即蒙特卡羅方法(Monte Carlo Simulation,簡(jiǎn)稱MC)和馬爾科夫鏈(Markov Chain ,也簡(jiǎn)稱MC)。要弄懂MCMC的原理我們首先得搞清楚蒙特卡羅方法和馬爾科夫鏈的原理。
蒙特卡羅方法
蒙特卡羅原來是一個(gè)賭場(chǎng)的名稱,用它作為名字大概是因?yàn)槊商乜_方法是一種隨機(jī)模擬的方法,這很像賭博場(chǎng)里面的扔骰子的過程。原理是通過大量隨機(jī)樣本,去了解一個(gè)系統(tǒng),進(jìn)而得到所要計(jì)算的值。
它非常強(qiáng)大和靈活,又相當(dāng)簡(jiǎn)單易懂,很容易實(shí)現(xiàn)。對(duì)于許多問題來說,它往往是最簡(jiǎn)單的計(jì)算方法,有時(shí)甚至是唯一可行的方法。
給定統(tǒng)計(jì)樣本集,如何估計(jì)產(chǎn)生這個(gè)樣本集的隨機(jī)變量概率密度函數(shù),是我們比較熟悉的概率密度估計(jì)問題。
求解概率密度估計(jì)問題的常用方法是最大似然估計(jì)、最大后驗(yàn)估計(jì)等。但是,我們思考概率密度估計(jì)問題的逆問題:給定一個(gè)概率分布p(x),如何讓計(jì)算機(jī)生成滿足這個(gè)概率分布的樣本。
這個(gè)問題就是統(tǒng)計(jì)模擬中研究的重要問題–采樣(Sampling)。
如果有一個(gè)我們很難求解出f(x)f(x)的原函數(shù)的函數(shù), 現(xiàn)要求其在定義域 [a, b] 上的積分, 如果這個(gè)函數(shù)是均勻分布, 那么我們可以采樣 [a,b] 區(qū)間的 n 個(gè)值:${x 0,x_1,…x {n-1}}$,用它們的均值來代表 [a,b] 區(qū)間上所有的 f(x) 的值。這樣我們上面的定積分的近似求解為:
如果不是均勻分布, 并 假設(shè)我們可以得到 xx 在[a,b][a,b]的概率分布函數(shù)p(x)p(x) ,那么我們的定積分求和可以這樣進(jìn)行:
上式最右邊的這個(gè)形式就是蒙特卡羅方法的一般形式。當(dāng)然這里是連續(xù)函數(shù)形式的蒙特卡羅方法,但是在離散時(shí)一樣成立。
可以看出,最上面我們假設(shè)x在[a,b]之間是均勻分布的時(shí)候,p(xi)=1/(b−a),帶入我們有概率分布的蒙特卡羅積分的上式,可以得到:
也就是說,我們最上面的均勻分布也可以作為一般概率分布函數(shù)p(x)p(x)在均勻分布時(shí)候的特例。那么我們現(xiàn)在的問題轉(zhuǎn)到了如何在已知分布求出 x 的分布p(x)p(x)對(duì)應(yīng)的若干個(gè)樣本上來。
采樣方法
主要有概率分布采樣及接受-拒絕采樣方法.
概率分布采樣
如果求出了xx的概率分布,我們可以基于概率分布去采樣基于這個(gè)概率分布的 n 個(gè)xx的樣本集,帶入蒙特卡羅求和的式子即可求解。但是還有一個(gè)關(guān)鍵的問題需要解決,即如何基于概率分布去采樣基于這個(gè)概率分布的 n 個(gè)xx的樣本集。
一般而言均勻分布 Uniform(0,1) 的樣本是相對(duì)容易生成的。 通過線性同余發(fā)生器可以生成偽隨機(jī)數(shù),我們用確定性算法生成[0,1]之間的偽隨機(jī)數(shù)序列后,
這些序列的各種統(tǒng)計(jì)指標(biāo)和均勻分布 Uniform(0,1) 的理論計(jì)算結(jié)果非常接近。這樣的偽隨機(jī)序列就有比較好的統(tǒng)計(jì)性質(zhì),可以被當(dāng)成真實(shí)的隨機(jī)數(shù)使用。線性同余隨機(jī)數(shù)生成器如下:
式中a,c,ma,c,m是數(shù)學(xué)推導(dǎo)出的合適的常數(shù)。這種算法產(chǎn)生的下一個(gè)隨機(jī)數(shù)完全依賴當(dāng)前的隨機(jī)數(shù),當(dāng)隨機(jī)數(shù)序列足夠大的時(shí)候,隨機(jī)數(shù)會(huì)出現(xiàn)重復(fù)子序列的情況。
當(dāng)然,也有很多更加先進(jìn)的隨機(jī)數(shù)產(chǎn)生算法出現(xiàn),比如 numpy 用的是 Mersenne Twister 等.
而其他常見的概率分布,無論是離散的分布還是連續(xù)的分布,它們的樣本都可以通過 Uniform(0,1) 的樣本轉(zhuǎn)換而得, 但是如何產(chǎn)生滿足其他分布下的隨機(jī)數(shù)呢?
比如二維正態(tài)分布的樣本(Z1,Z2)(Z1,Z2)可以通過通過獨(dú)立采樣得到的 uniform(0,1) 樣本對(duì)(X1,X2)(X1,X2)通過如下的式子轉(zhuǎn)換而得:
其他一些常見的連續(xù)分布,比如 t分布 , F分布 , Beta分布 , Gamma分布 等,都可以通過類似的方式從 uniform(0,1) 得到的采樣樣本轉(zhuǎn)化得到。在Python的 numpy , scikit-learn 等類庫(kù)中,都有生成這些常用分布樣本的函數(shù)可以使用。
不過很多時(shí)候,我們的x的概率分布不是常見的分布,這意味著我們沒法方便的得到這些非常見的概率分布的樣本集。那這個(gè)問題怎么解決呢?
接受-拒絕采樣
對(duì)于概率分布不是常見的分布,一個(gè)可行的辦法是采用接受-拒絕采樣來得到該分布的樣本。既然 p(x)p(x) 太復(fù)雜在程序中沒法直接采樣,那么我設(shè)定一個(gè)程序可采樣的分布 q(xq(x) 比如高斯分布,然后按照一定的方法拒絕某些樣本,以達(dá)到接近 p(x)p(x) 分布的目的,其中q(x)q(x)叫做 proposal distribution 。
具體采用過程如下,設(shè)定一個(gè)方便采樣的常用概率分布函數(shù) q(x),以及一個(gè)常量 k,使得 p(x) 總在 kq(x) 的下方。
首先,采樣得到 q(x) 的一個(gè)樣本$z 0$,采樣方法如上,使用 uniform(0,1) 轉(zhuǎn)換得到。然后,從均勻分布(0,kq(z0))(0,kq(z0))中采樣得到一個(gè)值uu。如果uu落在了上圖中的灰色區(qū)域,則拒絕這次抽樣,否則接受這個(gè)樣本z0z0。重復(fù)以上過程得到 n 個(gè)接受的樣本 $z_0,z_1,…z {n−1}$,則最后的蒙特卡羅方法求解結(jié)果為:
整個(gè)過程中,我們通過一系列的接受拒絕決策來達(dá)到用q(x)模擬p(x)概率分布的目的。
小結(jié)
使用接受-拒絕采樣,我們可以解決一些概率分布不是常見的分布的時(shí)候,得到其采樣集并用蒙特卡羅方法求和的目的。但是接受-拒絕采樣也只能部分滿足我們的需求,在很多時(shí)候我們還是很難得到我們的概率分布的樣本集。比如:
- 對(duì)于一些二維分布p(x,y)p(x,y),有時(shí)候我們只能得到條件分布p(x|y)p(x|y)和p(y|x)p(y|x),卻很難得到二維分布p(x,y)p(x,y)一般形式,這時(shí)我們無法用接受-拒絕采樣得到其樣本集。
- 對(duì)于一些高維的復(fù)雜非常見分布p(x1,x2,…,xn)p(x1,x2,…,xn),我們要找到一個(gè)合適的q(x)和q(x)和k$非常困難。
從上面可以看出,要想將蒙特卡羅方法作為一個(gè)通用的采樣模擬求和的方法,必須解決如何方便得到各種復(fù)雜概率分布的對(duì)應(yīng)的采樣樣本集的問題。
此時(shí)就需要使用一些更加復(fù)雜的隨機(jī)模擬的方法來生成樣本。比如馬爾科夫鏈蒙特卡羅方法,了解這個(gè)算法我們首先要對(duì)馬爾科夫鏈的平穩(wěn)分布的性質(zhì)有基本的認(rèn)識(shí)。






