カイ二乗検定に関するメモ

カイ二乗検定って何でああいう形の式
$\frac{(Y_1 - n p_1)^2}{n p_1} + \cdots + \frac{(Y_K - n p_K)^2}{n p_K}$
を使うのか?という疑問について。

検定について

カイ二乗検定自体の前に、検定一般について必要そうなことを書いておく(そもそも統計一般がすごくわかりにくい気がしてしょうがない。「帰無仮説」とか「有意」みたいな用語も、わざとわかりにくくするために使っているんじゃないかとすら思う)。
統計的検定は、科学における反証手続きに何となく似ているので、それと比較してみる。

反証

AならばB。
でもBでは無かった。
したがってAではない。

というのが反証の基本構造。例をあげると次のような感じ。

ニュートン力学が正しいならば、○○の時刻に××の位置に惑星が観測される。
でも違う位置で観測された。
したがってニュートン力学は間違っている。

もちろん実際には前提Aにあたる部分には「ニュートン力学が正しい」以外の様々な仮定が置かれているので、Bで無かったとしても前提のどの部分を否定するのかは一概には決められない。またBだった場合、論理的にはそこからAが正しいことは導けない。でもBが非常に正確な予言だったり意外な結果だったりする場合、BであることはAの正しさを支持するように思いたくなる。でも反証主義的には、Bだったという結果はAが正しいことを支持しない(けれど、信頼性は上がる)。
この「前提の何を否定するか」というのと「否定ではない結果だった場合をどう考えるか」という問題は、統計的検定でも同じように現れる問題。

統計的検定

統計的検定の構造は、一見すると反証に似ている。

Aならば、Bである確率が非常に高い。
Bでは無かった。なので「Bである確率が非常に高い」は間違いに違いない(間違いだと見なそう)。
したがってAではない。

あるいは

Aならば、Cである確率は非常に低い。
Cだった。なので「Cである確率は非常に低い」は間違いに違いない(間違いだと見なそう)。
したがってAではない。

というのが、統計的検定(における「棄却」)の基本構造。
反証の構造に似ているけど、反証が論理的な推論なのに対して棄却はそうではないという大きな違いがある。なぜなら「Bでない」「Cである」ことと「Bである確率が非常に高い」「Cである確率は非常に低い」ことに論理的な推論関係がないので。
むしろ、あったらおかしな事になる。宝クジを例に取ってみる。統計的検定を使って次のように考えたとする。

宝クジの当選番号決定が公正ならば、××××番が1等になる確率は非常に低い。
××××番が1等だった。だから「××××番が1等になる確率は非常に低い」は間違いに違いない。
だから、宝クジの当選番号決定は公正でない。

「××××番」をどんな番号にしても「××××番が1等になる確率は非常に低い」が成り立つ。だから、もしも棄却が論理的推論なら、どの番号が1等になっても「宝クジの当選番号決定が公正でない」という結果になってしまう。
何を確率の低いことと見なすかは、実践的な方法論とか態度によって決めるものであって、誰かが正しさを保証してくれるようなものではない。また、確率が低ければそちらの方が良いというわけでもない。「この100面体サイコロは公正である(各面の出る確率はそれぞれ1/100である)」というのを考えてみる。

1の目が出る確率は非常に低い(確率1パーセント)。
100回振った平均が50±7に入らない確率は非常に低い(確率2パーセント以下)。

確率だけなら1の方が低いけど、やり方としては1よりは2を使うのが良い気がする。でも統計的検定という手法そのものには、どちらを選んだ方が良いかに対する処方箋はない(と思う)。

本題: カイ二乗検定

それで問題は「何を確率の低いことと見なすか」ということ。
カイ二乗検定はいくつかの場合に使われるけど、一番基本的なのは次のような場合。

あり得る結果はK種類(それぞれ事象1、……、事象Kと呼ぶことにする)。それぞれが起こる確率は $p_1, p_2 , \ldots, p_K$ である。

この確率分布 $p_1, p_2 , \cdots, p_K$ の正しさを検定したいときに使う。上で例にあげた100面体サイコロの公正さの検定もこれに含まれる。
このような場合、確率空間がきっちり定まっているから、どんな結果の確率でも計算できる。だから検定のやり方として、結果の平均を調べてその値が取らなそうな範囲を決めて検定するとか、結果の分散を調べてその値が……とか色々なやり方が考えられる。にもかかわらず普通はカイ二乗検定が使われる。なぜ他のやり方ではなくカイ二乗検定を使うのか、というのがここでの問題。

カイ二乗検定ができるまで

まず平均や分散(標本の平均と標本の分散)を使いたくない理由はわかる。平均や分散を求めるためには結果が数でないといけない。でも結果である事象1から事象Kに対して番号や数値がふってあるとは限らないし、ふってあってもその番号や数値は、確率分布 $p_1, p_2 , \cdots, p_K$ の正しさには関係がない。だから、そのようなものには依存しないものによって検定をしたい。
もうひとつ重要なのは計算が(比較的)簡単にできるということ。直接計算しようとすると多項係数が出てきて計算が面倒そう。計算を楽にするということで、まず思いつくのは中心極限定理が使えないかということ。

中心極限定理

中心極限定理というのは、確率論で中心に位置する極限定理ってことで、ポーヤ・ジェルジ(ジョージ・ポリア)が命名。
例えばサイコロを何度も振ると、

6、2、1、4、5、3、3、2……

みたいな結果が得られる。次にサイコロを100回振るごとの平均を取ると

3.46、3.70、3.32、3.49、3.56、3.30、3.35、3.37……

のような結果が得られ、10000回ごとの平均を取ると

3.4932、3.5286、3.4888、3.4947、3.4816、3.4929、3.4881、3.5315……

のような結果が得られる。これを見るとnが大きい時ほど、サイコロの期待値3.5に近い値が多く出ていることがわかる。
この例に限らず一般的に、nが大きくなると「n回の平均と確率の期待値との差」が小さい確率は高くなる(正確な表現は省略)。この性質を大数の(弱)法則と呼ぶ。
でも、平均と期待値の差が小さいといっても、完全に一致するわけではなく誤差が生じる。その誤差がどんな分布をしているか見るために、平均と期待値との差を√n倍したものの確率分布を調べてみる(√n倍するのは、そうすると分散が一定になるから)。するとnを大きくしていったとき、元の分布がどうだったかに関係なく、平均と期待値との誤差(を√n倍したもの)の分布は正規分布に近づいていく。このことを中心極限定理と呼ぶ。
そうなる理由は特性関数を見るとわかる。
まず、確率変数Xの特性関数 $\phi(t)$ というのは、 $e^{\sqrt{-1}t X}$ の期待値 $E[e^{\sqrt{-1}t X}]$ のこと。特性関数については次のような性質がある。

独立な確率変数 $X_1$ と $X_2$ の特性関数がそれぞれ $\phi_1(t)$ と $\phi_2(t)$ のとき、 $X_1 + X_2$ の特性関数は $\phi_1(t) \phi_2(t)$ となる。
期待値が $\mu$ で標準偏差が $\sigma$ の確率変数の特性関数をテーラー展開すると $\phi(t) = 1 + \sqrt{-1}\mu t - \frac{1}{2}(\mu^2+\sigma^2) t^2 + \cdots$ となる。
期待値が $\mu$ で標準偏差が $\sigma$ の正規分布の特性関数は $e^{\sqrt{-1}\mu - \frac{1}{2}\sigma^2 t^2}$ である。

同じ事柄をn回繰りかえした結果を、独立な確率変数 $X_1,\ldots,X_n$ で表しておく(それぞれ期待値は $\mu$ で標準偏差は $\sigma$ とする)。n回の平均と期待値との差を√n倍したものは
$(\frac{X_1 + \cdots + X_n}{n} - \mu) \times \sqrt{n} = \frac{X_1 - \mu}{\sqrt{n}} + \ldots + \frac{X_n - \mu}{\sqrt{n}}$
と書き換えられる。 $\frac{X_1 - \mu}{\sqrt{n}}, \cdots , \frac{X_n - \mu}{\sqrt{n}}$ は互いに独立で、期待値はそれぞれ0で、標準偏差はそれぞれ $\frac{\sigma}{\sqrt{n}}$ なので、 $(\frac{X_1 + \cdots + X_n}{n} - \mu) \times \sqrt{n}$ の特性関数は $( 1 - \frac{\sigma^2}{2n} t^2 + \cdots)^n$ となる。説明は略すけど、nを大きくしたとき3次以上の項(式の「 $\cdots$ 」の部分)の影響は小さくなっていく(つまり、元の分布に応じて特性関数の3次以上の項は違っているのだけど、その影響はnが大きくなるにつれて小さくなっていく)。
そしてnを大きくしたとき $( 1 - \frac{\sigma^2}{2n} t^2 + \cdots)^n \to e^{- \frac{1}{2}\sigma^2 t^2}$ というように正規分布の特性関数に収束していく。
このように、nを大きくしていったときに、元の特性関数のうちの2次の項以外の影響が消えていくということが、元の確率変数の分布に関係なく中心極限定理が成り立つ理由。

中心極限定理を使ってみる

もとの問題に話を戻す。
とりあえず結果が事象1だったとき1で事象1以外だったとき0となる確率変数 $X_1$ を考える。n回繰りかえしたときに事象1が起こった回数は $Y_1 = X_1^{(1)} + \cdots + X_1^{(n)}$ と書ける。中心極限定理により、nを大きくすると、 $\frac{Y_1 - n p_1}{\sqrt{n}}$ の分布は正規分布に近づいていく。これで確率計算が簡単になる。でもこの値だけでは事象1しか見ていないので、他の事象についても、それが起こった回数 $Y_2,\cdots,Y_K$ と期待値との誤差を考える。
この場合、それぞれだけ考えた場合は「各誤差はそれぞれ正規分布になる」で済むのだけど、互いに独立ではないので、全体の確率分布についてはもう一度調べ直す必要がある(中心極限定理の場合と同じように特性関数を計算していく)。
計算してみると、
$(Z_1,\ldots,Z_K) = (\frac{Y_1 - n p_1}{\sqrt{n p_1}},\ldots, \frac{Y_K - n p_K}{\sqrt{n p_K}})$
の特性関数が
$\phi(t_1, \ldots, t_K) = \exp ( - \frac{1}{2} \sum_i t_i^2 + \frac{1}{2} \sum_{i,j} \sqrt{p_i p_j} t_i t_j + O(\frac{1}{\sqrt{n}}) )$
となり、「 $O(\frac{1}{\sqrt{n}})$ 」の部分を無視すると、多次元正規分布の特性関数になっている。
でもこれでもまだ $(Z_1,\ldots,Z_K)$ についての確率計算は面倒。コンピュータを使えば計算自体はそれほど大変ではない気もするけど、数表を用意するのは無理だし、どこを「確率の低い部分」に選ぶのかという問題もある。

変数変換する

$(Z_1,\ldots,Z_K)$ は互いに独立ではないので、互いに独立な新しい確率変数に変数変換できないかを考える。そうすると、説明は省略するけど、適当な直交行列によって、
$\left( \begin{array}{c}W_1 \\ \vdots \\ W_{K-1} \\ 0 \end{array} \right) = \left( \begin{array}{ccc}U_{11} & \cdots & U_{1K} \\ \vdots & \ddots & \vdots \\ U_{(K-1)1} & \cdots & U_{(K-1)K} \\ \sqrt{p_1} & \cdots & \sqrt{p_K} \end{array}\right) \left( \begin{array}{c}Z_1 \\ \vdots \\ Z_{K-1} \\ Z_K \end{array} \right)$
と変換できることがわかる。ただし $(W_1,\ldots,W_{K-1})$ は互いに独立で、それぞれが期待値0、標準偏差1の正規分布。
今度は、適当な直交行列を見つけて計算するという面倒な話が出てきた。

カイ二乗検定量

しかし直交行列を実際に計算しなくて済むやり方がある。
直交行列はノルムを保存するという性質により $W_1^2 + \cdots + W_{K-1}^2 = Z_1^2 + \cdots + Z_K^2$ となっている。左辺 $W_1^2 + \cdots + W_{K-1}^2$ は、独立な正規分布を二乗して足したものなので確率分布の計算は難しくない(具体的には自由度K-1のカイ二乗分布になる)。一方右辺は $Z_1^2 + \cdots + Z_{K}^2 = \frac{(Y_1 - n p_1)^2}{n p_1} + \cdots + \frac{(Y_K - n p_K)^2}{n p_K}$ なので、それぞれの事象が起こった回数から簡単に計算できる。
結局 $\frac{(Y_1 - n p_1)^2}{n p_1} + \cdots + \frac{(Y_K - n p_K)^2}{n p_K}$ の値(カイ二乗検定量と呼ばれる)について検定をすれば、比較的楽に計算がおこなえることがわかった。

カイ二乗検定を使う理由

これまでの数式の展開を見ると、カイ二乗検定の式と同じくらい計算が簡単で検定に使えそうな値は他には無い(少なくとも簡単には見つけられない)ことがわかる。これが、他でなくカイ二乗検定が使われる理由だと思う。

どの値を「めったに起こらないこと」に選ぶか。

$\frac{(Y_1 - n p_1)^2}{n p_1} + \cdots + \frac{(Y_K - n p_K)^2}{n p_K}$ を使うことに決めても、この値がとり得る範囲のうちどの領域をめったに起こらないこと(めったに取らない値)と見なすかはまだ決まっていない。
カイ二乗検定量は期待値とのズレの大きさのようなものなので、それが大きすぎる領域は当然選ばれるだろう。
問題はズレが小さい(0に近い)領域をどうするか。ズレが小さいのだから問題ないと考えることもできるし、ズレが小さすぎるのもおかしい(多少ズレがある方が自然、ズレがないのは確率モデルの立てかたに問題があるとか作為の介入が原因かもしれない)と考えることもできる。どちらにするかを決める標準的な方法はたぶんなくて、使用者が場合ごとに判断しかないのだろう。何かもっともらしい説明や処方箋があると便利そうだけど、それともやりたいようにやれるのが便利なのか。