カイ2乗検定

ポアソン分布からカイ2乗分布へ

個数 n1, n2, n3, が独立で,それぞれ平均 μ1, μ2, μ3, ポアソン分布をするとします。ポアソン分布ですから ni の分散はその平均 μi に等しく,ni が10程度以上であれば正規分布で近似できるので,個数から平均を引いて標準偏差で割ればほぼ標準正規分布 N(0,1) に従います。すると,

χ2=i=1N(niμi)2μi

は独立な標準正規分布の2乗の N 個の和です。一般に,独立な標準正規分布の2乗の N 個の和の確率分布を,自由度 N のカイ2乗分布と呼びます。

具体的な問題を考えましょう。毎年の事故発生件数が n1, n2, n3, , nN であるとき,1年あたりの平均事故件数は

m=1Ni=1Nni

になります。したがって,これらが独立にポアソン分布に従うという帰無仮説のもとで

χ2=i=1N(nim)2m

は自由度 N のカイ2乗分布に従いそうですが,ここでは m は実際のデータから計算していますので,正しくは自由度 N1 のカイ2乗分布になります。

個数の列が与えられたとき,上の χ2 の値を計算し,それがカイ2乗分布で予想される範囲から外れていれば,独立なポアソン分布であるという帰無仮説が疑わしくなります。このことを使った検定法は,ピアソン(Karl Pearson)のカイ2乗検定(χ2 検定)と呼ばれます。

例えばPISA「盗難事件」問題ほかで挙げた6年間の事故発生件数

> x = c(762,792,795,794,849,834)

が独立なポアソン分布に従うという帰無仮説についてカイ2乗検定をしてみましょう。Rでは次のようにして簡単に計算できます。

> chisq.test(x)

	Chi-squared test for given probabilities

data:  x 
X-squared = 6.2329, df = 5, p-value = 0.2842

p 値は 0.2842 ですので,帰無仮説と矛盾しません。

分割表のカイ2乗検定

2×2分割表の縦計・横計

abecdfghn

が与えられたとき,もし表の各行(各列)が独立ならば,

egnehnefgnfhnfghn

になるはずです。この予測値(期待度数)と実際の値(観測度数)との差の2乗を予測値で割ったものの合計

(aeg/n)2eg/n+(beh/n)2eh/n+(cfg/n)2fg/n+(dfh/n)2fh/n

は近似的に自由度1の χ2(カイ2乗)分布に従います。自由度が1であることは,縦横の合計を固定した 2×2 分割表の自由度が一つしかないことを考えれば納得できるでしょう。n1×n2 分割表では自由度は (n11)(n21) になります。

2×2 分割表の場合に限って,それぞれの「予測値と実際の値との差」を0に近い側に0.5ずつ修正するというイェーツの連続性の補正(Yates' continuity correction)をする場合があります。

Rで χ2 検定をする関数は chisq.test() です。単に

> chisq.test(matrix(c(228,863,284,851),nrow=2))

などとすれば,2×2 分割表の場合に連続性の補正を行います。結果は p=0.02377 です。連続性の補正を行いたくないなら,

> chisq.test(matrix(c(228,863,284,851),nrow=2), correct=FALSE)

のようにします。結果は p=0.02082 となります。ちなみに Fisherの正確検定 の結果は p=0.02335 です。

「カイ2乗…」と書いても「カイじじょう…」と読む人が多いのは,昔「2乗」のことを「自乗」(じじょう)と言ったことの名残りだろうと思います(私は「にじょう」と読みます)。

p 値は,連続性の補正をすれば大きめになり,しなければ小さめになります。 Fisherの正確検定かカイ2乗検定か も参考にしてください。

[この項要再検討] 観測度数を Oi,期待度数を Ei とすると,χ2 検定は

χ2=i(OiEi)2/Ei

を使いますが,類似のものとして,G検定(G-test)があります。これは

G=2iOilog(Oi/Ei)

の分布が同じ χ2 分布で近似できることを使います(Oi=0 の項は対数が計算できませんが,単に 0 として計算します)。分野によっては χ2 検定より人気があるようです。図は各セルに入る確率が [0.120.280.180.42] で与えられる総度数50の 2×2 分割表を1000個作って,横軸はFisherの方法の p 値,縦軸は連続性の補正なしの χ2 検定の p 値(○)とG検定の p 値(×)をプロットし,対応する○と×を線分で結んだものです。

少なくともこの場合にはG検定の p 値は小さすぎることがわかります。

このシミュレーションに使ったプログラムは次の通りです:

prob = c(0.12,0.18,0.28,0.42)
ni = 1000
pfis = pchi = pg = numeric(ni)
for (i in 1:ni) {
  k = sample(1:4, 50, replace=TRUE, prob=prob)
  na = sum(k == 1)
  nb = sum(k == 2)
  nc = sum(k == 3)
  nd = sum(k == 4)
  a = matrix(c(na, nb, nc, nd), nrow=2)
  pfis[i] = fisher.test(a)$p.value
  pchi[i] = chisq.test(a, correct=FALSE)$p.value
  pg[i]   = g.test(a)$p.value
}
plot(pfis, pchi, xlim=c(0,0.1), ylim=c(0,0.1), xlab="", ylab="")
abline(0,1)
for (i in 1:ni) {
  if (!is.na(pchi[i])) {
    points(pfis[i], pg[i], pch=4)
    lines(c(pfis[i],pfis[i]), c(pchi[i],pg[i]))
  }
}

度数分布データと仮説との整合性を調べる χ2 検定では,χ2(0) の値が大きいほど仮説からの外れが大きいので,通常は大きい側だけを調べればいいのですが,逆に 0 に近すぎる(仮説と合いすぎる)場合も,データ改竄が疑われるという意味で,興味があります。例えば有名な遺伝の法則を発見したメンデルのデータは,理論値 3:1 に合いすぎているという指摘がされています(R. A. Fisher, ``Has Mendel's work been rediscovered?'', Annals of Science 1, 115-137 (1936); C. R. ラオ『統計学とは何か―偶然を生かす』(藤越康祝ほか訳,丸善,1993年))。


Last modified: