FAQ/rHSD - CBU statistics Wiki
location: FAQ / rHSD

Using Tukey's test for main effects to correct for multiplicity in a repeated measures ANOVA using R

(There is also example in the post-hoc talk demo file available at the Graduate statistics talk page from here)

install.packages(c("MASS"))
install.packages(c("multcomp"))
library(multcomp)
library(MASS)
library(car)
library(nlme)

> oats

B

V

N

Y

1

I

Victory

0.0cwt

111

2

I

Victory

0.2cwt

130

3

I

Victory

0.4cwt

157

4

I

Victory

0.6cwt

174

5

I

Golden.rain

0.0cwt

117

6

I

Golden.rain

0.2cwt

114

7

I

Golden.rain

0.4cwt

161

8

I

Golden.rain

0.6cwt

141

9

I

Marvellous

0.0cwt

105

10

I

Marvellous

0.2cwt

140

11

I

Marvellous

0.4cwt

118

12

I

Marvellous

0.6cwt

156

13

II

Victory

0.0cwt

61

14

II

Victory

0.2cwt

91

15

II

Victory

0.4cwt

97

16

II

Victory

0.6cwt

100

17

II

Golden.rain

0.0cwt

70

18

II

Golden.rain

0.2cwt

108

19

II

Golden.rain

0.4cwt

126

20

II

Golden.rain

0.6cwt

149

21

II

Marvellous

0.0cwt

96

22

II

Marvellous

0.2cwt

124

23

II

Marvellous

0.4cwt

121

24

II

Marvellous

0.6cwt

144

25

III

Victory

0.0cwt

68

26

III

Victory

0.2cwt

64

27

III

Victory

0.4cwt

112

28

III

Victory

0.6cwt

86

29

III

Golden.rain

0.0cwt

60

30

III

Golden.rain

0.2cwt

102

31

III

Golden.rain

0.4cwt

89

32

III

Golden.rain

0.6cwt

96

33

III

Marvellous

0.0cwt

89

34

III

Marvellous

0.2cwt

129

35

III

Marvellous

0.4cwt

132

36

III

Marvellous

0.6cwt

124

37

IV

Victory

0.0cwt

74

38

IV

Victory

0.2cwt

89

39

IV

Victory

0.4cwt

81

40

IV

Victory

0.6cwt

122

41

IV

Golden.rain

0.0cwt

64

42

IV

Golden.rain

0.2cwt

103

43

IV

Golden.rain

0.4cwt

132

44

IV

Golden.rain

0.6cwt

133

45

IV

Marvellous

0.0cwt

70

46

IV

Marvellous

0.2cwt

89

47

IV

Marvellous

0.4cwt

104

48

IV

Marvellous

0.6cwt

117

49

V

Victory

0.0cwt

62

50

V

Victory

0.2cwt

90

51

V

Victory

0.4cwt

100

52

V

Victory

0.6cwt

116

53

V

Golden.rain

0.0cwt

80

54

V

Golden.rain

0.2cwt

82

55

V

Golden.rain

0.4cwt

94

56

V

Golden.rain

0.6cwt

126

57

V

Marvellous

0.0cwt

63

58

V

Marvellous

0.2cwt

70

59

V

Marvellous

0.4cwt

109

60

V

Marvellous

0.6cwt

99

61

VI

Victory

0.0cwt

53

62

VI

Victory

0.2cwt

74

63

VI

Victory

0.4cwt

118

64

VI

Victory

0.6cwt

113

65

VI

Golden.rain

0.0cwt

89

66

VI

Golden.rain

0.2cwt

82

67

VI

Golden.rain

0.4cwt

86

68

VI

Golden.rain

0.6cwt

104

69

VI

Marvellous

0.0cwt

97

70

VI

Marvellous

0.2cwt

99

71

VI

Marvellous

0.4cwt

119

72

VI

Marvellous

0.6cwt

121

Aov.mod <- aov(Y ~ N + V + Error(B/V), data = oats) 
Lme.mod <- lme(Y ~ N + V, random = ~1 | B/V, data = oats) 
summary(Aov.mod) 

Error = B

Df

Sum Sq

Mean Sq

F

Residuals

5

15875

3175.1

Error = B:V

Df

Sum Sq

Mean Sq

F

V

2

1786.4

893.18

1.4853

Residuals

10

6013.3

601.33

Error = Within

Df

Sum Sq

Mean Sq

F

N

3

20020.5

6673.5

41.053

Residuals

51

8290.5

162.6

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

anova(Lme.mod) 

numDF

denDF

F-value

p-value

(Intercept)

1

51

245.14566

<.0001

N

3

51

41.05276

<.0001

V

2

10

1.48534

0.2724

summary(Lme.mod) 

Linear mixed-effects model fit by REML

  • Data: oats

AIC

BIC

logLik

586.0688

605.7756

-284.0344

Random effects:

  • Formula: ~1 | B
    • (Intercept)

StdDev: 14.64486

  • Formula: ~1 | V %in% B
    • (Intercept) Residual

StdDev: 10.47344 12.74987

Fixed effects: Y ~ N + V

Value

Std.Error

DF

t-value

p-value

(Intercept)

79.91667

8.220345

51

9.721814

0.0000

N0.2cwt

19.50000

4.249956

51

4.588283

0.0000

N0.4cwt

34.83333

4.249956

51

8.196164

0.0000

N0.6cwt

44.00000

4.249956

51

10.353049

0.0000

VMarvellous

5.29167

7.078905

10

0.747526

0.4720

VVictory

-6.87500

7.078905

10

-0.971195

0.3544

Correlation:

(Intr)

N0.2cw

N0.4cw

N06.cw

VMrv11

-0.259

-

-

-

-

-0.259

0.500

-

-

-

-0.259

0.500

0.500

-

-

-0.431

0.000

0.000

0.000

-

-0.431

0.000

0.000

0.000

0.500

Standardized Within-Group Residuals:

Min

Q1

Med

Q3

Max

-1.84134117

-0.66279719

-0.06694235

0.63822499

1.66 066729

Number of Observations: 72

Number of Groups:

|||||| B || V %in% B ||

6

18

summary(glht(Lme.mod, linfct=mcp(V="Tukey"))) 
  • Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts

Fit: lme.formula(fixed = Y ~ N + V, data = oats, random = ~1 | B/V)

Linear Hypotheses:

Estimate

Std. Error

z value

Pr(>|z|)

Marvellous - Golden.rain == 0

5.292

7.079

0.748

0.735

Victory - Golden.rain == 0

-6.875

7.079

-0.971

0.595

Victory - Marvellous == 0

-12.167

7.079

-1.719

0.198

(Adjusted p values reported -- single-step method)

Similarly for the within subjects effect

summary(glht(Lme.mod, linfct=mcp(N="Tukey"))) 

None: FAQ/rHSD (last edited 2013-03-08 10:17:27 by localhost)