我們所說的抽樣，其實是指從一個概率分布中生成觀察值（observations）的方法。而這個分布通常是由其概率密度函數（PDF）來表示的。而且， 即使在已知PDF的情況下，讓計算機自動生成觀測值也不是一件容易的事情。從本質上來說，計算機只能實現對均勻分布（Uniform distribution）的采樣。 那如何實現計算機很好的采樣數據樣本呢？今天我們一起來看看實現方法。

在采樣問題上我們可能會面對這些問題：

• 計算機只能實現對均勻分布的采樣，但我們仍然可以在此基礎上對更為復雜的分布進行采樣，那具體該如何操作呢？
• 隨機分布的某些數字特征可能需要通過積分的形式來求解，但是某些積分可能沒有（或者很難求得）解析解，彼時我們該如何處理呢？
• 在貝葉斯推斷中，后驗概率的分布是正?于先驗和似然函數之積的，但是先驗和似然函數的乘積形式可能相對復雜，我們又該如何對這種形式復雜的分布進行采樣呢？

針對這些問題衍生出一系列求解的方法。

一、接受拒絕采樣（Acceptance-Rejection Sampling）

在數學中，拒絕抽樣是用來從分布產生觀測值的基本技術。它也被稱為接受拒絕方法或“接受 - 拒絕算法”，是一種蒙特卡羅方法，這種方法與Metropolis-Hastings算法也有一定關系。

1. 簡單認識

下圖是一個隨機變量的密度函數曲線，試問如何獲得這個隨機變量的樣本呢？

利用這個曲線的形狀來抽取樣本,用一個矩形將這個密度曲線套起來，把密度曲線框在一個矩形里，如下：

然后，向這個矩形里隨機投點。隨機投點意味著在矩形這塊區域內，這些點是滿足均勻分布的。投了大概10000個點,如下面這個樣子：

顯然，有的點落在了密度曲線下側，有的點落在了密度曲線的上側。上側的點用綠色來表示，下側的點用藍色來表示，如下圖：

只保留密度曲線下側的點，即藍色的點：

到這里，提一個問題，在密度曲線以下的這塊區域里，這些點滿足什么分布？均勻分布！這是拒絕采樣最關鍵的部分，搞個矩形、向矩形里投點等等，所做的一切都是為了獲得一個密度曲線所圍成區域的均勻分布。只要能獲得這樣一個在密度曲線下滿足均勻分布的樣本，我們就可以獲得與該密度曲線相匹配的隨機變量的采樣樣本。方法是，只需把每個藍點的橫坐標提取出來，這些橫坐標所構成的樣本就是我們的目標樣本。下圖左側，是按照以上方法獲得的一個樣本的直方圖以及核密度估計，下圖右側，是開始的密度曲線。

可見，采樣樣本的核密度估計與目標密度曲線基本一致，可以肯定這個樣本就是目標樣本。

最開始時候用到了一個矩形，這個矩形就是一個滿足均勻分布的建議分布，建議分布只是獲得目標密度函數曲線下均勻分布樣本的一個輔助工具。采用均勻分布作為建議分布有時效率很低，為什么這么說？從上例就可以看出，均勻分布的好多點（那些綠點）都被剔除了，造成了一種浪費。可以選擇一些其他曲線來把密度曲線框起來，效率會提高一點，如下圖：

數曲線為h(x), 對應于下圖中的藍線，建議分布密度曲線為g(x),我們把g(x)乘上一個常數因子c,然后用cg(x)這條曲線將目標密度曲線框起來。

我們假定滿足g(x)的隨機變量易于采樣，那么拒絕采樣的步驟如下：

• 從g(x)采到一個樣本數據，記為$x^{?}$，我們把它作為一個建議

• 要不要接受這個建議，作為滿足h(x)分布的一個樣本數據呢？我們定義一個接受概率:

也就是說，我們以$\alpha$的概率接受$x^{?}$作為h(x)分布的一個樣本數據。實際操作中，我們是取一個$U(0,1)$的隨機數$\mu$,如果$\mu <\alpha$,就接受$x^{?}$作為h(x)的一個樣本數據，否則，把它舍棄掉，回到1步繼續循環。最終可以獲得一個樣本。

• 文章開頭是一下子抽取10000個點，到后來怎么成了一個個抽了呢？其實它們是對應的，把藍點去掉的過程就相當于你做是否拒絕判斷的過程。
• 如果有密度曲線下的均勻分布樣本，就可以得到與密度曲線相匹配的分布的一個樣本。
• 如果建議分布的形狀和目標分布越接近，采樣的效率就越高。

2. Acceptance-Rejection Sampling過程

3. Acceptance-Rejection Sampling的直觀解釋

4. Acceptance-Rejection Sampling有效性證明(待)

5.python實現

2. 生成代碼如下：

import random import math import matplotlib.pyplot as plt import seaborn as sns import numpy as np%matplotlib inline sns.set_style('darkgrid') plt.rcParams['figure.figsize'] = (12, 8)def AceeptReject(split_val):global cglobal powerwhile True:x = random.uniform(0, 1)y = random.uniform(0, 1)if y*c <= math.pow(x - split_val, power):return xpower = 4 t = 0.4 sum_ = (math.pow(1-t, power + 1) - math.pow(-t, power + 1)) / (power + 1) #求積分 x = np.linspace(0, 1, 100) #常數值c c = 0.6**4/sum_ cc = [c for xi in x] plt.plot(x, cc, '--',label='c*f(x)') #目標概率密度函數的值f(x) y = [math.pow(xi - t, power)/sum_ for xi in x] plt.plot(x, y,label='f(x)') #采樣10000個點 samples = [] for i in range(10000):samples.append(AceeptReject(t)) plt.hist(samples, bins=50, normed=True,label='sampling') plt.legend() plt.show()

5. 小結

要想將蒙特卡羅方法作為一個通用的采樣模擬求和的方法，還的需馬爾科夫鏈的幫忙。

https://gaolei786.github.io/statistics/reject.html
https://zhuanlan.zhihu.com/p/75264565

總結

以上是生活随笔為你收集整理的接受拒绝采样（Acceptance-Rejection Sampling）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。