# Simultaneous Inferences in R

## Using package:multcomp

### mice example

The `multcomp` does not come with the standard R distribution. You need to download it from CRAN. See Install R Packages for example on the installation of R add-on packages.

> require(multcomp) # load package:multcomp

```Loading required package: multcomp
[1] TRUE
```

> mice <- data.frame(group=as.factor(rep(c("C","1","2"),rep(6,3))),

+ score=scan())

1: 58 32 59 64 55 49 73 70 68 71 60 62 53 74 72 62 58 61

19:

```Read 18 items
```

> levels(mice\$group) # show factor levels

```[1] "1" "2" "C"
```

> # We wish to use "C" as the reference level, so do:

> mice\$group <- relevel(mice\$group, ref="C")

> levels(mice\$group) # now "C" becomes the reference level

```[1] "C" "1" "2"
```

> mice.aov <- aov(score ~ group, data=mice) # fit aov model

> summary(mice.aov) # produce the summary of the fit

```            Df  Sum Sq Mean Sq F value  Pr(>F)
group        2     673   336.5   4.549 0.0286 *
Residuals   15    1110    74.0
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

>

>

> # The following shows Tukey's all pairwise comparisons

> # using glht (general linear hypothesis test)

> mice.glht <- glht(mice.aov, linfct = mcp(group="Tukey"))

> summary(mice.glht) # the summary of the tests

```         Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts

Fit: aov(formula = score ~ group, data = mice)

Linear Hypotheses:
Estimate Std. Error t value p value
1 - C == 0   14.500      4.965   2.920  0.0268 *
2 - C == 0   10.500      4.965   2.115  0.1204
2 - 1 == 0   -4.000      4.965  -0.806  0.7054
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Adjusted p values reported -- single-step method)
```

>

> # The results can be displayed graphically (function confint

> # was silently called to produce confidence intervals)

> windows(width=5,height=3,pointsize=10)

> plot(mice.glht)

> mtext("Tukey Honest Significant Differences",side=3,line=0.5)

>

>

> # The following demonstrates Dunnett's method

> mice.Dunnett <- glht(mice.aov, linfct=mcp(group="Dunnett"))

> summary(mice.Dunnett)

```       Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Dunnett Contrasts

Fit: aov(formula = score ~ group, data = mice)

Linear Hypotheses:
Estimate Std. Error t value p value
1 - C == 0   14.500      4.965   2.920  0.0195 *
2 - C == 0   10.500      4.965   2.115  0.0916 .
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(Adjusted p values reported -- single-step method)
```

> windows(width=5,height=3,pointsize=10)

> plot(mice.Dunnett,sub="Mice Data")

> mtext("Dunnet's Method",side=3,line=0.5)

> confint(mice.Dunnett) # show the confidence intervals

```         Simultaneous Confidence Intervals

Multiple Comparisons of Means: Dunnett Contrasts

Fit: aov(formula = score ~ group, data = mice)

Quantile = 2.4392
95% family-wise confidence level

Linear Hypotheses:
Estimate lwr     upr
1 - C == 0 14.5000   2.3884 26.6116
2 - C == 0 10.5000  -1.6116 22.6116
```

> confint(mice.Dunnett,level=0.90) # can change confidence level

```         Simultaneous Confidence Intervals

Multiple Comparisons of Means: Dunnett Contrasts

Fit: aov(formula = score ~ group, data = mice)

Quantile = 2.0663
90% family-wise confidence level

Linear Hypotheses:
Estimate lwr     upr
1 - C == 0 14.5000   4.2399 24.7601
2 - C == 0 10.5000   0.2399 20.7601
```

>

> # The plot below displays 90\% Dunnett's intervals

> plot(confint(mice.Dunnett,level=0.90),sub="Mice Data",

+ ylim=c(0.5,2.5))

> mtext("Dunnett's Method",side=3,line=0.5)

>

>

> # Can also do general contrasts

> contr <- rbind("avg1&2 - C"=c(-1,.5,.5),"2 - 1"=c(0,-1,1))

> confint(glht(mice.aov,linfct=mcp(group=contr)))

```       Simultaneous Confidence Intervals for General Linear Hypotheses

Multiple Comparisons of Means: User-defined Contrasts

Fit: aov(formula = score ~ group, data = mice)

Estimated Quantile = 2.4742

Linear Hypotheses:
Estimate lwr      upr
avg1&2 - C == 0  12.5000   1.8605  23.1395
2 - 1 == 0       -4.0000 -16.2854   8.2854

95% family-wise confidence level
```