[R] Random intercept model with time-dependent covariates, results different from SAS
Christoph Buser
buser at stat.math.ethz.ch
Wed Dec 22 08:42:15 CET 2004
Answering on a mail from
>From Keith Wong <keithw_at_med.usyd.edu.au>
Date Sun 04 Jul 2004 - 17:21:36 EST
Subject [R] Random intercept model with time-dependent
covariates, results different from SAS
Hi all
I've got a question about the degrees of freedom in a mixed model,
calculated with lme from the lme4 package.
Since I've no access to the original data of Keith Wong, I generated a
dataframe with the same structure (but balanced) by
> set.seed(1)
> dat <- data.frame(y = rnorm(50), group = rep(c("a","b"), each = 25),
time = factor(rep(1:5,10)), w = rnorm(50), z = rnorm(50),
id = factor(rep(1:10, each = 5)))
> str(dat)
The subject id is nested in group. We have 5 repeated measures for each
subject.
Surly there is no meaningfull interpretation for the results but for
demonstration it will do well.
I do the model without the covariates:
> library(nlme)
> options(contrasts = c("contr.SAS", "contr.poly"))
> g2 <- lme(y ~ time+group+time:group, data = dat, random = ~ 1 | id)
> anova(g2)
Analysis of Variance Table
numDF denDF F-value p-value
(Intercept) 1 32 0.6340507 0.4317
time 4 32 0.1103619 0.9780
group 1 8 0.2924309 0.6034
time:group 4 32 0.4614766 0.7634
I quit R and restart it and now use library(lme4)
> library(lme4)
> options(contrasts = c("contr.SAS", "contr.poly"))
> g2 <- lme(y ~ time+group+time:group, data = dat, random = ~ 1 | id)
> anova(g2)
Analysis of Variance Table
Df Sum Sq Mean Sq Denom F value Pr(>F)
time 4 0.351 0.088 40.000 0.1104 0.9782
group 1 0.233 0.233 40.000 0.2925 0.5916
time:group 4 1.468 0.367 40.000 0.4615 0.7635
I get other degrees of freedom for the denominator. How can I tell lme that
id is nested in the fixed factor group to get 8 as denominator DF for
group? In my example the difference is small (I've generated the data
randomly) but in other more meaningfull analysis the different DF can
change the results.
I Use R 2.0.1 under Linux.
Package: nlme
Version: 3.1-53
Package: lme4
Version: 0.6-11
In the original message the example was also calculated with SAS and there
was the problem of the degrees of freedom for group, too, but the syntax
wasn't correct. One must specify
> PROC MIXED;
> CLASS id time group;
> MODEL y = time group time*group /solution;
> RANDOM int /subject = id(group);
> RUN;
It is important to tell SAS that the subject id is nested in group and you
can do this by using
subject = id(group)
instead of
subject = id
Then you will get the correct degrees of freedom for the test of group:
Type 3 Tests of Fixed Effects
Num Den
Effect DF DF F Value Pr > F
time 4 32 0.11 0.9780
group 1 8 0.29 0.6033
time*group 4 32 0.46 0.7634
Thanks in advance for an answer. I'm interested if I must use the lme
syntax in another way for the lme4 package.
Regards,
Christoph Buser
More information about the R-help
mailing list