腑に落ちない人のための区間推定

区間推定は理解するのに色々とややこしいことが多く、私も計算方法とかは覚えたけど結局なんなんだっけ？となった記憶があります。

今回は、そんな区間推定の簡潔な理解に努めてみようと思います。

※初めて信頼区間を勉強している段階の人が読むと、知らなさそうな言葉が頻発しているかもと思ったので、そういう言葉には注釈を入れています。結構多くなってしまったので、雰囲気だけ知りたい方はあまり注釈を気にしなくても大丈夫です。

区間推定とは

区間推定というと、普通は信頼区間の話になるのですが、その意味をまどろっこしい日本語で説明することが多く、日本語の話をしているんだか数学の話をしているんだか分からなくなる時が多々あります。

ここでは、そういう細かい日本語のことは置いておいて、もっと前段のそもそも区間推定って何だっけというのを説明しようと思います。

区間推定とは、パラメータを区間で推定する方法です。

統計的推測にも種類¹があり、特にパラメータの値を推定することを点推定といいます。

例えば、$X=(X_1,X_2,...,X_n)$ を未知の平均 $\mu$ の正規分布から得られたランダム標本とした時、（未知の）真の平均 $\mu$ を標本平均 $\bar{X} = \frac{1}{n} \sum_{i=1}^{n}\,X_i$ によって推測するような行為は、まさに点推定です。

この標本平均は、色々と良い性質を持った推定量²ではありますが、$\mu$ を完全に正確に推定することができるわけではありません。³

そこで、真の平均を1点で推測するのではなく、範囲で推測した方が合理的なのではないか、ということで考えるのが区間推定です。

区間推定は、パラメータを区間で推定する方法なんだと理解できたら、今回の主目的はを果たしたことになります。ストックを押し、教科書の問題を解いて身につけましょう。お疲れさまでした。

勉強が進んだら、ストックからまた帰ってきて頂けると嬉しいです。

もうちょっと読んでやってもいいよという方は、続きを読んで頂けるとすんごく嬉しいです。。。

信頼区間

区間推定っていうのは、パラメータを点ではなく、範囲で推定する方法 だということが分かったと思います。ここでは、区間推定で一般的に使われている信頼区間について、定義と求め方を紹介します。

標本⁴ $X=(X_1,X_2,...,X_n)$ の分布が、統計モデル⁵ $M=\{P_{\theta}：\theta \in \Theta \subset R \,\}$⁶ に属するとします。統計量 $T_1(X)\,,\,T_2(X)$ が、

『 全ての $\theta \in \Theta$ で、$P_\theta\{\,T_1(X)≤\theta≤T_2(X)\,\}≥1-\alpha（0≤\alpha≤1）$ 』

を満たす時、閉区間 $[\,T_1(X)\,,\,T_2(X)\,]$を $\theta$ の $100×(1-\alpha)\,\%$ 信頼区間といいます。($1-\alpha$ のことを信頼係数といいます）

ポイントは $T_1(X)\,,\,T_2(X)$ が統計量、つまり確率変数だということです。

言い換えると、分布しているのは $\theta$ ではなく、区間 $[\,T_1(X)\,,\,T_2(X)\,]$ の方だということです。

手順で書いた方が分かりやすいかもしれません。⁷

＜区間推定の手順＞
①モデル $P_{\theta}$ を決める。
②信頼水準 $1-\alpha$ を決める。
③モデル $P_{\theta}$ と観測値 $x$ を使い、モデル $P_{\theta}$ のパラメータ $\theta$ が、確率 $1-\alpha$ で含まれる区間 $[\,T_1(X)\,,\,T_2(X)\,]$⁸ を推定する。

例：正規母集団を例に。

ここでは、人工的に生成したデータを使って、上に書いた手順①～③を実行していこうと思います。

現実の問題では、観測値 $x$ がどんな確率分布から生成されたものなのかを知ることはできないのですが、ここでは奇跡的に観測値 $x$ を生成している確率分布が、正規分布 $N(\mu\,,\,25)$ だと分かっていることにしましょう。（真の分散 $\sigma^2=25$ だということも分かっている状況です。⁹）

今、この正規母集団から独立に400個の観測値を得たとします。
（（未知の）真の平均だけは神様だけが知っていて、その値は $\mu=10$ です。）

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(seed=42)
X = np.random.normal(loc=10, scale=5, size=400)
plt.hist(X, bins=20)

この観測値を使って、真の平均 $\mu$ を区間推定します。これを先ほどの手順 ①～③ に沿って進めていきます。

① モデルを決める。

今回は神のお告げかなんかで、正規分布だと分かっちゃってるので、素直に $N(\mu, 25)$ だと決めましょう。

②信頼水準を決める。

次に ②信頼水準をきめるのですが、慣例的によく用いられる $1-\alpha=0.95$ にしましょうか。

③信頼区間を求める。

そして、いよいよ ③ 信頼区間を求めます。

まず、$X_1,X_2,...,X_n$ が独立に正規分布 $N(\mu,\sigma^2)$ に従う時、正規分布の再生成から、

\sum_{i=1}^n X_i \,\,\,～\,\,\,N(n\mu\,,\,n \sigma^2)

であるから、

\bar{X}=\frac{1}{n} \sum_{i=1}^n X_i \,\,\,～\,\,\,N(\mu\,,\,\frac{\sigma^2}{n})

となり、標本平均 $\bar{X}$ は $N(\mu, \frac{\sigma^2}{n})$ に従うということが分かります。

なので、この標本平均の標準化 $\frac{\bar{X}-\mu}{\frac{\sigma}{\sqrt{n}}}$ が $N(1,0)$ に従うことも分かります。

これが、

P_{N(0,1)}\,\left \{T_1≤\frac{\bar{X}-\mu}{\frac{\sigma}{\sqrt{n}}}≤T_2 \right \}=0.95

となればいいですね。

今、$\frac{\bar{X}-\mu}{\frac{\sigma}{\sqrt{n}}}$ の従う分布は $N(1,0)$ だと分かっているので、$T_1$ と $T_2$ はそれぞれ正規分布表¹⁰から読み取ればよく、$T_1=-1.960\,,\,T_2=1.960$ だということが分かります。

あとは愚直に式変形していけばよく、
（最初は、ほほーんって流さないで、手を動かした方がいいと思うので、冗長すぎる式変形をおいておきますね）

\begin{align}
P_{N(0,1)}\,\left \{-1.960≤\frac{\bar{X}-\mu}{\frac{\sigma}{\sqrt{n}}}≤1.960 \right \}=0.95 \\
\\
\\
P_{N(0,1)}\,\left \{-1.960\frac{\sigma}{\sqrt{n}}≤\bar{X}-\mu≤1.960\frac{\sigma}{\sqrt{n}} \right \}=0.95\\
\\
\\
P_{N(0,1)}\,\left \{-\bar{X}-1.960\frac{\sigma}{\sqrt{n}}≤-\mu≤-\bar{X}+1.960\frac{\sigma}{\sqrt{n}} \right \}=0.95\\
\\
\\
P_{N(0,1)}\,\left \{\bar{X}+1.960\frac{\sigma}{\sqrt{n}}≥\mu≥\bar{X}-1.960\frac{\sigma}{\sqrt{n}} \right \}=0.95\\
\\
\\
\end{align}

ということになり、信頼区間は $\left [\bar{X}-1.960\frac{\sigma}{\sqrt{n}}\,,\,\bar{X}+1.960\frac{\sigma}{\sqrt{n}} \right]$ だということが分かりました。

$n$ と $\sigma$ は分かっているので、 $\bar{X}$ を求めれば信頼区間も求まりますね。

mu = np.mean(X)
mu

=====実行結果=====
10.034189972943238
==================

というわけで、あとはさっきの式に代入して、

$\left [10.034-1.960×\frac{5}{\sqrt{400}}\,,\,10.034+1.960×\frac{5}{\sqrt{400}} \right]$

頑張って計算すると（機械にできることは、機械にやらせましょう）

平均 $\mu$ の $95\%$ 信頼区間が

$\left [\,9.544\,,\,10.524 \, \right]$

と推定されました。

計算自体は簡単ですね。

区間推定量の推移

信頼区間は、推定量なのでもちろん分布しています。ここでは、先ほどの正規母集団を使って、サンプルサイズごとの区間推定量の推移を見てみます。

まずは、サンプルサイズごとに信頼区間を計算しておきます。

N = 200 #サンプルサイズ最大 
ymin_list = []
ymax_list = []

for i in range(1, N):
    X = np.random.normal(loc=10, scale=5, size=i)
    # 標本平均
    mu = np.mean(X)
    # 上側信頼限界
    ymin_list.append(mu - 1.96 * 5 / np.sqrt(i))
    # 下側信頼限界
    ymax_list.append(mu + 1.96 * 5 / np.sqrt(i))

んでもって、図をカキカキしていきます。

SampleSize = np.arange(1, N, 1)

fig = plt.figure()
plt.xlim(0, N)
plt.ylim(min(ymin_list), max(ymax_list))

plt.plot(SampleSize, ymax_list, color="magenta", label="upper confidence limit")
plt.plot(SampleSize, ymin_list, color="Blue", label="lower confidence limit")
plt.fill_between(SampleSize, ymin_list, ymax_list, alpha=0.2, facecolor="yellow")

plt.hlines(y=10, xmin=1, xmax=N, color="red")
plt.xlabel("SampleSize")
plt.ylabel("X")
plt.annotate('True Mean', xy=(N/2, 10), xytext=(2/3*N, 2.5), fontsize=15, arrowprops={'width':0.4,'headwidth':7,'color':'black'})
plt.legend(loc='upper right', borderaxespad=1, fontsize=10)
plt.show()

マゼンタが上側信頼限界で、青が下側信頼限界です。
このマゼンタと青で囲まれた黄色の領域が、サンプルサイズごとの $\mu$ の推定値です。

また、真の平均を赤線で描き加えています。

この図を見ると、確かに信頼区間が確率的に変化していそうだということが分かると思います。

まとめ

くどいようですが、区間推定は、パラメータを区間で推定する方法ですということを覚えて頂ければ、私は満足です。

今回記事を書いていての発見

区間推定は発想自体は単純なのに、なんとなく分かりにくいのは何でだろうと思っていましたが、ググりまくって分かりました。たぶんこんな図をよく見るからなんじゃないかなと。

赤線は標本平均、青で塗りつぶされているのが信頼区間です。

これを見て、なんとなく違和感持ちませんか？

私は、「間推定の話をしているのに、なんでわざわざ点推定と区間推定を並べるの？」という気持ちになりました。区間推定の話をしているのであれば、これだけでいいような気がします。

線形モデルとかの話でも、親切にも線と信頼区間を一緒に書いてくれていることも多いですが、帯だけ書いてくれれば誤解も少ないのかなとも思ったりします。早い話が、なんの話をしているのか明確にしてくれ！ということなのですが。。。

そう思ったのが、私だけだったらごめんなさいです。

ちょっと注意しておきたいこと

今回は明示的に①モデルを決めると書きましたが、現実問題では、これが結構重要なんですよね。

なぜなら信頼区間を作る時に、①で決めたデータを生成していると考えられるモデル（確率分布）を使っているので、このモデルが違うと、得られる信頼区間も違ってきてしまうからです。そう考えると、確かにこの『①モデルを決める』という作業が極めて重要だということが分かるかなと思います。

今回の例では、データの生成分布が正規分布だと分かってしまっていたので、疑いもなく推測をしているわけですが、現実では真のモデルなんて知りようがありません。

なので、分析者は決めたモデルを明示的に伝える必要がありますし、分析した結果を見る側はどういうモデルを想定したのかを聞いておく必要があるということなんですね。

統計ってちょっとめんどくさいですね。

（んでもまぁ現実には、そんなことまで気にする人はすごく珍しいので、むしろそのことの方が問題なんじゃないかなとも思ったり、思わなかったり、ちょっぴりさびしかったり、、、）

★参考文献★

[1] 野田,宮岡：数理統計学の基礎(1992)
[2] 薩摩：確率統計(1989)

---

1. 統計的推測の種類：パラメータの値を推定することを点推定、区間で推定することを区間推定、パラメータの属しているパラメータ空間の部分集合を推定することを仮説検定といいます。 ↩

2. 良い推定量：一般に、統計的な推測において、推測方法は無数にあります。母平均を点推定するという方法だけでも、例えば、よく使う標本平均や中央値、刈り込み平均のようなロバストな推測法、最大値と最小値を足して2で割る(中点値)なんかも考えられます。こうやって色々な推測方法があると、どれを使っていいのかよく分からなくなってしまいますよね。そこで、どの推定方法が、どういう意味で良い性質を持っているのかというのを考えてあげるんですね。 ↩

3. 連続確率変数が1点をとる確率は0なので 、標本平均 $\bar{X}$ が真の平均 $\mu$ を完全に正確に推測する確率は $P(\bar{X}=\mu)=0$ です。ただし、これは$\bar{X}$ が真の平均 $\mu$ をピッタリと当てることが起こらないという意味ではありません。確率ってややこしいですね。 ↩

4. 標本：統計学では確率現象の考察を行います。そのために観測値 $x$ を確率変数 $X$ の実現値とみなし、$X$ の確率分布を推測するということをします。このように、観測値 $x$ を確率変数 $X$ の実現値だとみなす時、確率変数 $X$ のことを標本と呼びます。早い話が統計学では、確率変数のことを標本って呼んでいますと考えてください。 ↩

5. 統計モデル：標本 $X$ の真の確率分布は神のみぞ知ることで、私達は知ることができません。しかし、「これは真の確率分布じゃないかな～」という候補は考えることができます。この候補と考えられる確率分布の集まりのことを標本 $X$ の統計モデルと呼びます。 ↩

6. $P_\theta$：統計モデルに属する各確率分布を区別するために、各確率分布をパラメータ $\theta$ で特徴付けて$P_{\theta}$ と表します。また、すべてのパラメータの集まりを $\Theta$ と表し、これをパラメータ空間と呼びます。 ↩

7. もちろん、この手順はあくまでもよく使われる区間推定の方法だというだけで、パラメータを区間で推定していれば、区間推定といいます（よね？）。例えば、意味のないことですが、$(-\infty, \infty)$ としてもよいわけです。 ↩

8. ちなみに $[\,T_1(X)\,,\,T_2(X)\,]$ のことを区間推定量と呼ぶのですが、観測値から計算した $[\,T_1(x)\,,\,T_2(x)\,]$ のことは何て呼ぶんですかね？区間推定値？？ ↩

9. 教科書とかを読んでいて、真の分散が分かっている状況なんてそんなのありえないのに、何でそんな設定をするの？と思う方は多いと思います。私も勉強し始めた時にはそう思っていましたので、そう思うこと自体は別におかしいことではないと思います。たぶん推定の順番が平均→分散という流れなのに、先に分散が分かっているという状況が気持ち悪いのだと思います。しかし、考えやすい状況から出発して、徐々に制限を緩めていきながら難しいこと考えていくという流れは、すごく自然なことなのです。（そもそも正規母集団を仮定していることも、現実問題では気持ち悪いはずなのですが、、、） ↩

10. 正規分布表って、見たこと、、、ありますよね、、、？？ ↩