目次
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