機率與數理統計-蒙特卡洛之拒絕採樣

聲明

轉：> https://blog.csdn.net/u010159842/article/details/78959515

介紹

蒙特卡洛（Monte Carlo）方法是二十世紀四十年代中期因爲科學技術的發展和電子計算機的發明，而被提出的一種以機率統計理論爲基礎的數值計算方法。它的核心思想就是使用隨機數（或更常見的僞隨機數）來解決一些複雜的計算問題。

當所求解問題能夠轉化爲某種隨機分佈的特徵數（好比隨機事件出現的機率，或者隨機變量的指望值等）時，每每就能夠考慮使用蒙特卡洛方法。經過隨機抽樣的方法，以隨機事件出現的頻率估計其機率，或者以抽樣的數字特徵估算隨機變量的數字特徵，並將其做爲問題的解。這種方法多用於求解複雜的高維積分問題。

實際應用中，咱們所要面對的第一個問題就是如何抽樣？在統計學中， 抽樣（或稱採樣）是指從目標整體中抽取一部分個體做爲樣本的過程。

例如，咱們想知道一所大學裏全部男生的平均身高。可是由於學校裏的男生可能有上萬人之多，因此爲每一個人都測量一下身高可能存在困難，因而咱們從每一個學院隨機挑選出100名男生來做爲樣本，這個過程就是抽樣

可是在計算機模擬時，咱們所說的抽樣，實際上是指從一個機率分佈中生成觀察值（observations）的方法。而這個分佈一般是由其機率密度函數（PDF）來表示的。並且，即便在已知PDF的狀況下，讓計算機自動生成觀測值也不是一件容易的事情。從本質上來講，計算機只能實現對均勻分佈（Uniform distribution）的採樣。

具體來講，咱們可能要面對的問題包括：

• 計算機只能實現對均勻分佈的採樣，但咱們仍然能夠在此基礎上對更爲複雜的分佈進行採樣，那具體該如何操做呢？
• 隨機分佈的某些數字特徵可能須要經過積分的形式來求解，可是某些積分可能沒有（或者很難求得）解析解，彼時咱們該如何處理呢？
• 在貝葉斯推斷中，後驗機率的分佈是正⽐於先驗和似然函數之積的，可是先驗和似然函數的乘積形式可能相對複雜，咱們又該如何對這種形式複雜的分佈進行採樣呢？

Inverse CDF 方法

比較簡單的一種狀況是，咱們能夠經過PDF與CDF之間的關係，求出相應的CDF。或者咱們根本就不知道PDF，可是知道CDF。此時就可使用Inverse CDF的方法來進行採樣。這種方法又稱爲逆變換採樣（Inverse transform sampling）。

若是你對PDF和CDF的概念有點模糊，咱們不妨先來一塊兒回顧一下它們的定義。對於隨機變量$X$，以下定義的函數 $F$: $$F(x)=P{X<= x}, -\infty < x < \infty$$ 稱爲 $X$ 的累積分佈函數（CDF，Cumulative Distribution Function）。對於連續型隨機變量 $X$ 的累積分佈函數 $F(x)$，若是存在一個定義在實數軸上的非負函數 $f(x)$，使得對於任意實數 $x$，有下式成立： $$F(x)=\int_{-\infty}^{\infty}{f(t)}dt$$

則稱 $f(x)$ 爲 $X$ 的機率密度函數（PDF，Probability Density Function）。顯然，當機率密度函數存在的時候，累積分佈函數是機率密度函數的積分。

因此，一般咱們能夠經過對PDF（以下圖中的左圖所示爲正態分佈的PDF）進行積分來獲得機率分佈的CDF（以下圖中的右圖所示爲正態分佈的CDF）。而後咱們再獲得CDF的反函數 $F^{-1}(u)$，若是你想獲得 $m$ 個觀察值，則重複下面的步驟 $m$ 次:

• 從 Uniform(0,1) 中隨機生成一個值（前面已經說過，計算機能夠實現從均勻分佈中採樣），用 $u$ 表示。
• 計算 $F_{-1}(u)$ 的值 $x$，則 $x$ 就是從 $f(x)$ 中得出的一個採樣點。

如下圖爲例，若是從 Uniform(0,1) 中隨機生成的值 $u=0.8413$，則能夠算得$F_{-1}(u)=1$，則這次從正態分佈中生成的隨機數就是 1。

你可能會好奇，面對一個具備複雜表達式的函數， Inverse CDF 方法真的有效嗎？來看下面這個例子。假設如今咱們但願從下面這個PDF中抽樣：

$$ f(x)=\left{\begin{array}{ll}{8 x} & {, \text {if } 0 \leq x<0.25} \ {\frac{8}{3}-\frac{8}{3} x} & {, \text { if } 0.25 \leq x \leq 1} \ {0} & {, \text { otherwise }}\end{array}\right. $$

能夠算得相應的CDF爲

$$ F(x)=\left{\begin{array}{ll}{0} & {, \text { if } x<0} \ {4 x^{2}} & {, \text { if } 0 \leq x<0.25} \ {\frac{8}{3} x-\frac{4}{3} x^{2}-\frac{1}{3}} & {, \text { if } 0.25 \leq x \leq 1} \ {1} & {, \text { if } x>1}\end{array}\right. $$

對於 $u∈[0,1]$，它的反函數爲:

$$ F^{-1}(u)=\left{\begin{array}{ll}{\frac{\sqrt{u}}{2}} & {, \text { if } 0 \leq u<0.25} \ {1-\frac{\sqrt{3(1-u)}}{2}} & {, i f 0.25 \leq u \leq 1}\end{array}\right. $$

拒絕採樣

咱們已經看到 Inverse CDF 方法確實有效。但其實它的缺點也是很明顯的，那就是有些分佈的 CDF 可能很難經過對 PDF 的積分獲得，再或者 CDF 的反函數也很不容易求。這時咱們可能須要用到另一種採樣方法，這就是咱們即將要介紹的拒絕採樣。

下面這張圖很好地闡釋了拒絕採樣的基本思想。假設咱們想對 PDF 爲 $p(x)$ 的函數進行採樣，可是因爲種種緣由（例如這個函數很複雜），對其進行採樣是相對困難的。可是另一個 PDF 爲 $q(x)$ 的函數則相對容易採樣，例如採用 Inverse CDF 方法能夠很容易對對它進行採樣，甚至 $q(x)$ 就是一個均勻分佈（別忘了計算機能夠直接進行採樣的分佈就只有均勻分佈）。那麼，當咱們將 $q(x)$ 與一個常數 MM 相乘以後，能夠實現下圖所示之關係，即$M⋅q(x)$ 將 $p(x)$ 徹底「罩住」

你固然能夠採用嚴密的數學推導來證實Reject Sampling的可行性。但它的原理從直觀上來解釋也是至關容易理解的。你能夠想象一下在上圖的例子中，從哪些位置抽出的點會比較容易被接受。顯然，紅色曲線和綠色曲線所示之函數更加接近的地方接受機率較高，也便是更容易被接受，因此在這樣的地方採到的點就會比較多，而在接受機率較低（即兩個函數差距較大）的地方採到的點就會比較少，這也就保證了這個方法的有效性。

Reject Sampling方法確實能夠解決咱們的問題。可是它的一個不足涉及到其採樣效率的問題。 最理想的狀況下，參考分佈應該跟目標分佈越接近越好，從圖形上來看就是包裹的越緊實越好。可是這種狀況的參考分佈每每又不那麼容易獲得。在知足某些條件的時候也確實能夠採用所謂的改進方法，即Adaptive Rejection Sampling。

自適應的拒絕採樣（Adaptive Rejection Sampling）

前面咱們已經分析了，拒絕採樣的弱點在於當被拒絕的點不少時，採樣的效率會很是不理想。同時咱們也支持，若是可以找到一個跟目標分佈函數很是接近的參考函數，那麼就能夠保證被接受的點佔大多數（被拒絕的點不多）。這樣一來便克服了拒絕採樣效率不高的弱點。若是函數是 log-concave 的話，那麼咱們就能夠採樣自適應的拒絕採樣方法。什麼是 log-concave 呢？仍是回到咱們以前介紹過的 Beta 分佈的PDF，咱們用下面的代碼來繪製 Beta(2, 3) 的函數圖像，以及將 Beta(2, 3) 的函數取對數以後的圖形

(......更詳細的內容看原文)