[R-sig-ME] random effects specification
Sebastian P. Luque
spluque at gmail.com
Fri Apr 4 01:32:23 CEST 2008
On Thu, 3 Apr 2008 17:32:40 -0500,
"Douglas Bates" <bates at stat.wisc.edu> wrote:
[...]
> It seems that the observations are indexed by subject and treatment so
> the number of levels in the factor treatment:id equals the number of
> observations. You can't estimate a variance for such a term and also
> estimate a residual variance.
> I would start with
> n ~ treatment * community +(1|id)
Yes, the observations are indexed by subject and treatment in the sense
that id levels are the same within treatments of the same community, but
are different among communities. This is a subset of the data:
---<---------------cut here---------------start-------------->---
id community treatment n
1 A 1 13.93
2 A 1 14.42
3 A 1 13.56
1 A 2 14.61
2 A 2 14.74
3 A 2 15.59
1 A 3 13.95
2 A 3 15.21
3 A 3 14.51
1 A 4 13.61
2 A 4 14.99
3 A 4 15.13
4 B 1 14.79
5 B 1 13.41
6 B 1 14.71
4 B 2 14.69
5 B 2 13.46
6 B 2 14.28
4 B 3 14.30
5 B 3 13.18
6 B 3 13.58
4 B 4 14.54
5 B 4 13.25
6 B 4 14.09
---<---------------cut here---------------end---------------->---
Of course, there are many more individuals, but the levels of id differ
among communities, and are the same among treatments. lmer did converge
rapidly with your suggested formula though:
---<---------------cut here---------------start-------------->---
Linear mixed-effects model fit by REML
Formula: n ~ treatment * community + (1 | id)
Data: isotope.m.ph
AIC BIC logLik MLdeviance REMLdeviance
450 481 -216 410 432
Random effects:
Groups Name Variance Std.Dev.
id (Intercept) 0.193 0.439
Residual 0.232 0.481
number of obs: 240, groups: id, 61
Fixed effects:
Estimate Std. Error t value
(Intercept) 14.9748 0.1170 128.0
treatment2 0.0884 0.1222 0.7
treatment3 -0.2829 0.1222 -2.3
treatment4 0.3568 0.1222 2.9
communitysanikiluaq -0.5749 0.1678 -3.4
treatment2:communitysanikiluaq -0.2471 0.1763 -1.4
treatment3:communitysanikiluaq -0.7479 0.1763 -4.2
treatment4:communitysanikiluaq -0.6169 0.1763 -3.5
Correlation of Fixed Effects:
(Intr) trtmn2 trtmn3 trtmn4 cmmnty trtm2: trtm3:
treatment2 -0.522
treatment3 -0.522 0.500
treatment4 -0.522 0.500 0.500
cmmntysnklq -0.697 0.364 0.364 0.364
trtmnt2:cmm 0.362 -0.693 -0.347 -0.347 -0.524
trtmnt3:cmm 0.362 -0.347 -0.693 -0.347 -0.524 0.503
trtmnt4:cmm 0.362 -0.347 -0.347 -0.693 -0.524 0.503 0.503
---<---------------cut here---------------end---------------->---
However, I don't understand how (1 | id) accounts for treatment being
nested within community. Maybe it's time for me to re-read some more
examples from "Mixed-effects models in S and S-plus". Thanks.
Cheers,
--
Seb
More information about the R-sig-mixed-models
mailing list