如何在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个样品的机会，将有在140是至少一个样品 出的10000将小于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，因为中间结果中的浮点四舍五入。

结论

如果你想分布均值等于2.0，和重尾达到高达140上10000样本，需要从不同的泊松分布。你可以考虑帕累托分布，scipy.stats.pareto与参数b = 2

这里是10000个随机样本

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

和

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

的比较可以清楚地看到具有相同平均值的帕累托分布具有更重的尾部。

供参考的情节的代码如下

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