正態分布可以生成均勻分布嗎?

01-01

這是阿里的一道筆試題. 假如有一個服從正態分布的真隨機數生成器, 是否可以通過程序的處理, 得到一個服從均勻分布的隨機數生成器? 如果可以, 具體應該怎樣做呢? 謝謝!

可以啊。好像e^(x^2)積分不太好積，但是可以用拒絕採樣的方法。

------------------------------------------------------------------------------------------------

好吧，以上做法是錯誤的，想的簡單了。

正確的做法：生成兩個獨立的正太分布變數Z0,Z1,然後arctan(z0/z1)/(2pi)+0.5，可以生成0-1均勻分布的變數,已經通過程序驗證。

好吧。。那就稍微解釋一下 @王芊的做法為什麼對，其實不用計算：熟知生成二維標準正態分布的方法就是取兩個獨立的標準正態分布變數X和Y放在一起(X, Y)就行了，然後二維標準正態分布在直角坐標系裡有各向同性，也就是(X, Y)這個點所指的方向和X軸(或者任何一個給定方向)的夾角是均勻分布的。

怒答一記。這個題目非常有趣，如果你有一定的概率和編程背景其實非常容易。也潛在的考察了面試者對一些基本數值技巧的積累。

一個更常見的問題是，如果你有一個均勻分布的隨機數生成器，如何產生正態（ @王芊你的答案出賣了你，這不是正太分布。。）分布的變數？通常來說我們一般的程序都能更直接地生成均勻分布的隨機變數，而不是想本題所述。

一種業界廣泛使用的做法叫Box-Muller-Wiener演算法。首先生成一對[0,1]上的均勻分布的，獨立的隨機變數A和B。然後用A *2pi 作為角度，sqrt(-2*log(B)) 作為半徑，得到極坐標下的一個點。這個點的兩個坐標（X,Y）就是一個二維的標準正態分布。

那麼回到本題，如何利用正態分布的兩個變數X,Y生成均勻分布的隨機變數？現有的答案們提到了其中一種做法，即計算點（X,Y）在二維平面內與x軸（或任意其他固定方向）的夾角——這樣我們就得到了前述過程中的隨機變數A * 2pi。其實另一種方式是計算exp( -0.5 * (X^2 + Y^2))——這樣我們就得到了前述過程中的隨機變數B。某種意義上來說，後者更為穩定，因為不涉及除法運算（例如arctan(Y/X)），可以以很大概率避免overflow等等問題。

值得注意的是，點(X,Y)到原點距離的平方X^2 + Y^2是服從卡方分布的隨機變數。當變數數為2的時候，這是一個參數為2的指數分布。對於指數分布的隨機變數來說，其累積密度函數有解析的表達式，可以直接通過均勻分布來生成。同樣的辦法對一維正態分布則只能近似數值處理。

若 $X$ 服從正態分布，

則 $Y=Phi(X)$ 服從[0,1]間的均勻分布，其中 $Phi$ 是正態分布函數的cumulative distribution function。

支持一下 @王相及的答案，直接取Y=normcdf(X)。

他的答案評論裡面有個聲稱用matlab做了3年simulation的人說算normcdf需要積分。

關於如何計算erf函數隨便搜了一個：Abramowitz and Stegun. Page 297，顯然可以用級數算。

在matlab里實驗一下兩種方法的用時：

&>&> clear &>&> a = rand(10000); &>&> aa = rand(10000); &>&> tic; b = atan(a); toc; Elapsed time is 1.221234 seconds. &>&> clear b &>&> tic; b = atan(a./aa)/2/pi+0.5; toc; Elapsed time is 1.603967 seconds. &>&> clear b &>&> tic; b = normcdf(a); toc; Elapsed time is 4.962613 seconds. &>&> clear b &>&> tic; b = erfc(a); toc; Elapsed time is 0.565037 seconds. &>&> clear b &>&> tic; b = 0.5 * erfc(-a ./ sqrt(2)); toc; Elapsed time is 1.369548 seconds.

normcdf不是built-in函數，m文件里調用了erfc，所以解釋器會額外花一些時間。如果直接調用erfc計算結果的話，速度基本上是一樣快的，erfc反倒略快一點點。何況算夾角那個需要兩個正態分布才能算。

所以請不要自作聰明地認為計算arctan比計算 $Phi$ 函數要快。

對於任何的隨機變數，如果已知其分布，均可以生產均勻分布。

設隨機變數為 $X$ , 累積分布為 $F_X$ .因為隨機變數的函數也是隨機變數，所以 $F_X(X)$ 也是隨機變數。現在來求其分布， $P{F_X(X) leq a} = P{X leq F^{-1}_X(a)}=F_X(F^{-1}_X(a))=a$ ,考慮到累積函數的值域，其為[0, 1]上的均勻分布。

繼續解釋 @王相及的答案

其實這個是利用了 sklar『s 定理的引理（Inverse Sampling Theorem）只需利用正態分布函數的CDF，以及滿足正態分布的隨機變數X，即可生成一個滿足均勻分布的隨機變數。

方法如下：（就假設最簡單的標準正態分布的情形）

設標準正態分布函數的CDF為F（x），隨機變數X服從標準正態分布

那麼定義隨機變數 Y=F（X）， Y即服從[0,1]上的均勻分布

對於任何連續隨機變數，只要知道其CDF函數貌似都可以由這種方法產生均勻分布。

---------------------------------------------------------------------------------------------------------------------------

下面用matlab代碼給個小例子

X=randn（1000,1）;%產生服從標準正態分布的隨機變數X

Y=normcdf(X); %獲得服從[0,1]上均勻分布的隨機變數Y

---------------------------------------------------------------------------------------------------------------------------

其實這個的證明也不難詳見

https://en.wikipedia.org/wiki/Inverse_transform_sampling

@王芊與 @王相及都已經給出了正確的方法，當然在實際應用上@王芊的方法更好一些，因為只需要計算反三角函數。

下面我來給出另一個通用的方法（實際應用效果不一定好，但是構造具有普適性）：

假設我們手頭有一個分布未知的隨機數生成函數rnd()能夠生成一系列I.I.D. 的(非常數)隨機數，現在需要得到[0, 1]區間上（滿足某個精度）的均勻分布的隨機數函數。

我的方法分兩步(實際上可以推廣到有任一個未知隨機變數得到一個指定的隨機變數)：

第一步，由rnd()函數得到一個滿足 $P(X=0) = P(X=1) = 0.5$ 的隨機函數X
第二步，使用X來生成連續n個數 $x_1,x_2,cdots, x_n$ ，則二進位小數 $0.x_1cdots x_n$ 是滿足[0, 1]區間上均勻分布的一個精度為n的小數。

證明很簡單，下面上代碼：

其中rnd1為一個高斯分布，rnd2為一個伽馬分布，可以看到經過函數uniform變換，都能夠正確生成[0, 1]上的均勻分布。

(代碼在Python3.4下運行通過，Python2.x不一定可以運行)

In [1]: def X(rnd): ...: while True: ...: x, y = rnd(), rnd() ...: if x != y: ...: return 1.0 if x &> y else 0.0


In [2]: def uniform(rnd, precision=8):

   ...:     var = 0.0

   ...:     for n in range(precision):

   ...:         var += X(rnd) / 2**(n+1)

   ...:     return var

In [3]: import random In [4]: rnd1 = lambda : random.gauss(0, 1) In [5]: rnd2 = lambda : random.gammavariate(2, 3) In [6]: xs1 = [uniform(rnd1) for _ in range(100000)] In [7]: xs2 = [uniform(rnd2) for _ in range(100000)] In [8]: import matplotlib.pyplot as plt In [9]: %matplotlib inline In [10]: plt.hist(xs1)