网站开发天晟合益福建省建设厅网站职业资格

泊松分布

如果在有限时间$(0,1)$内进行$n$次伯努利实验，那么每次伯努利实验所占用的时间为$\frac{1}{n}$，按照自然规律，一件事情肯定是时间越长越容易发生，假定事件发生的概率与$\frac{1}{n}$成正比，记作$\frac{\lambda }{n}$，则二项分布变为

P { X = k } = ( n k ) ( λ n ) k ( 1 − λ n ) n − k P\{X=k\}=\binom{n}{k}(\frac{\lambda}{n})^k(1-\frac{\lambda}{n})^{n-k} P{X=k}=(kn​)(nλ​)k(1−nλ​)n−k

如果这个时间差$\frac{1}{n}$非常小，以至于事件变得近似连续，则有

lim ⁡ n → ∞ P { X = k } = lim ⁡ n → ∞ ( n k ) ( λ n ) k ( 1 − λ n ) n − k \lim_{n\to\infty} P\{X=k\}=\lim_{n\to\infty}\binom{n}{k}(\frac{\lambda}{n})^k(1-\frac{\lambda}{n})^{n-k} n→∞lim​P{X=k}=n→∞lim​(kn​)(nλ​)k(1−nλ​)n−k

其中

$$\underset{n\to \infty }{\lim }\left(\genfrac{}{}{0ex}{}{n}{k}\right)(\frac{\lambda }{n}{)}^k=\frac{{\lambda }^k}{k!}$$

而后面有一项更是传说中的重要极限

$$\underset{n\to \infty }{\lim }(1-\frac{\lambda }{n}{)}^n=e^{-\lambda }$$

综上就得到了一个新的分布

$$P(x=k)=\frac{e^{-\lambda }{\lambda }^k}{k!}$$

此即泊松分布，表示某个随机事件在连续时间内发生的概率，其中$\frac{\lambda }{n}$表示单位时间内某件事发生的概率。

和二项分布的关系

下面通过scipy中的stats模块，模拟二项分布和泊松分布之间的关联。

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as ssp,q = 0.2, 0.8
ns = [10, 500, 25000, 1250000]fig = plt.figure()
for i,n in enumerate(ns):lam = n * prs = ss.binom(n, p).rvs(50000)rv = ss.poisson(n*p)st, ed = rv.interval(0.999)xs = np.linspace(st, ed, 100)ax = fig.add_subplot(2,2,i+1)ax.hist(rs, density=True, bins='auto', alpha=0.2)ax.plot(xs, rv.pmf(xs))plt.title(f"n={n}")plt.show()

效果如下

和正态分布的关系

泊松分布是二项分布在有限时间内的极限情况，而根据中心极限定理，随着伯努利试验次数的增加，二项分布将逼近于${\sigma }^2=npq,\mu =np$的高斯分布。

在通过二项分布推导泊松分布的过程中，比较关键的一步是$\frac{\lambda }{n}$作为概率的一个假定，当$n\to \infty$时，可以发现$\frac{\lambda }{n}$将趋近于0，如果继续沿用二项分布的模型，那么$q\bar{1}-\frac{\lambda }{n}->0$，从而$npq\approx np$。

由此可以得到，当事件趋近于无穷多时，泊松分布将趋近于$\mu =\lambda ,{\sigma }^2=\lambda$的正态分布。

下面对这个关系进行测试

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as ssp = 0.1
#lam = 100
ns = [10, 500, 25000, 1250000]fig = plt.figure()
for i,n in enumerate(ns):lam = 0.1 * nrs = ss.poisson(lam).rvs(n)rv = ss.norm(lam, np.sqrt(lam))st, ed = rv.interval(0.999)xs = np.linspace(st, ed, 100)ax = fig.add_subplot(2,2,i+1)ax.hist(rs, density=True, bins='auto', alpha=0.2)ax.plot(xs, rv.pdf(xs))plt.title(f"n={n}")plt.show()

效果如下