はじめに
クロス集計表の検定にχ二乗検定というのがあって、それをやるとカテゴリー同士の関連（たとえば男はお茶を好むが女は水を好むなど）を統計的に検定することができる。χ二乗検定そのものについて詳しく知りたい人は上のリンク参照。
ただχ二乗検定を行うにあたっては２つほど注意点があり、たとえp値が有意水準を下回ったとしても
- カテゴリー別に見てもサンプル数が十分に確保されているか
- 実際にどのカテゴリーとどのカテゴリーが関連しているのか
を検討する必要がある。結局χ二乗検定で検定しているのはクロス集計表全体の偏りの有無であり、検定結果が有意になったとしてもすべてのカテゴリーの組み合わせに関連がみられたということではないのだ。
コクランの規則
本来はχ二乗検定をやる前に確かめておく必要がある基準である。χ二乗検定の適用基準としてクロス集計表の期待値が５未満のセルが全体の２０％以上になってはいけないというのがあるが、これがコクランの規則にあたる。「２０％以上」の部分には諸説あり「２５％以上」「２０％より大きく」など色々な表記をみかけます。
Pythonでχ二乗検定を行う場合はscipy.stat.chi2_contingencyを使うことが多いと思うので、この関数が返す期待値表を使ってコクランの規則を確かめる。
# χ二乗検定 crossはnumpy２次元配列です
x2, p, dof, expected = stats.chi2_contingency(cross)
expected = np.array(expected)
# コクランの規則
expected < 5
これでTrueのセルが全体の２０％未満であればコクランの規則を満たしている。
もしもデータがこの規則を満たさないことが分かった場合はフィッシャーの正確確率検定に移行するのがいいだろう。
残差分析
Rでχ二乗検定をすれば検定と同時に各セルの調整済み標準化残差を返してくれるので問題ないのだが、Pythonの場合は手計算でやる必要がありそうだ。
残差の定義は
残差 = 観測値 - 期待値
であるが、調整済み標準化残差を算出するために、新たに残差分散というものを定義する必要がある。
残差分散 = (1 - \frac{横周辺和}{総数})(1 - \frac{縦周辺和}{総数})
具体的にどういうものかっていうのは、下の参考サイトを参照してください。
とにかくこれをもとに、調整済み標準化残差が
調整済み標準化残差 = \frac{残差}{\sqrt{期待値*残差分散}}
のように計算できます。
ここまでの流れをPythonコードで書くとこんなかんじ。
# 残差
res = cross - expected
# 残差分散を求める
res_var = np.zeros(res.shape)
it = np.nditer(cross, flags=['multi_index'])
while not it.finished:
var = (1 - (cross[:,it.multi_index[1]].sum() / cross.sum()))*(1-(cross[it.multi_index[0],:].sum() / cross.sum()))
res_var[it.multi_index[0], it.multi_index[1]] = var
it.iternext()
# 調整済み標準化残差を求める
stdres = res / np.sqrt(expected * res_var)
# この調整済み標準化残差が絶対値1.96以上になれば有意な差を主張できます。ここでは正規分布表から値をp値に直して表示しています。
np.apply_along_axis(stats.norm.sf, 0, np.abs(stdres[0,:]))
お役に立てましたでしょうか？