中心極限定理は強力なツール

　不偏推定量の項で私たちが統計量を計算するのは、母集団のことを知るためだ、と述べた。これまでの項では、その母集団を正規分布として解析してきた。確かに、この世のデータのばらつき方には正規分布が現れることが多い。しかしながら、実際には正規分布を仮定できないような母集団を持つデータもたくさん存在している。こんなとき、今まで述べてきた手段が必ずしも役立たないことがある。例えば、t検定は母集団が正規分布を仮定できないときには、統計の精度が落ちてしまい、正しい判定を下すことができなくなってしまう。ゆえに、母集団の性質に合わせた統計ツールを使うように心がけなければならない。

　ところが、母集団の分布がどんなものであっても、平均値の推定を行える方法が存在する。中心極限定理と呼ばれる（中二ごころくすぐる名前の）強力なツールだ。本稿では、中心極限定理がどんなものかを解説し、どんな使い方をするのかを紹介しよう。この定理のおかげで、母集団のことがわからなくても、まずは平均と分散を確認することの妥当性が保証されているのである。

中心極限定理とは？

　中心極限定理の形にはいくつかあるが、もっとも利用しやすいものは、「母集団の平均がμ、標準偏差σのとき、この母集団から生じたサンプルサイズnの標本の平均値は平均がμ、標準偏差σ/√nの正規分布に従う」というものだ。標本平均の分布における標準偏差σ/√nは、標本平均の標準誤差と呼ばれる（母集団が正規分布の場合の証明は追記）。

ここで重要なのは、母集団が正規分布と限らずとも、標本平均の分布が正規分布となる点だ。正規分布は、世の中のデータでの出現頻度が高いだけでなく、その性質が良く知られている。例えば、平均μ、標準偏差σとわかっている正規分布の場合、平均μから±σ、2σ、3σの範囲に、それぞれ約68%、95%、99%のデータが含まれていることになる。

標本平均が平均がμ、標準偏差σ/√nの正規分布に従うということは、つまり、標本平均を計算すると95%の確率で、μ ± 2(σ/√n)の範囲に収まっているということになる（※突っ込んだ話をすると、頻度主義統計の観点からはこの表現は正しくないことを明記しておく。興味があれば、本稿の最後に目を通してほしい。が、ここでは理解を優先するためにこのように表現しておく）。逆に、標本平均 ± 2(σ/√n)の範囲に95%の確率で母平均μが存在するともいえる。実は、±2(σ/√n)の範囲だと95%よりもやや多いため、ぴったり95%の範囲にするために2ではなく1.96という値を使うこともある。このような、母平均が存在する範囲を信頼区間と呼ぶ。例えば、95%の確率の範囲なら、95%信頼区間という。

　例えば、平均20、標準偏差5の母集団（正規分布でなくてよい）からサンプルサイズ30のデータが得られたとする。この時、標本平均は平均20、標準偏差5/√30 = 0.912の正規分布に従い、20 - 1.96 × 0.912 = 18.21～20 + 1.96 × 0.912 = 21.79の範囲に、95%の標本平均が存在する、という具合だ。

※分散が未知ならt分布に従う

　実際には、母平均と母標準偏差はわからないので、標本から計算した標本平均と不偏分散の平方根で代用する。不偏分散で代用するときは標本平均は自由度n-1のt分布に従うことが知られている（分散の既知/未知で計算が変わる理由はz検定参）。標本平均が13、不偏分散の平方根が2.5、サンプルサイズが20であれば、自由度19のt分布の中心から95%の範囲を与える2.09を使うことで、13 - 2.09 × 2.5/√20 = 11.83～13 + 2.09 × 2.5/√20 = 14.17の範囲に、95%くらいの確率で母平均が存在すると言える。母集団のことが全く分からない状況でも、標本平均と不偏分散を計算することで、11.83～14.17の範囲に母平均が存在するんじゃないか、と言えてしまった。これは心強い主張だ。

　もう一つ、中心極限定理の重要な性質が、標本平均の分布が標準偏差σ/√nの正規分布に従う、ということなので、サンプルサイズnを大きくすることで標準偏差が小さくなる、つまり推定精度を上昇できるという直感的にも正しい、理論的保証を与えている点だ。例えば、先ほどの標本平均が13、不偏分散の平方根が2.5の場合、サンプルサイズを40にすると、13 - 2 × 2.5/√40 = 12.21～13 + 2 × 2.5/√40 = 13.79の範囲に狭まる。ただし、推定精度の上昇は平方根がかかっているため、推定精度を10倍にしようと思ったら、100倍のサンプリングが必要になる。どれくらいの精度が求められているかは、コストと相談しよう。

中心極限定理をシミュレーションする

　では、実際にRを使ったシミュレーションで、中心極限定理を示してみよう。まず初めに、母集団が正規分布の場合である。

------------------------------------------------------

library(plotn)

#このライブラリは私オリジナルの作図ライブラリである。

#もし、plotnライブラリが欲しい場合は、

#以下のようにコマンドを打つ（#は消す）。devtoolsライブラリが必要。

#library(devtools)

#install_github("bugplant/plotn")

n <- rnorm(10000, mean = 50.5, sd = 28.86751) #正規分布に従う母集団を生成

#rep(0:100, 100)の一様分布に合わせた

histn(n, xlab = "value", ylab = "頻度") #正規分布の図示

mean(n)

## [1] 50.67655

sd(n)

## [1] 28.50751

------------------------------------------------------

この正規分布から、サンプルサイズ20のサンプリングを行い、標本平均を算出することを10000回繰り返す。

------------------------------------------------------

c <- NULL

for (i in 1:10000){

  c <- c(c, mean(sample(n, 20))) #20個ランダムサンプリングし、平均値を算出することを10000回繰り返す。

}

histn(c, xlab = "標本平均", ylab = "頻度")

mean(c)

## [1] 50.78209

sd(c) #標本平均の分布の標準偏差

## [1] 6.414747

sd(n)/sqrt(20) #母集団から計算した標準誤差。sd(c)と近い値になることを確認

## [1] 6.374474

sum(c < mean(c) + 2 * sd(n)/sqrt(20) & c > mean(c) - 2 * sd(n)/sqrt(20)) #95%信頼区間に入っている標本平均の数

## [1] 9539

------------------------------------------------------

確かに、標本平均の分布は平均50.5を中心とした釣り鐘型分布となる。また、標本平均の分布の標準偏差は、シミュレートして作成した標本平均から標準偏差の関数を使って算出した場合と、母集団の標準偏差とサンプルサイズで計算した場合で近い値となることも確認できる。また、平均±2×標準誤差で表される95%信頼区間に含まれる標本平均の数が10000個中約9500個、つまり95%に相当することもわかるだろう。母集団が正規分布ならば、標本平均も正規分布に従うだろう、というのは想像に難くない。では、こんどは母集団が正規分布ではない場合だ。

　母集団が、一様分布ならどうだろう。一様分布は、乱数がすべて等しい確率で出現する分布だ。つまり、まっ平らな分布となり、釣り鐘型分布とはかけ離れた様相だ。コードを書けば以下のようになる。

------------------------------------------------------

u <- rep(1:100, 100) #一様分布に従う母集団を生成

histn(u, breaks = 20, xlab = "value", ylab = "頻度")

mean(u)

## [1] 50.5

sd(u)

## [1] 28.86751

------------------------------------------------------

この一様分布から、同様にサンプルサイズ20のサンプリングを行い、標本平均を算出することを10000回繰り返す。

------------------------------------------------------

c <- NULL

for (i in 1:10000){

  c <- c(c, mean(sample(u, 20))) #20個ランダムサンプリングし、平均値を算出することを10000回繰り返す。

}

histn(c, xlab = "標本平均", ylab = "頻度")

mean(c)

## [1] 50.59695

sd(c) #標本平均の分布の標準偏差

## [1] 6.422352

sd(n)/sqrt(20) #母集団から計算した標準誤差。sd(c)と近い値になることを確認

## [1] 6.454972

sum(c < mean(c) + 2 * sd(n)/sqrt(20) & c > mean(c) - 2 * sd(n)/sqrt(20)) #95%信頼区間に入っている標本平均の数

## [1] 9544

------------------------------------------------------

母集団が正規分布とはかけ離れた形であっても、標本平均の分布は、母集団が正規分布のとき時と同様、平均50.5を中心とした釣り鐘型分布となる。標準誤差や95%信頼区間の性質も同様だ。これはなかなか感動できるのではないだろうか？

　どんどんいこう。母集団が、二峰分布ならどうだろう。正規分布は一山形の分布出会ったが、今度は二つの山が存在するような分布ということだ。このような分布の場合、出やすい値が山近辺の2か所存在することになる。コードを書けば以下のようになる。今回は平均35、標準偏差11の正規分布と平均81.5、標準偏差28の正規分布を2:1で混ぜ合わせた分布にしてみた。つまり、この分布では35や81.5近辺の値が出やすいということになる。

------------------------------------------------------

m <- NULL

for(i in 1:10000){

  t <- rbinom(1,1,1/3)

  if(t == 0){

    m <- c(m, rnorm(1, 35, 11))

  } else {

    m <- c(m, rnorm(1, 81.5, 28))

  }

}

histn(m, xlab = "value", ylab = "頻度")

mean(m)

## [1] 50.29333

sd(m)

## [1] 28.39449

------------------------------------------------------

同様にサンプルサイズ20のサンプリングを行い、標本平均を算出することを10000回繰り返す。

------------------------------------------------------

c <- NULL

for (i in 1:10000){

  c <- c(c, mean(sample(m, 20))) #20個ランダムサンプリングし、平均値を算出することを10000回繰り返す。

}

histn(c, xlab = "標本平均", ylab = "頻度")

mean(c)

## [1] 50.31267

sd(c) #標本平均の分布の標準偏差

## [1] 6.260418

sd(n)/sqrt(20) #母集団から計算した標準誤差。sd(c)と近い値になることを確認

## [1] 6.349201

sum(c < mean(c) + 2 * sd(n)/sqrt(20) & c > mean(c) - 2 * sd(n)/sqrt(20)) #95%信頼区間に入っている標本平均の数

## [1] 9585

------------------------------------------------------

やはり同様の性質を示す。もうここまでくれば中心極限定理がかなり強い性質を持っていることがわかるだろう。

　では、最後に極めつけとして、母集団がたった3つの値しかとらない離散分布ならどうだろう。この場合、生じる値は3つしか存在せず、この3つの数字の組み合わせだけで標本平均が決定する。コードを書けば以下のようになる。今回は20、65、88だけしか存在しない分布にしてみた。それぞれが45%、30%、25%の確率で出現する。

------------------------------------------------------

d <- c(rep(20, 4500), rep(65, 3000), rep(88, 2500)) #離散分布に従う母集団を生成

histn(d, xlab = "value", ylab = "頻度")

mean(d)

## [1] 50.5

sd(d)

## [1] 28.86751

------------------------------------------------------

同様にサンプルサイズ20のサンプリングを行い、標本平均を算出することを10000回繰り返す。

------------------------------------------------------

c <- NULL

for (i in 1:10000){

  c <- c(c, mean(sample(d, 20))) #20個ランダムサンプリングし、平均値を算出することを10000回繰り返す。

}

histn(c, xlab = "標本平均", ylab = "頻度")

mean(c)

## [1] 50.59678

sd(c) #標本平均の分布の標準偏差

## [1] 6.484737

sd(n)/sqrt(20) #母集団から計算した標準誤差。sd(c)と近い値になることを確認

## [1] 6.454972

sum(c < mean(c) + 2 * sd(n)/sqrt(20) & c > mean(c) - 2 * sd(n)/sqrt(20)) #95%信頼区間に入っている標本平均の数

## [1] 9567

------------------------------------------------------

組み合わせ上、出現しない標本平均が存在するためか、今までの場合と異なり歯抜けのような分布となるが、やはり標本平均の分布は正規分布に近いものとなる。そして、標準誤差と95%信頼区間も同様の性質を示すことが確認できた。

中心極限定理が利用できない状況

　上記で示したように、母集団がほとんどどんな分布であっても中心極限定理は活用でき、平均値の推定を行うことができる。もちろん、最後に紹介した極端な離散分布にように、その平均値にどれほど意味があるかは、推定のうえで常に考えておく必要はある。注意しなければいけないのは、中心極限定理が有効に働かない場面があることだ。例えば、母集団がコーシー分布のような平均や分散が定義できない分布であるときは、標本平均の分布が正規分布に収束するわけではない。普通、標本は有限なサンプルサイズであり、必ず標本平均や不偏分散を計算できるので、形式上は中心極限定理を適用できる。しかし、外れ値が生じやすい標本の場合は、母集団がコーシー分布のようなものの可能性があるので要注意だ。例えば、xとyがともに正規分布に従うような変数の場合、x/yのような割り算値はコーシー分布に従う可能性があることが知られている（分母が0付近に近づくことで極端に大きい値が生じる）。

標準偏差か標準誤差か？

　話は変わるが、標準誤差は中心極限定理と密接にかかわる概念なので、ここで標準誤差に関する実用上の注意点を紹介する。標準誤差は、標準偏差とともにグラフのエラーバーとして表記することがある。

この時に標準偏差と標準誤差をはっきり区別せず、ともにデータのばらつきを示す指標くらいの気持ちで使ってしまうことがある。しかし、ここまで読んできた方なら、これらの2つが全く違うものを示していることはお分かりだろう。標準偏差は分散の平方根なので、標本のばらつきを示す指標だ。一方で、標準誤差は標本の平均値の信頼区間、精度を示している。標本そのものに興味があるか、平均値に興味があるかで使い分けよう。決して、標準誤差のほうがエラーバーが小さくて見栄えがする、などという理由で標準誤差を選んではいけない。

世の中に正規分布が多い理由の一つ

　本稿始めのほうで、中心極限定理を「母集団の平均がμ、標準偏差σのとき、この母集団から生じたサンプルサイズnの標本の平均値は平均がμ、標準偏差σ/√nの正規分布に従う」と紹介した。同値の内容として中心極限定理を「母集団の平均がμ、標準偏差σのとき、この母集団から生じたサンプルサイズnの標本の和は平均がnμ、標準偏差(√n)σの正規分布に従う」という、標本の平均ではなく、標本の和に注目して定義されることがある。世の中の現象は、多くの確率変数の効果の和として生じている。中心極限定理は、このような世の中の現象が正規分布に従って生じうることの一つの説明ともなっている。例えば、人間の身長はたくさんの遺伝子の効果と環境の効果の和で決まる。遺伝子のひとつひとつは、身長を高くする遺伝子かそうでないかの0/1の効果だし、環境の効果は温度のように連続的な効果かもしれないし、ある温度で閾値がある0/1の効果かもしれない。だが、これらが全部足されて生じた結果が身長なら、身長は正規分布に従う可能性があるのだ。事実、身長の分布はきれいな正規分布となることが知られている。

補足1：信頼区間の解釈

　信頼区間の説明で、本稿では簡単のため「95%の確率で母平均が収まっている範囲」と紹介した。頻度主義統計では、この表現は間違いだということを補足しておく。母平均は真のただ一つの値が存在している。そうであれば、信頼区間が決まった時に、その範囲には「母平均が存在するか、存在しないか、の2択」であるので、「95%の確率で存在する」という表現はおかしいはずだ。例えば、平均20、標準偏差5の母集団からサンプルサイズ30のデータが得られたと、わかっているとしよう。この時、95%信頼区間の範囲は± 1.96 × 0.912くらいだ。もし、標本平均が21だったら、19.21～22.79の範囲となり、母平均は含まれている。もし、標本平均が17だったら、15.21～18.79となり、母平均は含まれていない。このように、母平均は含まれるか含まれないかの2択だ。では、「95%」というのは何に対する確率なのかというと、「100回信頼区間を求めたら、95個くらいの信頼区間は母平均を含んでいる」ということなのだ。これが信頼区間の正しい解釈だ。コードを書いて確認してみよう。

------------------------------------------------------

n <- rnorm(10000, mean = 20, sd = 5) #母平均20、標準偏差5の正規分布の確認

histn(n, xlab = "value", ylab = "頻度")

mean(n)

## [1] 19.96774

sd(n)

## [1] 5.029822

------------------------------------------------------

------------------------------------------------------

ls <- NULL

for (i in 1:10000){ #信頼区間の判定を10000回行う

  s <- rnorm(30, mean = 20, sd = 5) #母平均20、標準偏差5の正規分布から30個のサンプルを生成

  m <- mean(s) #標本平均の計算

  d <- 5/sqrt(30) #標準誤差の計算

  upper <- m + 1.96 * d #上側95%信頼区間の計算

  lower <- m - 1.96 * d #下側95%信頼区間の計算

  e <- sign(upper - 20) * sign(lower - 20) == -1  #95%信頼区間に母平均が含まれるかの判定

#sign関数は符号判定の関数で、信頼区間から母平均を引いたときに信頼区間が母平均を含んでいれば、上側信頼区間と下側信頼区間の符号の積は負になる。

  ls <- rbind(ls, c(lower, m, upper, e))

}

ls <- as.data.frame(ls)

colnames(ls) <- c("Lower", "Mean", "Upper", "Exist")

sum(ls$Exist) #10000個中、95%信頼区間に母平均が含まれる個数

## [1] 9477

plotn(c(0,0), xlim = range(ls[,1:3]), ylim = c(1,100), xlab = "信頼区間", ylab = "試行回数")

for(i in 1:100){

  cols <- if(ls$Exist[i] == 1){

    "black"

  } else {

    "red"

  }

  overdraw(segments(ls$Lower[i],i,ls$Upper[i],i, col = cols))

}

overdraw(abline(v = 20, lty = 2))

------------------------------------------------------

確かに95%信頼区間を10000回求めると、約95%の信頼区間が母平均を含んでいた。さらにそれを図示したものが図10である。この図では赤で示された100個中3個の信頼区間が母平均を含んでいない。信頼区間の正確な解釈は上記の通りだが、通常、私たちは手元の標本のからただ一つの標本平均と標準誤差しか得ることはできない。そんなとき、95%信頼区間を求めたとして、「実際は母平均を含んでるかどうかわかりません」と解釈されても困ってしまうだろう。そこで、正確な解釈ではないと承知しつつも「95%くらいの確率で信頼区間には母平均は含まれてます」と解釈してしまう方が、まだ素直に受け取れよう。

補足2：母集団が正規分布の時の中心極限定理の証明

　中心極限定理は、母集団の分布に限らず成立する性質である。しかし、あらゆる分布について証明するのは骨が折れる（というか私にはできない）ので、母集団が正規分布の時に限って成立することを示すことにする。証明には、正規分布の再生性を利用する。

　下図のように、平均μ、分散σ^2の正規分布からサンプルサイズnでサンプリングすることを考える。このとき (データX1～Xnの和) = ΣXi が従う確率分布は、正規分布の和の再生性より、平均nμ、分散nσ^2である。

確率変数の和ΣXi をサンプルサイズnで除したものは、すなわち標本平均である。確率変数の和ΣXi が平均nμ、分散nσ^2の正規分布に従うのであれば、定数倍の再生性より、標本平均(1/n)ΣXi が従う正規分布は平均(1/n) nμ = μ、分散(1/n)^2 nσ^2 = σ^2/nの正規分布である。

基本的には上記の発想の下、母集団がポアソン分布、二項分布、ガンマ分布であっても、(1)これらの確率分布はパラメータによって正規分布への近似が可能になること、(2)それらの再生性を利用してゆくことで、各個の確率分布について中心極限定理が証明できる。