非心カイ2乗分布のシミュレーション
<岡本安晴「データを使いこなすための統計入門」22世紀アート>の付録 確率: E.7 非心カイ2乗分布において説明している非心カイ2乗分布のシミュレーションを行ってみた。非心カイ2乗分布は、平均が、標準偏差が1の正規分布の2乗和で表され、このときの非心度はである。シミュレーションのプログラム例をリスト1に示す。
リスト1 非心カイ2乗分布のシミュレーション
import numpy
as np
import scipy.stats
as ss
import matplotlib.pyplot
as plt
n_smpls = 10000
df = int(input('df = '))
lmbd = float(input('noncentrality
= '))
delta = (lmbd
/ df)**0.5
v = ss.norm.rvs(delta,
1, size = (df, n_smpls))
v = v**2
smpls = v.sum(axis
= 0)
plt.hist(smpls,
density = True, label = 'Simulation')
max_x = smpls.max()
x_coord = np.linspace(0,
max_x, 1000)
y_coord = ss.ncx2.pdf(x_coord,
df, lmbd)
plt.plot(x_coord, y_coord, label = 'Model')
plt.title(f'df = {df}')
plt.legend()
plt.show()
リスト1のスクリプトを実行すると、まず自由度と非心度の入力が求められる。
df = 5
noncentrality = 3
非心度を設定すると、計算が始まり、シミュレーションによるサンプルのヒストグラムと非心カイ2乗分布のグラフが表示される(図1)。
図1