Dear S users,
I ask a while ago a question about formula for lme (and glmmPQL) that I copy
below. Since nobody answered (is it too trivial?), I try a second time with
more details. I also give below the entire dataframe for possible tests. My
question is not an academic one, but I would like to use this type of models
for non balanced design and for non gaussian (e.g. binary) data, where aov does
not apply anymore.
Here are some of the model I tried that do not give the expected results or
give and error message:
My guess was this command:
> anova (lme(fixed = Reponse ~ A*B*C, random=~ 1|Sujet/(A*B), data =
> orienta,method="ML"))
Problem in getGroups.data.frame(dataMix, groups): Invalid formula for groups
> anova (lme(fixed = Reponse ~ A*B*C, random=~ 1|Sujet/(B+A), data =
> orienta,method="ML"))
Problem in getGroups.data.frame(dataMix, groups): Invalid formula for groups
> anova (lme(fixed = Reponse ~ A*B*C, random=~ 1|Sujet, data =
> orienta,method="ML"))
## gives the same results as
> summary (aov(formula = Reponse ~ A*B*C + Error(Sujet), data = orienta))
## which works but does not contain all the interactions wanted
> anova (lme(fixed = Reponse ~ A*B*C, random=~ 1|Sujet/(B), data =
> orienta,method="ML"))
## works also but does not contain all the interactions wanted
Olivier
At 09:51 13/05/02 +0200, I wrote:
Dear S users,
I'm afraid this is a trivial question, but I couldn't find the answer in
Pinheiro and Bates' book (which does not mean it is not there) and a long try
and error session was fruitless.
I would like to find the equivalent for lme and glmmPQL (with binomial link) of
the model that writes
aov(formula = Reponse ~ A*B*C + Error(Sujet/(A*B))
In the psychologists's dialect of statistics, "Sujet" is the subject, C is a
between- and A and B are within-subject factors.
In the engineer's dialect of statistics, Sujet is a random factor, A, B, and C
are fixed; Sujet is nested within C, and all two-way interactions are present,
as well as possible three way interactions ("A:B:C" and "A:B %in% Sujet", I
believe).
Below is the S+ output (equal to SPSS results). Note that if I have only one
experiment per condition (which is typically the case in psychoogy and which is
the case here), "A:B %in% Sujet" is counfounded with the error, which is a
problem when lme computes confidence interval. For aov, however, the above
model gives the same results as
aov(formula = Reponse ~ A*B*C + Error(Sujet/(A+B))
Thank you,
Olivier
"orienta" <-
structure(.Data = list("value" = structure(.Data = c(1, 1, 1, 1, 2, 2, 2, 2, 3,
3, 3, 3, 4, 4, 4, 4, 5, 5, 5, 5, 6, 6, 6, 6,
7, 7, 7, 7, 8, 8, 8, 8, 9, 9, 9, 9, 10, 10, 10, 10, 11, 11, 11, 11,
12, 12, 12, 12, 13, 13, 13, 13, 14, 14, 14, 14, 15, 15, 15, 15, 16,
16, 16, 16, 17, 17, 17, 17, 18, 18, 18, 18, 19, 19, 19, 19, 20, 20,
20, 20)
, levels = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11", "12",
"13", "14",
"15", "16", "17", "18", "19", "20")
, class = "factor"
)
, c(151.09999999999999, 172., 176., 172.90000000000001, 38.899999999999999,
63., 94., 147.09999999999999, 57.100000000000001, 146., 144.,
112.90000000000001, 30.899999999999999, 50., 46., 79.099999999999994,
48.899999999999999, 62., 100., 121.09999999999999, 1.1000000000000001,
24., 120., 146.90000000000001, 63.100000000000001, 126., 139.,
172.90000000000001, 11.1, 14., 20., 10.9, 119.09999999999999, 124.,
87., 37.899999999999999, 1.1000000000000001, 42., 40.,
82.099999999999994, 9.0999999999999996, 66., 128., 148.90000000000001,
26.899999999999999, 68., 110., 177.09999999999999, 3.1000000000000001,
4., 56., 30.899999999999999, 15.1, 48., 14., 35.100000000000001,
6.9000000000000004, 38., 64., 83.099999999999994, 67.099999999999994,
142., 148., 144.90000000000001, 80.900000000000006, 150., 168.,
162.90000000000001, 175.09999999999999, 124., 122., 135.09999999999999,
15.9, 48., 84., 158.09999999999999, 20.899999999999999, 30., 6.,
16.100000000000001)
, structure(.Data = c(1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2,
2, 1, 1, 2, 2,
1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1,
2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2,
1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2)
, levels = c("a1", "a2")
, class = "factor"
)
, structure(.Data = c(1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1,
2, 1, 2, 1, 2,
1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2,
1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2,
1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2)
, levels = c("b1", "b2")
, class = "factor"
)
, structure(.Data = c(2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1)
, levels = c("F", "H")
, class = "factor"
)
)
, names = c("Sujet", "Reponse", "A", "B", "C")
, 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", "43", "44", "45", "46", "47",
"48", "49", "50", "51", "52", "53", "54", "55", "56", "57", "58",
"59", "60", "61", "62", "63", "64", "65", "66", "67", "68", "69",
"70", "71", "72", "73", "74", "75", "76", "77", "78", "79", "80")
, class = "data.frame"
)
summary (aov(formula = Reponse ~ A*B*C + Error(Sujet/(A*B)), data = orienta))
Error: Sujet
Df Sum of Sq Mean Sq F Value Pr(F)
C 1 871.2 871.200 0.1034248 0.7514624
Residuals 18 151623.2 8423.513
Error: A %in% Sujet
Df Sum of Sq Mean Sq F Value Pr(F)
A 1 30326.47 30326.47 14.13652 0.0014341
A:C 1 276.77 276.77 0.12901 0.7236356
Residuals 18 38614.62 2145.26
Error: B %in% Sujet
Df Sum of Sq Mean Sq F Value Pr(F)
B 1 10296.72 10296.72 19.04355 0.0003742
B:C 1 62.66 62.66 0.11588 0.7374823
Residuals 18 9732.48 540.69
Error: A:B %in% Sujet
Df Sum of Sq Mean Sq F Value Pr(F)
A:B 1 1033.922 1033.922 2.069967 0.1673868
A:B:C 1 76.050 76.050 0.152256 0.7009675
Residuals 18 8990.768 499.487
|