どうでもよい覚え書き:McNemar検定と「対応のあるZ検定」

 これは仕事のための単なる覚え書きで、面白くも何ともない話であります。

 ごくごく教科書的な標準的知識でも、分野が違うとその分野に特化した表現となり、ちょっと面食らうことがある。そんなことありませんか? ありますよね?

 個人的な思い出話だけど、市場調査の会社に拾って頂いて、おまえ統計学できるんだろ、というのでいろいろ訊かれるようになったんだけど(おまえ日本人だからカラテできるだろうというの似ている)、ごく平凡な手法について訊ねられたのに、見たこともない式を示されて、アレッ、と面食らったことがあった。

 私は心理学の勉強をしてたんですけど、心理学専攻だったらきっとMcNemar検定っていうのを教えますよね。2×2のクロス表があり、セルの頻度が\(n_{11}, n_{12}, n_{21}, n_{22}\)だとして、行の割合と列の割合が違うかどうか調べたいとき、検定統計量は\(z = (n_{21} -n_{12}) / \sqrt{n_{21} + n_{12}} \)だ、っていうの。やりますよね。やりませんでしたっけ。
 いっぽうですね、市場調査の教科書をみると、二回の測定において観察された(たとえば)反応率をそれぞれ\(p_1, p_2\)、二回とも反応した人の割合を\(p_{12}\)として、検定統計量は\(z = \frac{p_1-p_2}{\sqrt{(p_1+p_2+2p_{12})/N}}\)、なんて書いてある。McNemar検定という呼び方はせず、たとえば「同じ集団内の比率の検定」とか、対応のあるZ検定などと呼んでいる。安易に外人の人名をつけられるのも困るけど、こういう曖昧な表現もそれはそれで困る。

 というわけで、平凡な検定手法のはずなのに見慣れぬ検定統計量が出てきて面食らい、しばし紙と鉛筆で考えこんだ覚えがある。頭の悪い哀しさで、この2つの検定統計量は等しいと納得するまでにしばらくかかった。

 それから幾星霜、当時からくたびれていた私はさらにくたびれたが、この話は数年に一回思い出して、いや待て、なんで同じなんだっけ… ああそうか… というのを繰り返している。このたびも元同僚の友人から急に電話がかかってきて、いや待て、なんでこうなるんだっけ… ああそうか… というのを繰り返した。
 というわけで、以下、覚え書きである。

 準備その1。
 変数\(X, Y\)の両方に変数\(Z\)を足したとき、その共分散はどうなるか。$$ cov[X+Z, Y+Z] $$ $$= E[((X-\mu_X)+(Z-\mu_Z))((Y-\mu_Y)+(Z-\mu_Z))]$$ $$=E[(X-\mu_X)(Y-\mu_Y)] + E[(X+Y-\mu_X-\mu_Y)(Z-\mu_Z)] + E[(Z-\mu_Z)^2]$$ $$ = Cov[X, Y] + Cov[X+Y, Z] + Var[Z]$$ つまり、\(Z\)の分散が上積みされるだけでなく、\(X+Y\)と\(Z\)の共分散も上積みされるわけね。

 準備その2。
 値\((0,1)\)をとる二値変数\(Y_{j}\)を考え(\(j\)は個体を表すインデクス)、\(E[Y_j] = \pi\)とする。分散は$$ var[Y_j] = E[ (Y_j – \pi)^2 ] = E[Y_j^2] – 2\pi E[Y_j] – \pi^2 = \pi – \pi^2 = \pi(1-\pi)$$ である。
 \(N\)個の個体のうち\(Y_j = 1\)であった個体の数を\(n\)とすれば、\(n=\sum_{j=1}^N Y_j\)だから、期待値は\(E[n] = N\pi\)、分散は\(var[n] = N\pi(1-\pi)\)である。
 \(N\)個の個体のうち\(Y_j = 1\)であった個体の割合を\(p = n/N\)とすれば、期待値は\(E[p] = \pi\)、分散は\(var[p] = \pi(1-\pi)/N\)である。
 ま、そりゃそうかな、って話ですね。

 準備その3。
 今度は、値\((0,1)\)をとるふたつの二値変数\(Y_{1j}, Y_{2j}\)を考え、\(E[Y_{1j}] = \pi_1, E[Y_{2j}] = \pi_2\)とする。共分散は$$ cov[Y_{1j}, Y_{2j}] = E[ (Y_{1j} – \pi_1)(Y_{2j} – \pi_2) ] = E[Y_{1j}Y_{2j}] – \pi_1 \pi_2$$である。
 さて、実はこの\(Y_{1j}, Y_{2j}\)は、\(c\)個のカテゴリを持つ多項変数のうち2つのカテゴリに注目し、値が1であることを\(Y_{1j}=1\), 2であることを\(Y_{2j}=1\)と表している変数だった、としよう。その場合、同時に1にはなれないから、\(Y_{1j}Y_{2j} = 0\), よって\( cov[Y_{1j}, Y_{2j}] = -\pi_1 \pi_2\)である。
 \(N\)個の個体のうち\(Y_{1j}=1\)であった個体の数を\(n_1\), \(Y_{1j}=2\)であった個体の割合を\(n_2\)とすれば、\( cov[n_1, n_2] = – N \pi_1 \pi_2 \)である。
 \(N\)個の個体のうち\(Y_{1j}=1\)であった個体の割合を\(p_1\), \(Y_{1j}=2\)であった個体の割合を\(p_2\)とすれば、\( cov[p_1, p_2] = – \pi_1 \pi_2 / N \)である。

 ここから本題。
 \(N\)個の個体を2回測定したとする。1回目の測定値を\(a\), 2回目の測定値を\(b\)とする。いずれもとりうる値は\(\{1, 2\}\)であるとする。
 1回目の測定で\(a\), 2回目の測定で\(b\)だった個体数を\(n_{ab}\)と書く。するとこういう2×2クロス表が得られる。$$ \begin{array}{c|cc|c} & b=1 & b=2 \\ \hline a=1 & n_{11} & n_{12} & n_{1+} \\ a=2 & n_{21} & n_{22} & n_{2+} \\ \hline & n_{+1} & n_{+2} & N \\ \end{array} $$ 割合を\(p_{ab} = n_{ab}/N\)と書く。データを多項分布からの標本と考え、ある人が\((a, b)\)となる確率を\(\pi_{ab}\)とする。
 いま関心があるのは、\(\pi_{a+}\)と\(\pi_{+a}\)の差である。

 割合の差\(\delta = \pi_{1+} – \pi_{+1} \)を推定することを考えよう。不偏推定量は\(d = p_{1+} – p_{+1}\)だ。では、\(d\)の分散は? $$ var[d] = var[p_{1+}]+var[p_{+1}] – 2cov[p_{+1}, p_{1+}]$$ 第1項、第2項はそれぞれ割合の分散だから、準備その2より$$ var[p_{1+}] = \frac{\pi_{1+} (1-\pi_{1+})}{N}$$ $$ var[p_{+1}] = \frac{\pi_{+1} (1-\pi_{+1})}{N}$$である。
 第3項の共分散は$$ cov[p_{+1}, p_{1+}] = cov[p_{11} + p_{21}, p_{11}+p_{11}] $$ 準備その1より $$ = cov[p_{21}, p_{11}] + cov[p_{21}+p_{12}, p_{11}] + var[p_{11}]$$ 準備その3より$$ = -\frac{\pi_{21}\pi_{12}}{N} -\frac{(\pi_{21}+\pi_{12}) \pi_{11}}{N} + \frac{\pi_{11}(1-\pi_{11})}{N}$$ $$ = -\frac{\pi_{21}\pi_{12}}{N} + \frac{\pi_{11}(1-\pi_{11}-\pi_{21}-\pi_{12})}{N}$$ $$ = \frac{\pi_{11}\pi_{22} -\pi_{21}\pi_{12}}{N} $$ というわけで、$$ var[d] = var[p_{1+}]+var[p_{+1}] – 2cov[p_{+1}, p_{1+}]$$ $$= \{ \pi_{1+} (1-\pi_{1+}) + \pi_{+1} (1-\pi_{+1}) – 2 (\pi_{11}\pi_{22} -\pi_{21}\pi_{12}) \} / N$$ である。その推定量は $$ \hat{var}[d] = \{ p_{1+} (1-p_{1+}) + p_{+1} (1-p_{+1}) – 2 (p_{11}p_{22} -p_{21}p_{12}) \} / N$$ である。

 やれやれ、推定量の分散の推定量が手に入った。
 標本が十分に大きければ\(d\)は正規分布に従うと考えられるので、\(\delta\)のWald信頼区間は $$ d\pm z_{\alpha/2} \sqrt{\hat{var}[d]}$$ である。\(z_{\alpha/2}\)ってのはなんか適当なアレね、1.96とか。

 ここから先は、\(\hat{var}[d]\)の中身をどうやって整理して書くか、というちがいである。

 まず、ひたすらセルの割合に書き換えて…$$ N \hat{var}[d] $$ $$ = p_{1+} (1-p_{1+}) + p_{+1} (1-p_{+1}) – 2 (p_{11}p_{22} -p_{21}p_{12})$$ $$ = (p_{11}+p_{12}) – (p_{11}+p_{12})^2 + (p_{11}+p_{21}) – (p_{11}+p_{21})^2 – 2 (p_{11}p_{22} -p_{21}p_{12})$$ $$ = (p_{11}+p_{12}) – (p_{11}^2 + p_{12}^2 +2p_{11}p_{12}) + (p_{11}+p_{21}) – (p_{11}^2+p_{21}^2+2p_{11}p_{21}) – 2 (p_{11}p_{22} -p_{21}p_{12}) $$ $$ = p_{11}(1 -p_{11} -2 p_{12} + 1 – p_{11} – 2p_{21} -2 p_{22}) +p_{12} – p_{12}^2 + p_{21} – p_{21}^2 + 2 p_{21}p_{12} $$ $$ = p_{11}(2 – 2p_{11} -2 p_{12} – 2p_{21} -2 p_{22}) +p_{12} + p_{21} -( p_{12}^2 + p_{21}^2 – 2 p_{21}p_{12}) $$ 第1項が消えて $$ = p_{12} + p_{21} -( p_{12} – p_{21})^2 $$ \(H_0: \delta = 0\)の下では第3項が消えて $$ \hat{var}_0[d] = \frac{p_{12}+p_{21}}{N} = \frac{n_{12}+n_{21}}{N^2}$$ 検定統計量は $$ z_0 = \frac{d}{\sqrt{\hat{var}_0[d]}} = \frac{(p_{12}-p_{21}) N}{\sqrt{n_{12}+n_{21}}} = \frac{n_{12}-n_{21}}{\sqrt{n_{12}+n_{21}}}$$ というわけで、ご存じMcNamar検定となる。

 上記ではひたすらセルの割合に書き換えたが、途中で周辺割合\(p_{1+}, p_{+1}\)を復活させるために\(p_{11}\) を二回足してやれば$$ N \hat{var}[d] = p_{12} + p_{21} – (p_{12} – p_{21})^2 = p_{1+} + p_{+1} – 2p_{11} – (p_{12} – p_{21})^2 $$ \(H_0: \delta = 0\)の下では第4項が消えて $$ \hat{var}_0[d] = \frac{p_{1+}+p_{+1}-2p_{11}}{N}$$ 検定統計量は $$ z_0 = \frac{d}{\sqrt{\hat{var}_0[d]}} = \frac{p_{12}-p_{21}}{\sqrt{\frac{p_{1+}+p_{+1}-2p_{11}}{N}}}$$ というわけで、市場調査の教科書に書いてある式になる。

 いうまでもなく、どっちの書き方が良いという話ではない。
 いまとなっては、市場調査の教科書で後者の表現がされている理由も、なんとなくわかる気がする。調査実務家からいえば、対象者に二値変数を2回反復して測定することは多いが(現行製品と新製品への購入意向とか)、二元クロス表をみるのは稀であり、むしろふたつの反応率が載っている集計表のほうがなじみ深いだろう。