In [55]:

    

# ポアソン分布に従う乱数を累積分布関数から生成
import numpy as np
from scipy.stats import uniform, poisson
import matplotlib.pyplot as plt

# ポアソン分布のパラメータ
lam = 4

# この値まで確率値を計算
K = 20

# サンプリング数
N = 10000

rv = poisson(mu=lam)

# cdfの表を計算
t = np.arange(K)
prob = rv.cdf(t)

X = []
for i in range(N):
    u = uniform.rvs(loc=0, scale=1, size=1)
    # prob < uはcdfがuより小さいときTRUEを返す
    # TRUEは1と解釈されるのでsum()でTRUEの数をカウントして
    # インデックスを求めている
    X.append(np.sum(prob < u))

# ポアソン分布に従う乱数の分布を描画
# hist()のnormed=Trueはバーの積分が1になる確率密度関数になるため離散分布では使えない
# 離散分布ではバーの高さの合計が1になる確率質量関数にする必要がある
# http://stackoverflow.com/questions/3866520/plotting-histograms-whose-bar-heights-sum-to-1-in-matplotlib
plt.figure(1)
nbins = np.arange(-0.5, K, 1.0)
weights = np.ones_like(X) / float(len(X))
plt.hist(X, nbins, weights=weights)
plt.plot(t, rv.pmf(t), 'ro-', lw=1)
plt.xlim((0, K))

plt.show()


    

In [56]:

    

# 二項分布に従う乱数を累積分布関数から生成
import numpy as np
from scipy.stats import uniform, binom
import matplotlib.pyplot as plt

# 二項分布のパラメータ
n = 20
p = 0.5

# この値まで確率値を計算
K = 40

# サンプリング数
N = 10000

rv = binom(n=n, p=p)

# cdfの表を計算
t = np.arange(K)
prob = rv.cdf(t)

X = []
for i in range(N):
    u = uniform.rvs(loc=0, scale=1, size=1)
    # prob < uはcdfがuより小さいときTRUEを返す
    # TRUEは1と解釈されるのでsum()でTRUEの数をカウントして
    # インデックスを求めている
    X.append(np.sum(prob < u))

# 二項分布に従う乱数の分布を描画
# hist()のnormed=Trueはバーの積分が1になる確率密度関数になるため離散分布では使えない
# 離散分布ではバーの高さの合計が1になる確率質量関数にする必要がある
# http://stackoverflow.com/questions/3866520/plotting-histograms-whose-bar-heights-sum-to-1-in-matplotlib
plt.figure(1)
nbins = np.arange(-0.5, K, 1.0)
weights = np.ones_like(X) / float(len(X))
plt.hist(X, nbins, weights=weights)
plt.plot(t, rv.pmf(t), 'ro-', lw=1)
ｐｌ
plt.show()