如何編程產生泊松分布的隨機數？

我是一名大二學生，概率論只學了一半，剛學完連續型隨機變數聯合分布函數。我們VBA老師留了一項作業：生成泊松分布的隨機數，手頭工具只能使用Excel自帶的rnd()函數。請問生成的原理是啥？想了好久都想不明白該怎麼做，希望有人能指點一下下，要的是思路和原理，不是代碼，謝謝！

已經百度到的信息：

1. rnd()是均勻分布函數；

2. 有人給的原理貌似是：

x1*x2*...*xn&<=exp^(-lambda)&<=x1*x2*...*x(n+1)

貌似n還是啥就是生成的隨機數，但是這個式子我完全看不懂。

求現階段我能理解的答案。。。非VBA亦可，之前也學過C++，只要是原理就成。

---

你說的這個問題叫做Poisson過程，其產生的思路如下:

1，從的均勻分布出發，用Inverse CDF 方法產生一系列獨立的指數分布（參數為）隨機數；

2，記。如果，則停止，輸出；若否，則繼續生成，直到為止；

3，重複過程2。

容易證明，輸出的一系列整數就滿足服從參數為的Poisson分布。

---

提供一個生成離散型隨機數的一般思路，以泊松分布為例：

非常感謝 @王贇 Maigo的提醒，我把原來的廢話都刪了，思路很簡潔。

只要生成一個0到1之間的隨機數，然後看泊松分布的前幾項和剛好大於這個隨機數就行了。

下面是vba代碼，這個可以擴展到其他的離散型隨機變數，但運行效率不理想：

Function PoissonRand(lambda As Double) As Integer
Dim rand As Single
Dim k As Integer
Randomize
rand = Rnd
k = 0
Do While rand &> WorksheetFunction.POISSON(k, lambda, True)
k = k + 1
Loop
PoissonRand = k
End Function


這裡直接用了POISSON函數可以直接獲得累積概率。

根據 @王贇 Maigo優化版，去除階乘計算提升效率，可讀性稍受影響，且其他離散型隨機變數很難復用：

Function PoissonRand(lambda As Double) As Integer
Dim rand As Single
Dim k As Integer
Dim p As Single
Dim sump As Single

Randomize
rand = Rnd
k = 0
"p(0)化簡
p = 1 / Exp(lambda)
Do While rand &> sump
k = k + 1
"從p(k)轉換到p(k+1)
p = p * lambda / k
sump = sump + p
Loop
PoissonRand = k
End Function


另外，由於浮點數精度、rnd偽隨機性和試驗次數較少等原因，隨機數和理論值存在偏差。

我隨機了10000組PoissonRand(20)：

數據源：

---

Poisson distribution

Wikipedia 上有 Knuth 給出的演算法，其中λ是參數

algorithm poisson random number (Knuth):
init:
Let L ← e?λ, k ← 0 and p ← 1.
do:
k ← k + 1.
Generate uniform random number u in [0,1] and let p ← p × u.
while p &> L.
return k ? 1.


不要問我為什麼， 凡是Knuth說的都對！

－－－－－－－－－－－－－－－－－－－－－－－－－－－

如果將 L都取一個對數，就變成了

algorithm poisson random number (Knuth):
init:
Let L ← ?λ, k ← 0 and p ← 0.
do:
k ← k + 1.
Generate uniform random number u in [0,1] and let p ← p + ln(u).
while p &> L.
return k ? 1.


而 ln(u) 正是產生指數分布。

當然，我們喜歡看正數一些，對 L加個負號，就和頂樓的方法一樣啦。

algorithm poisson random number (Knuth):
init:
Let L ← λ, k ← 0 and p ← 0.
do:
k ← k + 1.
Generate uniform random number u in [0,1] and let p ← p - ln(u).
while p &< L. return k ? 1.

---

Ross, Sheldon M.

Simulation

---

一個比較理想的辦法是構造一個函數Y=f(X)，如果輸入的隨機變數X服從一個分布(例如[0,1)的均勻分布)，那麼輸出的隨機變數Y服從我們需要的分布。但是這個問題里泊松分布是離散的，數學上這個函數比較難寫出來，所以可以用另一個辦法。

考慮泊松分布的物理意義，是一個泊松過程在一個定長時間內的事件發生次數。

假如說你在一個公交車站，一輛公交車到下一輛公交車之間的間隔如果是參數為lambda的指數分布，那一個單位時間內到達公交車的數量服從泊松分布。

餘下你就模擬一下就好了。

---

最快的是不是還是分段rejective sampling

---

Inverse transform sampling

---

看scipy的源碼，不只播送分布，各種分布都有，想看什麼看什麼

---

推薦閱讀：