产生幂律分布的随机数发生器？

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))

其中y是均匀变量，n是分布幂，x0和x1定义了分布的范围，x是幂律分布变量。

Answer 2

如果您知道您想要的分布(称为概率分布函数)并对其进行适当的归一化，则可以对其进行积分以获得累积分布函数，然后对累积分布函数进行倒置(如果可能)，以获得从均匀[0,1]分布到所需的转换。

所以你可以从定义你想要的分布开始。

P = F(x)

(对于0，1中的x)然后进行积分以给出

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

如果这个可以反转，你就会得到

y = F^{-1}(C)

因此，调用rand()并将结果作为C插入到最后一行中，然后使用y。

这个结果被称为抽样基本定理。这是一个麻烦，因为标准化需求和分析反转函数的需要。

或者，您可以使用拒绝技术:在所需的范围内均匀抛出一个数字，然后抛出另一个数字，并在第一次抛出的位置与PDF进行比较。如果第二次投掷超过PDF，则拒绝。对于有很多低概率区域的PDF，比如那些有长尾巴的PDF，效率往往很低。

一种中间方法涉及通过暴力反转CDF :将CDF存储为查找表，然后进行反向查找以获得结果。

这里真正糟糕的是，简单的x^-n分布在范围[0,1]上是不可归一化的，因此您不能使用抽样定理。尝试(x+1)^-n ...

Answer 3

我只是想进行一个实际的模拟，作为(正确的)公认答案的补充。虽然在R中，代码是如此简单，以至于是(伪)伪代码。

公认答案中的Wolfram MathWorld formula与其他可能更常见的方程之间的一个微小差异是，幂律指数 n (通常表示为α)没有显式的负号。因此，选择的alpha值必须为负值，通常介于2和3之间。

x0和x1代表分布的下限和上限。

所以就是这样：

set.seed(0)
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(1e7)   # Number of samples
x  = ((x1^(alpha+1) - x0^(alpha+1))*y + x0^(alpha+1))^(1/(alpha+1))
plot(density(x), ylab="log density x", col=2)

​

​

或以对数刻度绘制：

plot(density(x), log="xy", ylab="log density x", col=2)

​

​

以下是数据的摘要：

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