[R] Possible bug in lmer nested analysis with factors

[Yan Wong]

On 16 Sep 2005, at 17:21, Sundar Dorai-Raj wrote:

> My guess is he wants this:
>
> c1 <- factor(c)
> d1 <- factor(d)
> m <- lmer(a ~ b + (1|c1:d1)+(1|c1))
>
> which assumes d1 is nested within c1.
>
> Take a look at Section 3 in the "MlmSoftRev" vignette:
>
> library(mlmRev)
> vignette("MlmSoftRev")

Ah, that vignette is extremely useful: it deserves to be more widely  
known.
I mainly intended this reply to be a thank you to yourself and Harold.

Unfortunately, there is (as always), one last thing that is still  
puzzling me:
the df for fixed factors seems to vary between what I currently  
understand to
be equivalent calls to lme and lmer, e.g:

-------

a<-rnorm(36);
b<-factor(rep(1:3,each=12))
c<-factor(rep(1:2,each=6,3))
d<-factor(rep(1:3,each=2,6))
c <- evalq(b:c)[,drop=T] #make unique factor levels
d <- evalq(c:d)[,drop=T] #make unique factor levels

summary(lme(a ~ b, random=~1|c/d))
#  output includes estimates for fixed effects such as
#                    Value Std.Error DF    t-value p-value
#  (Intercept)  0.06908901 0.3318330 18  0.2082041  0.8374
#  b2           0.13279084 0.4692828  3  0.2829655  0.7956
#  b3          -0.26146698 0.4692828  3 -0.5571630  0.6163

# I understand the above lme model to be equivalent to
summary(lmer(a ~ b + (1|c) +(1|c:d))
#but this gives fixed effects estimates with differing DF:
#               Estimate Std. Error DF t value Pr(>|t|)
#  (Intercept)  0.069089   0.331724 33  0.2083   0.8363
#  b2           0.132791   0.469128 33  0.2831   0.7789
#  b3          -0.261467   0.469128 33 -0.5573   0.5811

Again, many thanks for your help: even more so if you or anyone
else can answer this last little niggle of mine.

Yan