python - 如何在 scipy.stats 中指定泊松分布的尾部值？

Question

我想创建一个均值为 2 的泊松分布，元素数量为 10000，最小值为 1，尾值为 140，到目前为止我只能指定最小值

 stats.poisson.rvs( 2, loc = 1,size=10000)

并生成以下分布

如何让它以 140 而不是 11 结束

Answer 1

泊松分布中均值等于 2 的样本大于 140 的可能性很小，以至于您不会从 10000 个样本中得到一个。

事实上，泊松分布只有一个参数 λ 和定义的概率质量函数

P(x=k) = λ^k * exp(-λ) / k!

平均值也等于λ。如果 λ = 2 那么

P(x=140) = 7.7e-199

所以如果有 10000 个样本，则至少有一个样本的可能性 10000 中的 140 将小于 7.7e-195。这是一个如此之小的数字，您不可能期望在一生中发生这种情况。

计算 λ=2 的泊松分布样本位于 140 以上的概率有点困难。您可以使用 scipy.stats.poisson.cdf 来查看

P(x>=22) = 1 - scipy.stats.poisson.cdf(21,2) = 5.5e-16

因此，即使 10000 个样本中有一个样本高于 21 的概率也小于 5.5e-12。以相同的方式计算 P(x>=140) 将返回 0，因为中间结果中有 float 舍入。

结论

如果您希望分布均值等于 2.0，并且重尾在 10000 个样本上达到 140，则您需要不同于泊松的分布。您可以考虑参数 b = 2 的帕累托分布，scipy.stats.pareto。

这里是对来自

的10000个随机样本的比较

scipy.stats.poisson.rvs(2,size=10000)

和

numpy.rint(scipy.stats.pareto.rvs(2,size=10000))

很明显，具有相同均值的 Pareto 分布具有更重的尾部。

下面是该图的代码供引用

import matplotlib.pyplot as plt
import scipy.stats
import numpy as np
pareto_x = np.rint(scipy.stats.pareto.rvs(2,size=10000))
poisson_x = scipy.stats.poisson.rvs(2,size=10000)
plt.figure(figsize=(8,4))
plt.subplot(121)
plt.title("Poisson distribution, a = 2")
plt.xlabel("sample number")
plt.ylabel("sample value")
plt.axis([0,10000,0,180])
plt.plot(range(0,10000),poisson_x,"o")
plt.subplot(122)
plt.axis([0,10000,0,180])
plt.title("Pareto distribution, b = 2")
plt.xlabel("sample number")
plt.plot(range(0,10000),pareto_x,"o")
plt.subplots_adjust(hspace=0.4,bottom=0.2)
plt.savefig("poisson_pareto.png")