ポアソン分布について、python のライブラリでの計算方法を調べてみました。
結果を記載します。
ポアソン分布について
そもそもポアソン分布とは何かの理解が曖昧なので、書籍、Webで調べた結果をざっくり記載します。
- 滅多に起こり得ない希少な事象の発生数の確率が取り得る分布
-
ポアソン分布では，平均と分散が等しくなる
-
ポアソン分布となる現実での例 (統計学入門 p115 からの引用)
- 交通事故件数
- 大量生産の不良品数
- 破産件数
- 火災件数
- 砲弾命中数
- 遺伝子の突然変異数
-
二項分布と、ポアソン分布の関係
p = 0.002の場合等、二項分布の成功確率pが十分に小さい場合、二項分布は、ポアソン分布に近似できる。
ポアソン分布の計算ができるライブラリについて
以下、ライブラリでポアソン分布の計算ができるようです。
statsmodels.families.Poisson というのもありましたが、こちらはポワソン回帰を計算する際に利用するライブラリのようで、分布の計算ができるわけではなさそうでしたので、ここでは記載しません。
- numpy.random.poisson
- scipy.stats.poisson
各ライブラリを使ってみる
numpy.random.poisson
numpy.random.poisson は、ポアソン分布に従う、ndarray を生成することができます。
import numpy as np
# λ=518/365 試行回数は100000 * 2回
s = np.random.poisson(518./365., [10000,2])
# ヒストグラム描画
%matplotlib inline
import matplotlib.pyplot as plt
count, bins, ignored = plt.hist(s, 14, normed=True)
plt.show()
# -------------------
# 発生確率を2倍にする
# -------------------
import numpy as np
# λ=518 * 2/365 試行回数は10000 * 2 回
s = np.random.poisson(518. * 2./365., [10000,2])
# ヒストグラム描画
%matplotlib inline
import matplotlib.pyplot as plt
count, bins, ignored = plt.hist(s, 14, normed=True)
plt.show()
scipy.stats.poisson
scipy.stats.poisson を使い、確率質量関数 の描画、ある事象が起きる確率を計算することができます。
- 確率質量関数 を描画する
from scipy.stats import poisson
import matplotlib.pyplot as plt
fig, ax = plt.subplots(1, 1)
# λ を指定する
mu = 518 / 365
# mements は 文字列 mvsk を指定すると、 m > 平均 v > 分散 > s > 歪度（わいど) k > 尖度(せんど) が返る
mean, var, skew, kurt = poisson.stats(mu, moments='mvsk')
print("平均", mean)
print("分散", var)
print("歪度", skew)
print("尖度", kurt)
# グラフ描画
x = np.arange(poisson.ppf(0.01, mu),
poisson.ppf(0.99, mu))
ax.plot(x, poisson.pmf(x, mu), 'bo', ms=8, label='poisson pmf')
ax.vlines(x, 0, poisson.pmf(x, mu), colors='b', lw=5, alpha=0.5)
rv = poisson(mu)
ax.vlines(x, 0, rv.pmf(x), colors='k', linestyles='-', lw=1,
label='frozen pmf')
ax.legend(loc='best', frameon=False)
plt.show()
平均 1.4191780821917808
分散 1.4191780821917808
歪度 0.8394243293074157
尖度 0.7046332046332047
-
ヒストグラムの見方について
個人的な理解では以下の通りです。- ある期間内に、事象が、0 回起きる確率(発生しない確率)は、0.24
- ある期間内に、事象が、1 回起きる確率は、0.35
- ある期間内に、事象が、2 回起きる確率は、0.24
- ある期間内に、事象が、3 回起きる確率は、0.12
- ある期間内に、事象が、4 回起きる確率は、0.04
-
確率の計算をする
# scipy.statsから、 poissonをインポートします
from scipy.stats import poisson
# λ を指定する
mu = 518 / 365
# 事象が発生しない確率を求める
zero = poisson.pmf(0,mu)
print('事象が発生しない確率:', zero)
# 事象が1回発生する確率を求める
one = poisson.pmf(1,mu)
print('事象が発生する確率:', one)
# 事象が5回発生する確率を求める
five = poisson.pmf(5,mu)
print('事象が発生する確率:', five)
事象が発生しない確率: 0.241912767619
事象が発生する確率: 0.343317297608
事象が発生する確率: 0.011605450952
参考
以下、参考になりました。
以上です。
コメント