目次
Python で統計学
標準正規分布
正規分布 \(\mathcal{N}\left(\mu, \sigma^2\right)\) は統計学の基本的な分布である．正規分布の密度関数は \begin{eqnarray*} f(x) = \frac{1}{\sqrt{2\pi} \sigma} \exp \left( - \frac{(x-\mu)^2}{2\sigma^2}\right) \end{eqnarray*} であり，ここで \(\mu\) が平均，\(\sigma^2\) が分散である．さらに \(\mu = 0\)，\(\sigma^2 = 1\) とした正規分布は標準正規分布と呼ばれる． \begin{eqnarray*} f(x) = \frac{1}{\sqrt{2\pi}} e^{-x^2/2} \end{eqnarray*}
標準正規分布の密度関数をプロットすると次のようになる．
ライブラリのインポート
Python では Scipy をインポートすることで，正規分布関する様々な値を取得できる．まずは，ライブラリをインポートする．
import numpy as np
from scipy.stats import norm
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')
確率密度関数
確率密度関数の値は norm.pdf( ) で求めることができる．例として，\(x = 0, 1, 2\) での確率密度関数の値を求めてみる．なお，norm.pdf(0) でも norm.pdf(x = 0) でもよい．
norm.pdf(0)
0.3989422804014327
norm.pdf(x = 0)
0.3989422804014327
norm.pdf(x = 1)
0.24197072451914337
norm.pdf(x = 2)
0.05399096651318806
確率密度関数をプロットしてみよう．
# 横軸用のnumpy配列を準備
x_step = 0.05
x_min = -5
x_max = x_step + 5
x = np.arange(x_min, x_max, x_step)
# 標準正規分布
y = norm.pdf(x)
# グラフを描く
fig, ax = plt.subplots(1, 1, figsize=(6, 4))
ax.plot(x, y)
ax.set_xlim(-4.5, 4.5)
ax.set_ylim(0, 0.45)
ax.set_xticks(np.arange(-4, 5, 1))
ax.set_yticks(np.arange(0.0, 0.5, 0.1))
ax.grid()
# plt.savefig('normal.png', dpi=300, facecolor='white')
plt.show()
分布関数
標準正規分布の分布関数の値を求めるには norm.cdf( ) を利用すると良い．なお，分布関数は
\begin{eqnarray*}
\Phi(x) = \frac{1}{\sqrt{2\pi}} \int_{-\infty}^{x} e^{-t^2/2}dt
\end{eqnarray*}
で定義されるので，\(-\infty\) から \(x\) までの面積を求めていることになる．
norm.cdf(x = 0)
0.5
norm.cdf(x = 1)
0.8413447460685429
norm.cdf(x = 2)
0.9772498680518208
標準正規分布の範囲
正規分布に従うときにある値がある範囲にある確率を求めることは頻繁に行われる．ある値が \(a\) から \(b\) の範囲になる確率は次の図で示すことができる．
これは分布関数値の差によって求めることができる．例えば，\((-1, 1)\) の範囲に含まれる確率は次のように求めることができる．
a = -1
b = 1
norm.cdf(x = b) - norm.cdf(x = a)
0.6826894921370859
同様に，\((-2, 2)\) や \((-3, 3)\) の範囲に含まれる確率を求めてみよう．
a = -2
b = 2
norm.cdf(x = b) - norm.cdf(x = a)
0.9544997361036416
a = -3
b = 3
norm.cdf(x = b) - norm.cdf(x = a)
0.9973002039367398
標準正規分布表の値
上の方法を応用することで統計学テキストの付録に必ず掲載されている「標準正規分布表」の値を求めることができる．この表に掲載される値の数式は \begin{eqnarray*} I(z) = \frac{1}{\sqrt{2\pi}} \int_{0}^{z} e^{-x^2/2}dt \end{eqnarray*} であり，次の図の面積を意味する．
つまり，これは上の標準正規分布の範囲において，\(a = 0\)，\(b = z\) と置いたものに他ならない．よって \(z = 1.00\) のときに，\(I(1.00) = 0.34134\) であるという標準正規分表に掲載されている値は次のように求めることができる．
norm.cdf(x = 1) - 0.5
0.3413447460685429
なお，上の 0.5 という値は，次の結果を用いていることに注意する．
norm.cdf(x = 0)
0.5
上側 \(\alpha\) 点の値
一般に \(0 < \alpha < 1 \) に対して， \begin{eqnarray*} \Pr\{Z \geq z\} = \frac{1}{\sqrt{2\pi}} \int_{z}^{\infty} e^{-x^2/2}dt = \alpha \end{eqnarray*} となるような \(z = z_\alpha \) 標準正規分布の上側 \(\alpha\) 点という．これは下図の面積が \(\alpha \) になる横軸の値 \(z_\alpha \) になる．
統計学の教科書には次の値がよく登場する．
| \(\alpha \) | 10% | 5% | 2.5% | 1% | 0.5% |
|---|---|---|---|---|---|
| 上側 \(\alpha \) 点 \((z_\alpha) \) | 1.2816 | 1.6449 | 1.9600 | 2.3263 | 2.5758 |
| 両側 \(\alpha \) 点 \((z_\alpha/2) \) | 1.6449 | 1.9600 | 2.2414 | 2.5758 | 2.8070 |
例えば，\(\alpha = 5\%\) の上側 \(\alpha \) 点は norm.isf( ) を使って次のように求めることができる．
norm.isf(q = 0.05)
1.6448536269514729
さらに 両側 \(\alpha \) 点 は次のようにすれば良い．
norm.isf(q = 0.05/2)
1.9599639845400545