[R] unbalanced mixed effects models for fully factorial designs
Murray Logan
murray.logan at sci.monash.edu.au
Thu Aug 3 00:25:03 CEST 2006
Does anyone know of a way of dealing with unbalanced mixed effects
(fixed and random factors) for fully factorial designs.
An example of such data is given below;
The response variable is SQRTRECRUITS
SEASON is a random factor
DENSITY is a fixed factor
Thus DENSITY:SEASON is a fixed factor.
Therefore, whereas the effects of SEASON and DENSITY:SEASON should be
tested against the overall residual (error) term, the effect of DENSITY
should be tested against the DENSITY:SEASON interaction.
To complicate matters, the data are unbalanced, and thus Type III SS are
preferable
quinn <-
structure(list(SEASON = structure(as.integer(c(2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 4, 4, 4, 4, 4, 4, 4, 4, 4)), .Label = c("Autumn",
"Spring", "Summer", "Winter"), class = "factor", contrasts = "contr.sum"),
DENSITY = structure(as.integer(c(2, 2, 2, 2, 2, 1, 1, 1,
1, 1, 1, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2,
1, 1, 1, 1, 1, 1, 2, 2, 2, 1, 1, 1, 1, 1, 1)), .Label = c("High",
"Low"), class = "factor"), RECRUITS = as.integer(c(15, 10,
13, 13, 5, 11, 10, 15, 10, 13, 1, 21, 31, 21, 18, 14, 27,
34, 49, 69, 55, 28, 54, 14, 18, 20, 21, 4, 22, 30, 36, 13,
13, 8, 0, 0, 10, 1, 5, 9, 4, 5)), SQRTRECRUITS = c(3.872983,
3.162278, 3.605551, 3.605551, 2.236068, 3.316625, 3.162278,
3.872983, 3.162278, 3.605551, 1, 4.582576, 5.567764, 4.582576,
4.242641, 3.741657, 5.196152, 5.830952, 7, 8.306624, 7.416198,
5.291503, 7.348469, 3.741657, 4.242641, 4.472136, 4.582576,
2, 4.690416, 5.477226, 6, 3.605551, 3.605551, 2.828427, 0,
0, 3.162278, 1, 2.236068, 3, 2, 2.236068), GROUP =
structure(as.integer(c(4,
4, 4, 4, 4, 3, 3, 3, 3, 3, 3, 6, 6, 6, 6, 6, 6, 5, 5, 5,
5, 5, 5, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 8, 8, 8, 7, 7, 7,
7, 7, 7)), .Label = c("AutumnHigh", "AutumnLow", "SpringHigh",
"SpringLow", "SummerHigh", "SummerLow", "WinterHigh", "WinterLow"
), class = "factor")), .Names = c("SEASON", "DENSITY", "RECRUITS",
"SQRTRECRUITS", "GROUP"), row.names = c("1", "2", "3", "4", "5",
"6", "7", "8", "9", "10", "11", "12", "13", "14", "15", "16",
"17", "18", "19", "20", "21", "22", "23", "24", "25", "26", "27",
"28", "29", "30", "31", "32", "33", "34", "35", "36", "37", "38",
"39", "40", "41", "42"), class = "data.frame")
I realise that Anova (car package) calculated Type III SS (given the
correct contrasts), however, this does not permit mixed models.
Conversely, if I was to specify a aov model such as;
summary(aov(SQRTRECRUITS ~ SEASON+DENSITY+Error(DENSITY:SEASON),
data=quinn))
purely to obtain a test for DENSITY (ignoring the test for SEASON), the
SS are Type I.
Although it is possible to calculate out the F-ratio (and p-value)
calculations manually and substitute them into the anova tables, I cant
help think that there must be a better solution.
Is there any expectation that there will be a summary routine that
provides Type II and Type II SS, and or is aov ever likely to
accommodate non-hierarchical mixed models?
Regards
Murray
More information about the R-help
mailing list