対照群と処理群の比較
ダネットの方法　　　　　　　　　Last modified: Aug 26, 2015

　各処理群ごとに，対照群と平均値に差があるかどうかの検定のみに関心がある場合には，ダネットの方法（Dunnett 1955，1964）が有効である。しかし，あらゆる平均値の比較を行う場合には，最初から対比較・線形比較を用いるべきである。

例題：

　「5 種類の飼料で飼育した魚の体長の平均値が表 1 のようであった。飼料 A 〜飼料 D のそれぞれを対照群と比較して平均値に差があるといえいるか，有意水準 5% で検定しなさい。」

表 1．飼料の違いによる魚の体長

　匹数　 　平均値　 　不偏分散　
　対照群　 20　 25.7　 2.5　
　飼料 A　 20　 35.6　 5.6　
　飼料 B　 20　 32.7　 4.6　
　飼料 C　 20　 29.4　 4.9　
　飼料 D　 20　 27.9　 3.7　
　全体　　 100　 30.26　 16.53　

検定手順:

  1. 前提

  2. 1 個の対照群と，$p$ 個の実験群の，合計 $k = 1 + p$ 群のデータがあるとする。

  3. ある実験群の平均値 $\bar{X}_t$ を，対照群の平均値 $\bar{X}_c$ と比較する。

  4. 群内平均平方（群内分散）を $V_w$，その自由度を $\nu$ とする。

  5. 両側検定用のダネットの表（$\alpha= 0.05$，$\alpha= 0.01$，あるいは片側検定用の表 $\alpha= 0.05$，$\alpha= 0.01$）から，$\nu$，$k$，$\alpha$ に対応する値を読みとり，これを $d$ とする。
    なお，対応する自由度が表にない場合には，自由度の逆数で比例配分して求める。

  6. 実験群の平均値 $\bar{X}_t$ と対照群の平均値 $\bar{X}_c$ の差の $(1-\alpha) 100 \%$ 信頼限界を次式で求める。

    \[ \left | \bar{X}_t - \bar{X}_c \right | \pm d \sqrt{V_w}\ \sqrt{\frac{1}{n_t}+\frac{1}{n_c}} \]

  7. 帰無仮説の採否を決める。

    $P$ 値を直接求める方法は，[R による計算プログラム] を参照のこと。
全ての実験群と対照群の差についての検定結果は，指定された有意水準のもとで正しい。

例題の解：

  1. 一元配置分散分析の結果，全体として平均値の差があることがわかる。
    $F_0 = 72.26761$，自由度 $= (4，95)$，有意確率 $P < 0.0001$。

    > n.i <- rep(20, 5)
> Mean.i <- c(25.7, 35.6, 32.7, 29.4, 27.9)
> Var.i <- c(2.5, 5.6, 4.6, 4.9, 3.7)
> k <- length(n.i)
> n <- sum(n.i)
> Mean <- sum(n.i*Mean.i)/n
> Sb <- sum(n.i*(Mean.i-Mean)^2)
> Sw <- sum((n.i-1)*Var.i)
> St <- (Sb+Sw)/(n-1)
> dfb <- k-1
> dfw <- n-k
> Vb <- Sb/dfb
> Vw <- Sw/dfw
> F <- Vb/Vw
> P <- pf(F, dfb, dfw, lower.tail=FALSE)
> Ps <- format.pval(P)
> cat(sprintf("F = %.5f(%i, %i), P %s%s\n", F, dfb, dfw, ifelse(substring(Ps, 1, 1)=="<", "", "="), Ps))
F = 72.26761(4, 95), P < 2.22e-16
  2. 群内分散 $V_w = 4.26$，その自由度は $95$ である。

    > Vw
[1] 4.26
> dfw
[1] 95
  3. 両側検定用のダネットの表（$\alpha$= 0.05）から（補間法により），$d = 2.48052632$ を得る。 
    > nu.a <- 60
> a <- 2.51
> nu.c <- 120
> c <- 2.47
> nu.b <- 95
> d <- a-(a-c)*(1/nu.a-1/nu.b)/(1/nu.a-1/nu.c)
> d
[1] 2.480526
  4. 多重比較の結果は以下のようになる。 
    > diff <- abs(Mean.i[2:5]-Mean.i[1])
> SE <- sqrt(Vw)*sqrt(1/n.i[2:5]+1/n.i[1])
> lower <- diff-d*SE
> upper <- diff+d*SE
> data.frame(diff, lower, upper)
  diff     lower     upper
1  9.9 8.2809933 11.519007
2  7.0 5.3809933  8.619007
3  3.7 2.0809933  5.319007
4  2.2 0.5809933  3.819007
    　比較する群　 　平均値の差　 95% 信頼区間 判定
    　飼料 A　 9.9 ［ 8.280993 ，11.51901 ］ 　有意な差である　
    　飼料 B　 7.0 ［ 5.380993 ，8.619007 ］ 　有意な差である　
    　飼料 C　 3.7 ［ 2.080993 ，5.319007 ］ 　有意な差である　
    　飼料 D　 2.2 ［ 0.580993 ，3.819007 ］ 　有意な差である　

    [R による計算プログラム]によれば， 

    > # データを再現する（本当のデータではないが，平均値，不偏分散は表 1 に示されたものと同じになる）
> n.i <- rep(20, 5)
> Mean.i <- c(25.7, 35.6, 32.7, 29.4, 27.9)
> Var.i <- c(2.5, 5.6, 4.6, 4.9, 3.7)
> data <- mapply(function(n, v, m) c(scale(rnorm(n))*sqrt(v)+m), n.i, Var.i, Mean.i)
> data
          [,1]     [,2]     [,3]     [,4]     [,5]
 [1,] 27.44489 36.03663 32.32625 28.82010 29.34327
 [2,] 25.17521 38.90546 33.38308 30.33822 28.04991
 [3,] 24.39793 34.59178 30.53276 26.72185 27.96165
 [4,] 28.55847 38.04829 29.87004 30.91803 28.49038
 [5,] 26.01307 33.39172 31.55247 29.00199 26.73114
 [6,] 24.07109 33.36879 34.25772 29.05738 30.69118
 [7,] 24.78360 33.65421 33.11500 29.63991 29.62892
 [8,] 22.75709 31.76294 36.17742 28.82967 25.76821
 [9,] 26.01198 34.00233 30.27184 26.70152 27.92164
[10,] 25.99516 40.20536 32.26772 31.04629 27.66308
[11,] 27.39708 34.78678 37.57392 29.97032 27.74809
[12,] 26.00918 38.86186 32.11273 24.23982 30.63102
[13,] 26.96409 35.71035 33.68013 31.98496 26.00846
[14,] 24.64677 34.75720 33.60452 28.63074 29.41766
[15,] 24.48767 36.27121 34.97363 30.28508 28.07175
[16,] 25.24359 33.39627 32.17666 30.46506 27.44650
[17,] 28.83124 38.18133 34.24452 27.95919 30.48866
[18,] 25.24484 32.64976 28.93742 27.10572 27.96583
[19,] 26.30222 37.41242 30.72537 32.77763 24.26720
[20,] 23.66484 36.00534 32.21679 33.50652 23.70544

> apply(data, 2, mean)
[1] 25.7 35.6 32.7 29.4 27.9

> apply(data, 2, var)
[1] 2.5 5.6 4.6 4.9 3.7

> data2 <- c(data)
> data2
  [1] 27.44489 25.17521 24.39793 28.55847 26.01307 24.07109 24.78360 22.75709 26.01198
 [10] 25.99516 27.39708 26.00918 26.96409 24.64677 24.48767 25.24359 28.83124 25.24484
 [19] 26.30222 23.66484 36.03663 38.90546 34.59178 38.04829 33.39172 33.36879 33.65421
 [28] 31.76294 34.00233 40.20536 34.78678 38.86186 35.71035 34.75720 36.27121 33.39627
 [37] 38.18133 32.64976 37.41242 36.00534 32.32625 33.38308 30.53276 29.87004 31.55247
 [46] 34.25772 33.11500 36.17742 30.27184 32.26772 37.57392 32.11273 33.68013 33.60452
 [55] 34.97363 32.17666 34.24452 28.93742 30.72537 32.21679 28.82010 30.33822 26.72185
 [64] 30.91803 29.00199 29.05738 29.63991 28.82967 26.70152 31.04629 29.97032 24.23982
 [73] 31.98496 28.63074 30.28508 30.46506 27.95919 27.10572 32.77763 33.50652 29.34327
 [82] 28.04991 27.96165 28.49038 26.73114 30.69118 29.62892 25.76821 27.92164 27.66308
 [91] 27.74809 30.63102 26.00846 29.41766 28.07175 27.44650 30.48866 27.96583 24.26720
[100] 23.70544

> group <- rep(1:5, each=20)
> group
  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3
 [42] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 5 5
 [83] 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5

> library(mvtnorm)
> dunnett(data2, group) # P 値は実行のたびに少し変化する
            t            p
1:2 15.168072 0.000000e+00
1:3 10.724900 0.000000e+00
1:4  5.668876 1.010806e-06
1:5  3.370683 4.062457e-03

    multcomp の glht を使う場合 

    > library(multcomp)
> factor.group <- factor(group)
> ans <- glht(aov(data2~factor.group), linfct=mcp(factor.group="Dunnett"), abseps=0.0001)
> summary(ans) # P 値は実行のたびに少し変化する

	 Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Dunnett Contrasts


Fit: aov(formula = data2 ~ factor.group)

Linear Hypotheses:
           Estimate Std. Error t value Pr(>|t|)
2 - 1 == 0   9.9000     0.6527  15.168  < 0.001
3 - 1 == 0   7.0000     0.6527  10.725  < 0.001
4 - 1 == 0   3.7000     0.6527   5.669  < 0.001
5 - 1 == 0   2.2000     0.6527   3.371  0.00411
(Adjusted p values reported -- single-step method)

演習問題：

応用問題：

・ 計算プログラム　[R]
・ 直前のページへ戻る　　・ E-mail to Shigenobu AOKI