区間推定の理論と実践

統計学における推定の考え方は仮説検定と双璧を成す重要な理論であり、初学者がまず目指すべき最高到達点である。推定は強力な手法であるが故に、恐らく大抵の教科書では記載があるだろうし、ググればいくらでもWebページがヒットする。

そんな状況下ではあるが、推定の中でも特に区間推定は個人的に興味をそそられる内容であり、ぜひともブログを書き起こしてみたくなった。そこで、本稿ではできるだけ二番煎じにならないよう、私が得た区間推定に対する理解について説明してみたいと思う。

本稿では正規母集団から抽出した標本データを用いた母平均の区間推定のみをターゲットとする。また、適当なサンプルデータを用いて実際に区間推定を行うことで、区間推定の理論に対して確かな実感を得てみたいと思う。

なお、本稿を執筆するに辺り、全体的に[1][2]を参考にした。

母平均の区間推定

母分散が既知の場合

正規母集団 N(\mu,\ \sigma^2)に従う独立な確率変数 X_1,\ X_2,\ \cdots .\ X_nの標本平均 \overline{X}を考える。このとき、 \overline{X}自体も確率変数となっており、これが従う確率分布は N(\mu,\ \sigma^2/n)となる。また、これを標準化した確率変数 Z = \frac{\overline{X} - \mu}{\sigma/\sqrt{n}}は、当然だが標準正規分布 N(0,\ 1)に従う。

今、母分散 \sigma^2は分かっているが、母平均 \muの値が不明な状況を考える。もちろん、実際にはこんなことはまずあり得ないのだが、このケースを考えることは後々の議論において意味がある。

この時、 \mu = E(\overline{X})という関係が成り立つので、標本平均をそのまま母分散だと考えてもあながち間違いではない。しかし、当然だが \overline{X}の値にはぶれがあり、どんなサンプルを抽出するかによって値が変わってしまう。しかも、極端なケースでは \overline{X}はとんでもなく大きくなったり小さくなったりすることがあり得る。そのため、 \overline{X}の値だけをみても、母平均について何も確実なことは言えない。

では、何か確率的な議論はできないだろうか?例えば「95%の確率で \muはAからBの間に入ります」というようなことが言えないだろうか?

結論から言うとこれは可能である。標準正規分布において、それより大きい値を取る確率が \alpha\ (0 < \alpha < 1)となるような点 ( 100 \alphaパーセント点) を Z_\alphaと表すと、先ほど導入した確率変数Zは確率 1 - \alphaで以下の不等式を満たす。

 \displaystyle{
 -Z_{\alpha/2} < Z < Z_{\alpha/2}
}

2点注意事項がある。まず、これは正規分布が左右対称だから言えることである。さもなければ、  -Z_{\alpha/2}の部分はそれより小さい値を取る確率が \alpha/2となるような点を別途求めて、それと置き換えてやる必要がある。

次に、 \alphaを2で割っている理由だが、これは上側の \alpha/2と下側の \alpha/2を両方排除し、残りが 1 - \alphaとなるようにするためである。

この不等式を変形すれば、以下のように \muが確率 1 - \alphaで取る値の範囲が分かる。

 \displaystyle{
\overline{X} - \frac{Z_{\alpha/2} \sigma}{\sqrt{n}} <  \mu < \overline{X} + \frac{Z_{\alpha/2}\sigma}{\sqrt{n}}
}

母分散が未知の場合

先ほども述べた通り、母平均は未知なのに母分散は分かっているということは普通はあり得ない。実際に得られるのは不偏分散 s^2だけというケースがほとんどである。不偏分散の定義式を以下に示す。

 \displaystyle{
s^2 = \frac{1}{n-1} \sum_{i=1}^n (X_i - \overline{X})^2
}

不偏分散に対しては \sigma^2 = E(s^2)という関係式が成り立つので、これを母分散の代わりに用いることはそれほど不自然ではない。しかし、やはりサンプル数が少ない場合には母分散とかけ離れた値を取ることがあり得る。そのため、先ほど考えた確率変数Zの式に対して単に \sigma = sを代入すると、これはもはや正規分布には従わない。

では、これはどのような分布になるのだろうか?それを調べるために、Zの式に \sigma = sを代入した式を以下のように変形してみよう。

 \displaystyle{
\begin{eqnarray}
t &=& \frac{\overline{X} - \mu}{\frac{s}{\sqrt{n}}} \\
&=& \frac{\overline{X} - \mu}{\frac{\sigma}{\sqrt{n}} \cdot \frac{s}{\sigma}} \\
&=& \frac{Z}{\frac{\sqrt{\frac{1}{n-1} \sum_{i=1}^n (X_i - \overline{X})^2}}{\sigma}} \\
&=& \frac{Z}{\sqrt{\frac{1}{n-1} \sum_{i=1}^n \left(\frac{X_i - \overline{X}}{\sigma}\right)^2}}
\end{eqnarray}
}

ここで、分母に \sum_{i=1}^n \left(\frac{X_i - \overline{X}}{\sigma}\right)^2という式が現れた。この式自体も確率変数であり、これは標準正規分布に従う確率変数の二乗和を表す。このような確率変数が従う確率分布は自由度n-1の \chi^2分布と呼ばれる。

まとめると、標準正規分布に従う確率変数Zと自由度n-1の \chi^2分布に従う確率変数Wによって、元々の確率変数tは以下のように表される。

 \displaystyle{
t = \frac{Z}{\sqrt{\frac{W}{n-1}}}
}

このような確率変数tが従う確率分布は自由度n-1のt分布と呼ばれる。結論として、母分散未知のケースにおける母平均の推定にはt分布を使用する必要があることが分かった。

自由度n-1のt分布を t(n-1)と書くと、母分散既知のケースと同様の計算を行うことで、 \muは確率 1 - \alphaで以下の不等式を満たすことが分かる。

 \displaystyle{
\overline{X} - \frac{t_{\alpha/2}(n-1) s}{\sqrt{n}} <  \mu < \overline{X} + \frac{t_{\alpha/2}(n-1) s}{\sqrt{n}}
}

ただし、 t_{\alpha/2}(n-1)は自由度n-1のt分布の 100 \alpha/2パーセント点を表す。ここでもt分布が左右対称であることを利用している。

サンプル数に関する注意

t分布は自由度が大きくなると正規分布に収束していく。そのため、t分布を用いる必要があるのはサンプル数が少ない場合に限る。実用上はサンプル数が30以上あれば正規分布を使用しても差し支えないようである[1]。

区間推定の実践

ここまでに述べたことを実践してみよう。統計処理をするにはpythonかRが簡単そうであるが、Android上で環境を構築できたのでpythonを使うことにする[3]。また、データセットとしてIrisは見飽きたので、ここではAbarone (アワビ) を使用している[4]。

Abaroneデータセットには、1列目に性別、3列目に直径の情報が入っている。今回はオスのアワビの直径データについて、母平均の区間推定を行ってみた。サンプル数は敢えて少なめに10とし、t分布を使った場合と正規分布を使った場合で比較を行った。もちろん、正しいのはt分布を使う方である。

以下に今回作成したソースコードを示す。作成の際には[5]を参考にさせて頂いた。日本語コメントは本稿のために後から付したものであり、このままでは実行できないので注意されたい。

import pandas as pd
from scipy import stats
import math

# 全部の列に名前を付けるのは諦めた...
abalone = pd.read_csv('file:///storage/emulated/0/Download/abalone.data', 
                 header=None, names=['sex', 'length', 'diameter', 'height', 'a', 'b', 'c', 'd', 'e'])

# オスの情報だけ抽出
male = abalone[abalone['sex'] == 'M']

n = 10
# 自由度n-1のt分布の95%信頼区間を取得
tleft, tright = stats.t.interval(alpha=0.95, df=n-1)
# 正規分布の95%信頼区間を取得
nleft, nright = stats.norm.interval(alpha=0.95)

mean = male.iloc[:n]['diameter'].mean()
# デフォルトで不偏分散の平方根が計算される
stddev = male.iloc[:n]['diameter'].std()

# 本文中の式と若干符号が異なるので注意
tleft = mean + tleft*stddev/math.sqrt(n)
tright = mean + tright*stddev/math.sqrt(n)

nleft = mean + nleft*stddev/math.sqrt(n)
nright = mean + nright*stddev/math.sqrt(n)

print("mean = {:.3f}, stddev = {:.3f}".format(mean, stddev))
print("t-dist: left = {:.3f}, right = {:.3f}".format(tleft, tright))
print("n-dist: left = {:.3f}, right = {:.3f}".format(nleft, nright))

実行結果を以下に示す。

mean = 0.339, stddev = 0.046
t-dist: left = 0.306, right = 0.372
n-dist: left = 0.310, right = 0.368

実行結果から分かるように、正規分布を使った場合は区間が狭めに出力されてしまう。そのため、正規分布での結果を元に考察を行うと判断を誤る恐れがある。

最後に特大の注意事項がある。今回はオスのアワビの直径データが正規分布に従うと仮定して分析を行ったが、本来であればこの仮定の妥当性は良く確かめる必要がある。誤った仮定を元に一生懸命計算したとしても、その結果には何の意味もないどころか、誤った帰結へと導かれる危険性がある。

データが正規分布に良く当てはまっているかどうかを判断する方法について私には知見がない。軽くググってみたらいろいろと情報がありそうだったので、興味のある方は調べてみると面白いかもしれない。

その他の話題

本稿では正規母集団から抽出したサンプルから母平均の区間推定を行った。この他にも、例えば母分散の区間推定を行う方法などは有名である[1]。また、実用上のタスクとしては正規母集団以外の母集団を考えることもあるだろう。しかし、どんな場合であれ、基本的な区間推定の考え方は本稿に述べたものと大きく異なることはないであろう。

まとめ

以上、区間推定に関する理論と実践について述べた。統計学は応用先が非常にたくさんある分野なので、理解を深めておいて損はない。また機会があれば、多変量解析などの進んだ内容も勉強し直したい。