Python：从幂律分布生成随机数

Question

我希望从负指数（a = -2）的幂律分布中随机生成介于2到15之间的随机变量。我找到了以下内容：

r = scipy.stats.powerlaw.rvs(a, loc = 2, scale = 13, size = 1000)

但是它不接受负数作为a的值。 有人知道出路吗？

- Panchi

---

尝试使用(1/2)代替。实际上应该是一样的。 - Rob Foley

这里有人回答了：https://dev59.com/z3TYa4cB1Zd3GeqP0fDY，如果你能从中理解些什么的话。 - Faflok

@RobFoley，我需要对不同的负指数进行检查...-1.5、-2、-3等等..因此我希望能得到更一般的建议:(@Faflok，我已经阅读了那篇文章..但是它并没有帮助我摆脱'a>0'的条件。 - Panchi

Answer 1

numpy.random和scipy.stats中定义的幂律分布在数学上不适用于负的a，因为在零点存在奇异性导致它们不可归一化，如这个问题的回答所解释的那样。因此，遗憾的是，数学上说“不行”。

如果您需要在一个不包含零的区间内使pdf与x^{g-1}成比例并且g < 0，你可以定义一个分布。

对于a <= x <= b，有pdf(x) = const * x**(g-1)，其中通过使用均匀变量(np.random.random)进行转换。

In [3]: def rndm(a, b, g, size=1):
    """Power-law gen for pdf(x)\propto x^{g-1} for a<=x<=b"""
   ...:     r = np.random.random(size=size)
   ...:     ag, bg = a**g, b**g
   ...:     return (ag + (bg - ag)*r)**(1./g)

那么你可以做例如以下事情：

In [4]: xx = rndm(1, 2, g=-2, size=10000)

等等等等。

为了完整起见，这里是pdf文件：

In [5]: def pdf(x, a, b, g):
    ag, bg = a**g, b**g
   ....:     return g * x**(g-1) / (bg - ag)

这一切都假定 a < b 且 g != 0。对于 a=0，b=1 和 g > 0，这些公式应与 numpy.power 和 scipy.stats.powerlaw 的结果相符。

- ev-br