Dunnett の多重比較を R で実行する方法
中澤 港 [nminato@med.gunma-u.ac.jp]

rev 1.00 / 30-Nov-2009

1 はじめに
多重比較とは，3 群以上の間で平均値を比べる場合に，2 群ずつの検定を繰り返すと，全体としての第一種の過誤が，個々の検定で設定する有意水準よりも大きくなってしまうため，そこを調整して検定する方法の総称である。
総当たりで比較するならば，もっとも単純な方法は，全体としての第一種の過誤を有意水準（例えば 0.05）以下にするために，個々の検定の有意水準を，検定の繰り返し回数で割っておくことである。Bonferroni の不等式の応用であり，Bonferroni の方法と呼ばれる。これだと調整しすぎなので，そうならないように配慮した Holm の方法が広く実用されている。

R では，pairwise.t.test() 関数や pairwise.wilcox.test() 関数の p.adjust.method="bon"や p.adjust.method="holm"で，これらが実行できる。
もしデータが正規分布しているなら，Tukey の HSD という方法の適用も考えられる。これは R
では TukeyHSD(aov(X~C)) という形で実行できる。
しかし，総当たりの比較ではなく，対照群と個々の実験群の比較をしたい場合もある。この場合も比較回数に基づいて単純に Bonferroni や Holm の方法を用い，検定の多重性の調整をすることは可能だが，それよりも Dunnett の方法や Williams の方法という，その目的に特化した方法を使う方が検出力が失われない。

2 multcomp ライブラリを使う
R をインストールした後，もしまだ multcomp ライブラリがインストールされていないなら，インターネットに接続した環境で，install.packages("multcomp",dep=T) としてインストールする必要がある。R 本体をインストールする場合と同様，Windows では管理者権限がないとライブラリインストールもできないことに注意が必要である。なお，作業途中で repository を尋ねる
ダイアログが出たら，日本では Japan(Tsukuba) か Japan(Tokyo) を選べばいいと思われる。
インストールが完了したのち，実際に multcomp ライブラリをロードして使うには，R のプロンプトに対して，require(multcomp) とするか，または library(multcomp) とする。これで準備

1

は完了である。
さて，いよいよ実際に Dunnett の多重比較を行うわけだが，R-2.6.2 より前（multcomp ライ
ブラリ 0.993 より前）は，simtest() という関数で，簡単に Dunnett の多重比較が実行できていた。かなり前から，ヘルプには「simtest() は古いユーザインターフェースなので，glht() とsummary.glht() を使え」と書かれていたが，simtest() の方が使いやすかったと思う。しかし，
R-2.6.2（multcomp-0.993) 以降，そのインターフェースは廃止されたので，以下の通り，新しいインターフェースを使うしかない。
対照群と複数の実験群の比較というデザインで使う Dunnett の多重比較の実行方法を手っ取り早く知るには，とりあえず，example(recovery) として表示される結果をみることをお薦めする。
検定の多重性を調整した有意確率を出すには glht() の結果のオブジェクトの summary() をとればいいことがわかる。

3 例示
この文書では，ごく簡単な実例で説明する。5 人ずつ 3 群にランダムに分けた高血圧患者がいて，他の条件（食事療法，運動療法など）には差をつけずに，プラセボを 1 ヶ月服用した群の収縮期血圧 (mmHg 単位) の低下が 5, 8, 3, 10, 15 で，代表的な薬を 1 ヶ月服用した群の低下は 20, 12,30, 16, 24 で，新薬を 1 ヶ月服用した群の低下が 31, 25, 17, 40, 23 だったとしよう（下表）。
グループ
プラセボ群

5831015

代表薬群

2012301624

新薬群

3125174023

このとき，プラセボ群を対照として，代表的な薬での治療及び新薬での治療に有意な血圧降下作用の差が出るかどうかを見たい（悪くなるかもしれないので両側検定で）という場合に，Dunnettの多重比較を使う。

R でこのデータを bpdown というデータフレームに入力して Dunnett の多重比較をするためには，次のコードを実行すればよい。要点は，一元配置分散分析を行う関数 aov(X~C) の結果のオブジェクトを，multcomp ライブラリをロードした後で glht() 関数の第１引数として渡すことと，その関数の第２引数として，linfct=mcp(群分け変数名="Dunnett") と指定することである。

2

¶

³

bpdown <- data.frame(
medicine=factor(c(rep(1,5),rep(2,5),rep(3,5)), labels=c("プ ラ セ ボ","代 表 薬","新薬")),sbpchange=c(5, 8, 3, 10, 15, 20, 12, 30, 16, 24, 31, 25, 17, 40, 23))
plot(sbpchange ~ medicine, data=bpdown)
res1 <- aov(sbpchange ~ medicine, data=bpdown)summary(res1)
library(multcomp)res2 <- glht(res1, linfct = mcp(medicine = "Dunnett"))
confint(res2, level=0.95)
summary(res2)µ´得られる結果を次の枠内に示す。

¶summary(res1) の結果³

Df Sum Sq Mean Sq F value
Pr(>F)
2 926.8
463.4 9.4765 0.003395 **12 586.8
48.9

medicine
Residuals--Signif. codes:

0***0.001**0.01*0.05.

0.11µ´

この一元配置分散分析の結果から，薬の種類が血圧降下に対して有意な影響を与えていることがわかる。

3¶confint(res2) の結果³

Simultaneous Confidence Intervals
Multiple Comparisons of Means: Dunnett Contrasts

Fit: aov(formula = sbpchange ~ medicine, data = bpdown)
Estimated Quantile = 2.5023
95% family-wise confidence level

Linear Hypotheses:
Estimate lwr
upr
1.1330 23.2670
19.0000
7.9330 30.0670

代表薬 - プラセボ == 0 12.2000
新薬 - プラセボ == 0µ¶summary(res2) の結果´³

Simultaneous Tests for General Linear Hypotheses
Multiple Comparisons of Means: Dunnett Contrasts

Fit: aov(formula = sbpchange ~ medicine, data = bpdown)
Linear Hypotheses:
代表薬 - プラセボ == 0
新薬 - プラセボ == 0

Estimate Std. Error t value Pr(>|t|)12.2004.423
2.759 0.03151 *19.0004.423
4.296 0.00196 **--Signif. codes: 0 *** 0.001 ** 0.01 * 0.05
(Adjusted p values reported -- single-step method).

0.11µ´

この結果から，代表的な薬も新薬も有意水準 5% でプラセボと血圧降下作用に有意差があったと言える。

4