產生冪律分佈的隨機數發生器？

Question

我正在爲C++命令行Linux應用程序編寫一些測試。我想用冪律/長尾分佈生成一堆整數。意思是，我經常收到一些數字，但其中大多數並不常見。

理想情況下，我只能用rand（）或stdlib隨機函數之一來使用一些魔術方程。如果沒有，一個簡單易用的C/C++將會非常棒。

謝謝！

Answer 1

這個page at Wolfram MathWorld討論瞭如何從均勻分佈（這是大多數隨機數發生器提供的）得到冪律分佈。

簡短的回答（在上面的鏈接派生）：

x = [(x1^(n+1) - x0^(n+1))*y + x0^(n+1)]^(1/(n+1)) 

其中ý是均勻的變量，Ñ是分佈功率，X0和X1限定的範圍內分佈，而x是你的冪律分佈變量。

Answer 2

如果您知道您想要的分佈（稱爲概率分佈函數（PDF））並將其正確化，則可以將其整合以獲得累積分佈函數（CDF），然後將CDF（如果可能）從統一的[0,1]分配到你想要的轉換。

所以，你首先定義你想要的發行版。

P = F(x) 

（對於x在[0,1]），然後積分得到

C(y) = \int_0^y F(x) dx 

如果能倒你

y = F^{-1}(C) 

於是呼rand()和堵塞結果作爲C在最後一行並使用y。

這個結果被稱爲抽樣的基本定理。這是一個麻煩，因爲規範化要求和需要分析反轉功能。

或者，您可以使用拒絕技術：在所需範圍內均勻地拋出一個數字，然後拋出另一個數字並與第一次拋出的位置處的PDF進行比較。如果第二次投擲超過PDF，則拒絕。對於具有很多低概率區域的PDF，傾向於效率低下，如長尾巴的那些......

中間方法涉及通過強力反轉CDF：將CDF作爲查找表存儲，並執行反轉查找以獲得結果。

---

這裏真正臭氣熏天就是這麼簡單x^-n分佈都在範圍[0,1]非normalizable，所以你不能使用採樣定理。嘗試（x + 1）^ - n改爲...

Answer 3

我無法評論產生冪律分佈所需的數學（其他職位有建議），但我建議您熟悉<random>中的TR1 C++標準庫隨機數設施。這些提供比std::rand和std::srand更多的功能。新系統爲發電機，發動機和配電系統指定了一個模塊化API，並提供了一堆預設。

所包含的分配預設有：

• uniform_int
• bernoulli_distribution
• geometric_distribution
• poisson_distribution
• binomial_distribution
• uniform_real
• exponential_distribution
• normal_distribution
• gamma_distribution

當你定義你的冪律分佈，你應該能夠與現有的發電機和引擎插入。本書Pete Becker的C++標準庫擴展對<random>有很大的幫助。

Here is an article有關如何創建其他分佈（與柯西，卡方，學生噸和費雪˚F例子）

Answer 4

我只是想進行實際模擬作爲補充，（理所當然）接受的答案。儘管在R中，代碼非常簡單，可以成爲（僞） - 僞代碼。在接受的答案和其他的Wolfram MathWorld formula之間

一個微小的差異，也許更常見的，方程是一個事實，即冪指數n（通常稱爲阿爾法）不進行明確的負號。所以選擇的alpha值必須是負數，通常在2和3之間。

x0和x1表示分佈的上限和下限。

所以在這裏，它是：

x1 = 5   # Maximum value 
x0 = 0.1   # It can't be zero; otherwise X^0^(neg) is 1/0. 
alpha = -2.5  # It has to be negative. 
y = runif(1e5) # Number of samples 
x = ((x1^(alpha+1) - x0^(alpha+1))*y + x0^(alpha+1))^(1/(alpha+1)) 
hist(x, prob = T, breaks=40, ylim=c(0,10), xlim=c(0,1.2), border=F, 
col="yellowgreen", main="Power law density") 
lines(density(x), col="chocolate", lwd=1) 
lines(density(x, adjust=2), lty="dotted", col="darkblue", lwd=2) 

或對數刻度繪製：

h = hist(x, prob=T, breaks=40, plot=F) 
    plot(h$count, log="xy", type='l', lwd=1, lend=2, 
    xlab="", ylab="", main="Density in logarithmic scale") 

下面是數據彙總：

> summary(x) 
    Min. 1st Qu. Median Mean 3rd Qu. Max. 
    0.1000 0.1208 0.1584 0.2590 0.2511 4.9388