Rによるボンフェローニ補正
検定を複数回繰り返すことには問題があるため (統計検定の多重性の問題),様々な回避法が考案されている.チューキー・クレーマー検定,ダネット検定,スティール・ドゥワス検定,スティール検定等のアドホックな検定法が解決したい問題に適している場合はそれらを使うのが良いが,それ以外の場合においてはより簡易的な補正法を選択することもできる.最も有名な方法はボンフェローニ補正であるが,R ではこれを含めたいくつかの補正法が代表値の差の検定に関してはデフォルトで用意されている.平均値の差の検定を行いたい場合は,pairwise.t.test() が,中央値の検定を行いたい場合は,pairwise.wilcox.test() が使用できる.どちらの方法も,対応がない場合,対応がある場合の検定を実行することができる.
対応がない場合のボンフェローニ補正をした中央値の検定
以下のサンプルサイズがそれぞれ,9,10,8,8からなる互いに対応がないデータAからDが得られたときの各データ間における中央値の値に差があるかどうかをボンフェローニ補正をしたウィルコクソンの順位和検定で検定する.帰無仮説 (H0) は各群の中央値はそれぞれの群間で等しいこととなる.
データA | 301, 311, 320, 291, 388, 412, 325, 361, 287 |
データB | 197, 180, 247, 260, 248, 199, 179, 134, 163, 200 |
データC | 209, 302, 187, 166, 234, 290, 175, 116 |
データD | 342, 216, 316, 386, 324, 145, 254, 228 |
最初に,以下のコマンドにて,上のデータAからDをまとめて変数 vx に格納する.
$ vx=c(301, 311, 320, 291, 388, 412, 325, 361, 287, 197, 180, 247, 260, 248, 199, 179, 134, 163, 200, 209, 302, 187, 166, 234, 290, 175, 116, 342, 216, 316, 386, 324, 145, 254, 228)
次に,以下のコマンドにて,読み込んだ各データのグループ名を指定するための変数 vg を作成する.データAからDの標本サイズはそれぞれ,9,10,8,8であるが,それらが以下のコマンドの c(...) に対応する.データの読み込み,およびその後の解析は,分散分析を行うときと同様にエクセルやテキストファイル形式のデータをデータフレームとしてそのまま読み込む方法でも良い.
$ vg=factor(rep(c("A", "B", "C", "D"), c(9, 10, 8, 8)))
これらを用いて,対応がない場合のボンフェローニ補正をした中央値の検定は以下のように行う.paired= には T か F を指定し,これに F を指定した場合は対応がない検定が行われる.p.adjust.method= には使用したい補正法の名称を文字列で指定する.
$ pairwise.wilcox.test(x=vx,g=vg,paired=F,p.adjust.method="bonferroni")
これを実行した結果は以下のようになる.行列形式で表される値がそれぞれの群間における p値である.有意水準を5%と設定した場合だと,このとき,データAとデータBの間にのみ中央値の差があると結論する.
Pairwise comparisons using Wilcoxon rank sum test data: vx and vg A B C B 0.00013 - - C 0.00592 1.00000 - D 1.00000 0.20568 0.49790 P value adjustment method: bonferroni
平均値の差の検定を行いたい場合,すなわちボンフェローニ補正したt検定を行いたい場合は,pairwise.wilcox.test() でなくて,pairwise.t.test() を使用する.使用方法は全く同じ.
対応がある場合のボンフェローニ補正をした中央値の検定
以下のサンプルサイズがそれぞれ,8,8,8,8からなる互いに対応があるデータAからDが得られたときの各データ間における中央値の値に差があるかどうかをボンフェローニ補正をしたウィルコクソンの符号順位検定で検定する.帰無仮説 (H0) は各群の中央値はそれぞれの群間で等しいこととなる.
データA | 301, 311, 320, 291, 388, 412, 325, 361 |
データB | 197, 180, 247, 248, 199, 179, 134, 163 |
データC | 209, 302, 187, 166, 234, 290, 175, 116 |
データD | 342, 216, 316, 186, 324, 145, 254, 228 |
最初に,以下のコマンドにて,上のデータAからDをまとめて変数 vx に格納する.
$ vx=c(301, 311, 320, 291, 388, 412, 325, 361, 197, 180, 247, 248, 199, 179, 134, 163, 209, 302, 187, 166, 234, 290, 175, 116, 342, 216, 316, 186, 324, 145, 254, 228)
次に,以下のコマンドにて,読み込んだ各データのグループ名を指定するための変数 vg を作成する.データAからDの標本サイズはそれぞれ,8,8,8,8であるが,それらが以下のコマンドの c(...) に対応する.これらの操作は,上の対応がない場合と完全に一致する.
$ vg=factor(rep(c("A", "B", "C", "D"), c(8, 8, 8, 8)))
これらを用いて,対応がある場合のボンフェローニ補正をした中央値の検定は以下のように行う.この場合は,データ間に対応があるので paired= には T を指定しする.
$ pairwise.wilcox.test(x=vx,g=vg,paired=T,p.adjust.method="bonferroni")
これを実行した結果は以下のようになる.行列形式で表される値がそれぞれの群間における p値である.有意水準を5%と設定した場合だと,このとき,データAとB,データAとCのそれぞれの間に中央値の差があると結論する.
Pairwise comparisons using Wilcoxon signed rank test data: vx and vg A B C B 0.047 - - C 0.047 1.000 - D 0.141 0.328 1.000 P value adjustment method: bonferroni
この場合も,t検定を行いたい場合は,pairwise.wilcox.test() でなくて,pairwise.t.test() を使用する.
オプション一覧
補正の方法は p.adjust.method= にて指定するが,これには bonferroni 以外にも以下で示す多くのものが使える.ボンフェローニの方法はかなり保守的な方法で有意差は出にくくなる.以下のオプションは全て,pairwise.t.test() と pairwise.wilcox.test() に共通して指定できる.
オプション | 詳細 |
---|---|
x | 検定の対象としたいデータをベクトル形式で指定. |
g | 入力するデータに対するグループの情報を持つベクトルデータを指定. |
paired | T または F で指定.T の場合,対応のある検定,F の場合,対応のない検定が実行される. |
p.adjust.method | 使用した補正法を表す文字列を指定.ボンフェローニ補正をしたい場合はbonferroni,ホルムの方法 (1979年) を使いたい場合は holm,ホックバーグの方法 (1988年) を使いたい場合は hochberg,ホンメルの方法 (1988年) を使いたい場合は hommel,ベンジャミニ・ホックバーグの方法を使いたいときは BH または fdr,ベンジャミニ・イェクティエリの方法を使用したいときは BY を指定する. |