目次
Python で統計学
カイ2乗分布
推定や仮説検定で重要になるカイ2乗分布の密度関数は， \begin{eqnarray*} f_n(x) = \frac{1}{2^{n/2}\Gamma\left(\frac{n}{2}\right)} x^{n/2-1}e^{-x/2}, \quad x > 0 \end{eqnarray*} である．ここで \(n~(=1,2,\cdots)\) は自由度である．つまりデータサイズが \(n+1\) のときに，自由度が \(n\) になる．さらに，\(\Gamma(\alpha)\) はガンマ関数であり， \begin{eqnarray*} \Gamma(\alpha) = \int_{0}^{\infty} x^{\alpha - 1}e^{-x}dx, \quad \alpha > 0 \end{eqnarray*} で定義される．カイ2乗分布の平均は \(n\) であり，分散は \(2n\) である．
カイ2乗分布に関する重要な定理として次のようなものが知られている．
- 定理
- 標準正規分布 \(\mathcal{N}\left(0, 1\right)\) に従う確率変数列 \(Z_1, Z_2, \cdots, Z_n\) が独立ならば，\(X_n = Z_1^2 + \cdots + Z_n^2\) は自由度 \(n\) のカイ2乗分布に従う．
ライブラリのインポート
Python でカイ2乗分布を利用するためには Scipy ライブラリをインポートするとよい．
import numpy as np
from scipy.stats import chi2
import matplotlib.pyplot as plt
from IPython.display import set_matplotlib_formats
# from matplotlib_inline.backend_inline import set_matplotlib_formats # バージョンによってはこちらを有効に
set_matplotlib_formats('retina')
確率密度関数
確率密度関数の値を取得するには chi2.pdf(x, df) を使用すると良い．例えば，自由度 5 のカイ2乗分布について，いくつか密度関数の値を計算してみよう．
chi2.pdf(x = 0, df = 5)
0.0
chi2.pdf(x = 1, df = 5)
0.0806569081730478
chi2.pdf(x = 2, df = 5)
0.1383691658068649
自由度 5 の確率密度関数をプロットしてみよう．
# 横軸用のnumpy配列を準備
x_step = 0.05
x_min = 0
x_max = x_step + 12
x = np.arange(x_min, x_max, x_step)
# 自由度5のカイ2乗分布を計算
y5 = chi2.pdf(x, 5)
# グラフを描く
fig, ax = plt.subplots(1, 1, figsize=(6, 4))
ax.plot(x, y5)
ax.set_xlim(x_min, x_max)
ax.set_ylim(0, 0.2)
ax.set_xticks(np.arange(0, 13, 2))
ax.set_yticks(np.arange(0.0, 0.25, 0.05))
ax.grid()
# plt.savefig('stat-chi2-5.png', dpi=300, facecolor='white')
plt.show()
様々な自由度でカイ2乗分布がどのような形状になるのかを確かめてみよう．
# 自由度1-4のカイ2乗分布を計算
y1 = chi2.pdf(x, 1)
y2 = chi2.pdf(x, 2)
y3 = chi2.pdf(x, 3)
y4 = chi2.pdf(x, 4)
# グラフを描く
fig, ax = plt.subplots(1, 1, figsize=(6, 4))
ax.plot(x, y1, label="n=1")
ax.plot(x, y2, label="n=2")
ax.plot(x, y3, label="n=3")
ax.plot(x, y4, label="n=4")
ax.plot(x, y5, label="n=5")
ax.legend() # 凡例を表示
ax.set_xlim(x_min, x_max)
ax.set_ylim(0, 0.6)
ax.set_xticks(np.arange(0, 13, 2))
ax.set_yticks(np.arange(0.0, 0.7, 0.1))
ax.grid()
# plt.savefig('stat-chi2-12345.png', dpi=300, facecolor='white')
plt.show()
上側 \(\alpha\) 点
統計学テキストの付録によく掲載されている自由度 \(n\) のカイ2乗分布の上側 \(\alpha\) 点の値を求めていみよう．この 上側 \(\alpha\) 点の値は下図の面積がちょうど \(\alpha\) になるような横軸の値 \(\chi_{n}^{2}(\alpha)\) を意味する．
上側 \(\alpha\) 点は chi2.isf(q, df) で求めることができる．例えばテキストの表にある「\(\alpha = 0.05\) の列」と「\(n = 5\) の行」にある「11.0705」という値は次のように求めることができる．
chi2.isf(q = 0.050, df = 5)
11.070497693516355