1
$\begingroup$

First off, I have read a plethora of answers on this site but can't seem to find something that satisfies my situation below:

library(lme4)

> dput(test2)
structure(list(Subject = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 
5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 
7L, 7L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 9L, 9L, 9L, 9L, 9L, 
9L, 9L, 9L, 9L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 
11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 12L, 12L, 12L, 12L, 
12L, 12L, 12L, 12L, 12L, 13L, 13L, 13L, 13L, 13L, 13L, 13L, 13L, 
13L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 15L, 15L, 15L, 15L, 
15L, 15L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 17L, 17L, 
17L, 17L, 17L, 17L, 17L, 18L, 18L, 18L, 18L, 18L, 18L, 18L, 18L, 
18L, 19L, 19L, 19L, 19L, 19L, 19L, 19L, 19L, 19L, 20L, 20L, 20L, 
20L, 20L, 20L, 20L, 20L, 20L, 21L, 21L, 21L, 21L, 21L, 21L, 21L, 
21L, 21L, 22L, 22L, 22L, 22L, 22L, 22L, 22L, 22L, 22L, 23L, 23L, 
23L, 23L, 23L, 23L, 23L, 23L, 23L), .Label = c("153", "204", 
"212", "790", "1422", "1427", "1430", "1507", "1508", "1511", 
"1512", "1602", "1609", "1618", "5327", "5466", "5753", "5935", 
"6424", "7004", "8339", "10806", "12802"), class = "factor"), 
    Status = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
    2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L), .Label = c("Con", "Sub", "Cl"), class = "factor"), 
    Stim = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 
    1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 
    3L, 3L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 
    2L, 3L, 3L, 3L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 
    2L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 
    2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 
    1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 
    3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 2L, 2L, 2L, 
    3L, 3L, 3L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 
    3L, 3L, 3L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 
    2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 
    1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 
    3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 
    2L, 3L, 3L, 3L), .Label = c("B", "C", "D"), class = "factor"), 
    Treatment = structure(c(3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 
    3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 
    3L, 2L, 1L, 3L, 2L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 3L, 2L, 
    1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 
    1L, 3L, 2L, 1L, 3L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 2L, 1L, 3L, 
    2L, 1L, 3L, 2L, 1L), .Label = c("None", "V2", "V1"), class = "factor"), 
    value = c(2.056, 5.819, 5.371, 15.264, 8.942, 7.934, 9.616, 
    7.752, 6.362, 1.416, 4.82, 5.394, 6.053, 3.23, 3.866, 1.334, 
    1.152, 1.152, 2.165, 2.267, 2.662, 6.26, 3.153, 5.022, 2.455, 
    1.023, 1.239, 4.024, 7.428, 41.86, 9.878, 9.313, 10.056, 
    4.729, 4.776, 1.514, 3.085, 3.845, 5.605, 5.411, 4.991, 1.771, 
    1.418, 2.934, 19.516, 11.404, 11.081, 30.123, 13.1, 13.57, 
    1.612, 0.903, 1.622, 3.293, 1.23, 1.277, 2.185, 4.405, 5.247, 
    22.581, 8.429, 7.529, 7.651, 6.508, 8.511, 1.427, 1.719, 
    1.884, 1.502, 1.436, 2.026, 0.376, 0.631, 0.284, 2.653, 2.307, 
    1.416, 15.086, 5.152, 2.666, 10.176, 3.574, 3.349, 6.896, 
    4.051, 5.826, 14.228, 6.891, 5.277, 10.309, 4.732, 7.404, 
    2.054, 2.602, 1.856, 8.099, 4.793, 5.258, 1.948, 2.016, 1.72, 
    3.523, 5.142, 5.03, 17.116, 7.487, 7.235, 8.407, 3.585, 3.04, 
    4.506, 4.026, 10.015, 3.807, 8.342, 4.873, 4.21, 5.79, 20.031, 
    3.12, 4.239, 42.889, 18.702, 19.384, 19.981, 23.454, 18.292, 
    62.969, 19.677, 25.142, 55.025, 30.228, 20.731, 18.017, 22.562, 
    8.579, 8.663, 17.277, 7.495, 7.351, 38.814, 24.896, 30.968, 
    39.639, 39.108, 18.909, 65.557, 56.319, 32.419, 7.28, 2.734, 
    1.494, 5.23, 4.252, 1.923, 9.488, 2.662, 1.644, 6.813, 8.478, 
    6.799, 33.731, 13.292, 10.854, 12.287, 3.597, 4.815, 33.283, 
    7.684, 8.271, 37.239, 9.026, 9.107, 38.36, 28.966, 32.334, 
    2.691, 2.737, 2.949, 3.163, 4.345, 5.337, 5.419, 9.01, 8.911, 
    1.155, 0.485, 0.863, 5.354, 0.697, 1.304, 4.537, 2.098, 1.517
    )), row.names = c(4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 
16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L, 28L, 29L, 30L, 31L, 
32L, 33L, 34L, 35L, 36L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 
51L, 52L, 53L, 54L, 55L, 56L, 57L, 58L, 59L, 66L, 67L, 68L, 69L, 
70L, 71L, 75L, 76L, 77L, 78L, 79L, 80L, 87L, 88L, 89L, 90L, 91L, 
92L, 93L, 94L, 95L, 99L, 100L, 101L, 102L, 103L, 104L, 105L, 
106L, 107L, 111L, 112L, 113L, 114L, 115L, 116L, 117L, 118L, 119L, 
123L, 124L, 125L, 126L, 127L, 128L, 129L, 130L, 131L, 135L, 136L, 
137L, 138L, 139L, 140L, 141L, 142L, 143L, 147L, 148L, 149L, 150L, 
151L, 152L, 153L, 154L, 155L, 159L, 160L, 162L, 163L, 164L, 165L, 
166L, 167L, 174L, 175L, 176L, 177L, 178L, 179L, 183L, 184L, 185L, 
186L, 187L, 188L, 189L, 190L, 191L, 195L, 198L, 199L, 200L, 201L, 
202L, 203L, 207L, 208L, 209L, 210L, 211L, 212L, 213L, 214L, 215L, 
219L, 220L, 221L, 222L, 223L, 224L, 225L, 226L, 227L, 231L, 232L, 
233L, 234L, 235L, 236L, 237L, 238L, 239L, 243L, 244L, 245L, 246L, 
247L, 248L, 249L, 250L, 251L, 255L, 256L, 257L, 258L, 259L, 260L, 
261L, 262L, 263L, 267L, 268L, 269L, 270L, 271L, 272L, 273L, 274L, 
275L), class = "data.frame")

Model:

mod.lmer <- lmer(value ~ Status + Treatment + Stim + Treatment*Status*Stim + (1|Subject),
                  data = test2)

Notably, the following things do not produce similar results:

  • anova() type 3 and Anova() type 3
  • summary() and anova() type 3

I have tried options(contrasts = c("contr.sum", "contr.poly")) and it does not change any results. I have run my pre-determined contrasts of interest using emmeans, which had high significance, so I'm having a hard time understanding why the Anova() type 3 and summary() aren't showing any significance/similar results. I was initially going to use type 3 Anova for my data due to the contrasts I'm interested in and mildly unbalanced data, but now I'm not sure which anova to use for my report.

> anova(mod.lmer, type = 1)
Type I Analysis of Variance Table with Satterthwaite's method
                       Sum Sq Mean Sq NumDF   DenDF F value    Pr(>F)    
Status                 213.07  106.53     2  19.999  2.7098 0.0909098 .  
Treatment             1796.00  898.00     2 147.061 22.8417 2.298e-09 ***
Stim                   622.26  311.13     2 147.647  7.9139 0.0005433 ***
Status:Treatment       549.26  137.32     4 147.067  3.4928 0.0093286 ** 
Treatment:Stim         584.18  146.04     4 147.072  3.7148 0.0065363 ** 
Status:Stim            385.45   96.36     4 147.755  2.4511 0.0486018 *  
Status:Treatment:Stim   59.95    7.49     8 147.088  0.1906 0.9918876    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

> anova(mod.lmer, type = 3)
Type III Analysis of Variance Table with Satterthwaite's method
                       Sum Sq Mean Sq NumDF   DenDF F value    Pr(>F)    
Status                 217.10  108.55     2  20.041  2.7611  0.087273 .  
Treatment             1837.07  918.53     2 147.097 23.3640 1.543e-09 ***
Stim                   523.71  261.86     2 147.852  6.6606  0.001699 ** 
Status:Treatment       564.47  141.12     4 147.093  3.5895  0.007990 ** 
Treatment:Stim         624.52  156.13     4 147.093  3.9713  0.004331 ** 
Status:Stim            387.16   96.79     4 147.777  2.4620  0.047779 *  
Status:Treatment:Stim   59.95    7.49     8 147.088  0.1906  0.991888    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

> Anova(mod.lmer, type = 3)
Analysis of Deviance Table (Type III Wald chisquare tests)

Response: value
                       Chisq Df Pr(>Chisq)
(Intercept)           1.7023  1     0.1920
Status                0.8426  2     0.6562
Treatment             0.2184  2     0.8966
Stim                  0.4451  2     0.8005
Status:Treatment      3.2780  4     0.5124
Treatment:Stim        3.2279  4     0.5204
Status:Stim           2.6306  4     0.6214
Status:Treatment:Stim 1.5248  8     0.9923

Summary:

> summary(mod.lmer)
    Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
    Formula: value ~ Status + Treatment + Stim + Treatment * Status * Stim +      (1 | Subject)
       Data: test2

REML criterion at convergence: 1196.2

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.0551 -0.3745  0.0005  0.2780  3.6174 

Random effects:
 Groups   Name        Variance Std.Dev.
 Subject  (Intercept) 74.54    8.634   
 Residual             39.31    6.270   
Number of obs: 194, groups:  Subject, 23

Fixed effects:
                             Estimate Std. Error        df t value Pr(>|t|)
(Intercept)                   4.89297    3.75025  48.48368   1.305    0.198
StatusSub                     2.39119    5.92696  56.30030   0.403    0.688
StatusCl                      5.05503    5.50718  43.78039   0.918    0.364
TreatmentV2                  -0.15767    3.25774 147.17807  -0.048    0.961
TreatmentV1                  -1.36205    3.25774 147.17807  -0.418    0.676
StimC                         1.70903    3.18599 147.39224   0.536    0.592
StimD                        -0.06519    3.18599 147.39224  -0.020    0.984
StatusSub:TreatmentV2        -0.81671    5.34245 147.25141  -0.153    0.879
StatusCl:TreatmentV2          0.27653    4.67391 147.10108   0.059    0.953
StatusSub:TreatmentV1         1.22546    5.23640 147.42567   0.234    0.815
StatusCl:TreatmentV1          7.13076    4.67391 147.10108   1.526    0.129
TreatmentV2:StimC            -0.34855    4.39879 147.11046  -0.079    0.937
TreatmentV1:StimC             6.53316    4.39879 147.11046   1.485    0.140
TreatmentV2:StimD            -0.19866    4.39879 147.11046  -0.045    0.964
TreatmentV1:StimD             4.00316    4.39879 147.11046   0.910    0.364
StatusSub:StimC              -3.95019    5.11692 147.69916  -0.772    0.441
StatusCl:StimC               -1.28903    4.62419 147.20128  -0.279    0.781
StatusSub:StimD              -0.97558    5.24777 148.10277  -0.186    0.853
StatusCl:StimD                4.74162    4.62419 147.20128   1.025    0.307
StatusSub:TreatmentV2:StimC   1.71993    6.96497 147.15964   0.247    0.805
StatusCl:TreatmentV2:StimC    2.77555    6.46642 147.06631   0.429    0.668
StatusSub:TreatmentV1:StimC   6.14800    6.88397 147.25849   0.893    0.373
StatusCl:TreatmentV1:StimC    4.09084    6.46642 147.06631   0.633    0.528
StatusSub:TreatmentV2:StimD   0.66771    7.09809 147.15476   0.094    0.925
StatusCl:TreatmentV2:StimD    4.43823    6.46642 147.06631   0.686    0.494
StatusSub:TreatmentV1:StimD   4.23759    7.01862 147.24976   0.604    0.547
StatusCl:TreatmentV1:StimD    2.84270    6.46642 147.06631   0.440    0.661

Correlation matrix not shown by default, as p = 27 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
$\endgroup$
1
  • $\begingroup$ It's not enough to just change that contrasts option. You have to re-fit the model afterward. Did you do that? $\endgroup$
    – Russ Lenth
    Commented Mar 28, 2021 at 12:47

1 Answer 1

3
$\begingroup$

I don't get the same results you do from car::Anova. I have:

> options(contrasts = c("contr.sum", "contr.poly"))
> mod.lmer <- lmer(value ~ Status + Treatment + Stim + Treatment*Status*Stim + (1|Subject),
+                  data = test2)

> car::Anova(mod.lmer, type=3)
Analysis of Deviance Table (Type III Wald chisquare tests)

Response: value
                        Chisq Df Pr(>Chisq)
(Intercept)           29.8855  1  4.583e-08
Status                 5.5222  2   0.063224
Treatment             46.7280  2  7.131e-11
Stim                  13.3213  2   0.001280
Status:Treatment      14.3580  4   0.006236
Treatment:Stim        15.8853  4   0.003177
Status:Stim            9.8479  4   0.043068
Status:Treatment:Stim  1.5248  8   0.992289

... for which the P values are in the same ballpark as what you show with anova(mod.lmer, type = 3). [BTW, if I try anova(mod.lmer, type = 3), I get a warning message that there is no type argument; so you must have some package loaded that I don't have loaded. So please tell us where you are getting that anova() function.]

I suspect that the contrasts being used may not be what you think they are. The contrasts option requires two arguments, as in my test. Check the contrasts actually used in fitting the model as follows:

> attr(mod.lmer@pp$X, "contrasts")
$Status
[1] "contr.sum"

$Treatment
[1] "contr.sum"

$Stim
[1] "contr.sum"

Another way of getting a type III anova that does not depend on what contrast coding is used is via joint_tests():

> emmeans::joint_tests(mod.lmer)
 model term            df1    df2 F.ratio p.value
 Status                  2  19.89   2.761 0.0875 
 Treatment               2 147.06  23.363 <.0001 
 Stim                    2 147.75   6.654 0.0017 
 Status:Treatment        4 147.03   3.589 0.0080 
 Status:Stim             4 147.39   2.461 0.0478 
 Treatment:Stim          4 147.03   3.971 0.0043 
 Status:Treatment:Stim   8 147.02   0.191 0.9919 

You can't directly compare the results of summary() with an ANOVA, because the summary provides tests of individual coefficients, not omnibus tests.

I don't think type III anovas are particularly useful. For model selection, type II anovas make sense; but type III tests consider things like main effects when two- and three-way interactions involving them are still in the model, and that's like testing the performance of your VW Rabbit, assuming that it is a Maserati.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.