Hi Patrick,

Here is the "classic approach," as done with the R aov() function:

> mma <- aov(Y ~ Fac1 + Error(Fac2), data = mm)
> summary(mma)

Error: Fac2
          Df Sum Sq Mean Sq F value Pr(>F)
Residuals  1    534     534               

Error: Within
          Df Sum Sq Mean Sq F value Pr(>F)
Fac1       2   49.8   24.89    9.92 0.0021
Residuals 14   35.1    2.51               


compare with

> (mmm <- lmer(Y ~ Fac1 + (1 | Fac2), data = mm))
Linear mixed model fit by REML 
Formula: Y ~ Fac1 + (1 | Fac2) 
   Data: mm 
  AIC  BIC logLik deviance REMLdev
 77.1 81.5  -33.5     75.1    67.1
Random effects:
 Groups   Name        Variance Std.Dev.
 Fac2     (Intercept) 59.01    7.68    
 Residual              2.51    1.58    
Number of obs: 18, groups: Fac2, 2

Fixed effects:
             Estimate Std. Error t value
(Intercept) -6.97e-14   5.47e+00    0.00
Fac1B        2.67e+00   9.14e-01    2.92
Fac1C       -1.33e+00   9.14e-01   -1.46

Correlation of Fixed Effects:
      (Intr) Fac1B 
Fac1B -0.084       
Fac1C -0.084  0.500
> anova(mmm)
Analysis of Variance Table
     Df Sum Sq Mean Sq F value
Fac1  2   49.8    24.9    9.92

For more examples of the "classic" approach in R see http://www.personality-project.org/r/r.anova.html and, regarding the absence of p-values in lmer: http://rwiki.sciviews.org/doku.php?id=guides:lmer-tests

MK

On Sep 10, 2013, at 8:57 PM, Patrick Ward <p.ward1@uq.edu.au> wrote:

> I teach statistics at Masters level. I would like to teach mixed models without the Expected Mean Squares (EMS) approach from the old days . But I cannot get "lmer" to give me the right F-statistics.
> 
> 
> Here is the data :
> 
> Y Fac1 Fac2
> 4 A 1
> 6 A 1
> 5 A 1
> 8 B 1
> 10 B 1
> 7 B 1
> 2 C 1
> 5 C 1
> 6 C 1
> -6 A 2
> -5 A 2
> -4 A 2
> 0 B 2
> -4 B 2
> -5 B 2
> -8 C 2
> -7 C 2
> -6 C 2
> 
> 
> Fac1 has three levels (A, B and C) and Fac2 has two levels (1 and 2). I would like to treat Fac1 as fixed and Fac2 as random. The old "EMS formula" technique produces this ANOVA table:
> 
>         Df     Sum Sq       Mean Sq    EMS               F value
> Fac1        2     49.78      24.89      6V1 + 3V1x2 + Ve    24.89/0.89 =  27.97
> Fac2        1    533.56       533.56    9V2 + Ve           533.56/2.78 = 191.93
> Fac1:Fac2   2      1.78     0.89        3V1x2 + Ve            0.8/2.78 =   0.32
> Residuals   12    33.33     2.78        Ve
> 
> 
> The fixed effect uses 0.89 as its denominator.
> 
> I cannot get lmer to give me that F-statistic of 27.97 for the fixed effect (Fac1).
> 
> Here is what I did with lmer:
> 
> fit <- lmer(Y ~ Fac1 + (Fac1 | Fac2))
>    anova(fit)
> 
> 
> Analysis of Variance Table
> 
>     Df Sum Sq Mean Sq F value
> Fac1  2  49.02   24.51  10.294
> 
> 
> This F-statistic of 10.29 is using the denominator 2.38, as can be seen from this summary of the random effects:
> 
> summary(fit)
> 
> Random effects:
> Groups   Name        Variance Std.Dev. Corr
> Fac2     (Intercept) 49.77980 7.05548
>          Fac1B        0.88499 0.94074  1.000
>          Fac1C        0.88500 0.94074  1.000 1.000
> Residual              2.38094 1.54303
> 
> 
> The F-statistic should be using the 0.885 as the denominator, and not 2.38.
> 
> So what am I doing wrong? Where is my thinking wrong? How do I get lmer to use 0.89 as the denominator for the fixed effect F-statistic? Why does lmer gives different F-statistics from the classical approach to mixed modelling? As far as I can tell, it effectively uses the residual Mean Square as the denominator for the F-statistics for the fixed effects, something that was ruled out classically as incorrect.
> 
> 
> My data set is deliberately simple, just to remove arithmetic complexity from the story ; I can produce the classical ANOVA table by hand, and so I know that it is arithmetically correct. So please no questions like "why do you want a factor with two levels to be random?". Thanks.
> 
> 
> Pat Ward
> 
> 	[[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-mixed-models@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models


	[[alternative HTML version deleted]]

