I asked the following question 10 days ago, and received no answer. If you
don't have time for a lengthy answer, a pointer to a reference would be
appreciated.
Olivier
Hi,
I'm not exactly sure how to call the design that we use, but it may contain the
keywords fractional design, repeated measure, BIB, MLM, and maybe more! Here is
the setting: We have three types of route (a1, a2, a3) and three conditions
(b1, b2, b3). Since each subject (S) can only do each route once, each subject
was asked to do only 3 experiments:
for 1/6 of subjects, it was the 3 combinations (a1,b1), (a2,b2), (a3,b3)
for 1/6 of subjects, it was the 3 combinations (a1,b1), (a2,b3), (a3,b2)
for 1/6 of subjects, it was the 3 combinations (a1,b2), (a2,b1), (a3,b3)
for 1/6 of subjects, it was the 3 combinations (a1,b2), (a2,b3), (a3,b1)
for 1/6 of subjects, it was the 3 combinations (a1,b3), (a2,b1), (a3,b2)
for 1/6 of subjects, it was the 3 combinations (a1,b3), (a2,b2), (a3,b1)
I guess the formulas for aov and lme are (on simulated data)
> summary (aov (resp ~ A+B+Error(S/(A+B)), data=esfr5))
> anova (lme (resp ~ A+B, random=~1|S, data=esfr5), type="marginal")
# See Output and Data below
with a * instead of the + if one includes the interactions. I have the
following questions:
1) Why is the error in aov divided only in two pieces. If the design were
complete (each subject having done 9 experiments), it would be divided in four?
Moreover, why is the 2nd Error term called "A %in% S"?
2) Is the interaction A:B estimable or is it confounded with something? As you
noticed, it appears twice in the tests. What do these two tests mean? This is a
simulated example, but what if they are very different?
3) Up to now, with _real_ data, lme always find singularity ('Problem in
.C("mixed_EM",: Singularity in backsolve, while calling subroutine mixed_EM' ).
Are they some parameters in lmeControl that can be tuned to solve this ? (aov
has no problem with these data, see output of one example of real data and the
data below)
4) Are the adjusted models different between lme and aov (especially when
including interaction) ?
5) Can we analyse these data with SPSS ?
Olivier
#Output of simulated data
> summary(aov(resp~A+B+Error(S/(A+B)),data=esfr5))
Error: S
Df Sum of Sq Mean Sq F Value Pr(F)
Residuals 29 3.114277 0.1073889
Error: A %in% S
Df Sum of Sq Mean Sq F Value Pr(F)
A 2 0.071461 0.0357305 0.458006 0.6348918
B 2 0.206622 0.1033109 1.324277 0.2741954
Residuals 56 4.368734 0.0780131
> summary(aov(resp~A*B+Error(S/(A*B)),data=esfr5))
Error: S
Df Sum of Sq Mean Sq F Value Pr(F)
A:B 4 0.491633 0.1229084 1.171608 0.347148
Residuals 25 2.622643 0.1049057
Error: A %in% S
Df Sum of Sq Mean Sq F Value Pr(F)
A 2 0.071461 0.0357305 0.463168 0.6318575
B 2 0.206622 0.1033109 1.339203 0.2709411
A:B 4 0.357267 0.0893167 1.157798 0.3401417
Residuals 52 4.011467 0.0771436
> anova(lme(resp~A+B, random=~1|S, data=esfr5), type="marginal")
numDF denDF F-value p-value
(Intercept) 1 56 190.1546 <.0001
A 2 56 0.4580 0.6349
B 2 56 1.3243 0.2742
> anova(lme(resp~A*B, random=~1|S, data=esfr5), type="marginal")
numDF denDF F-value p-value
(Intercept) 1 52 190.4814 <.0001
A 2 52 0.4597 0.6340
B 2 52 1.3292 0.2735
A:B 4 52 1.0644 0.3835
# Output of real data:
> summary(aov(Reponse~Lien+Repetition+Error(Sujet/(Lien+Repetition)),data=SPenSLR1))
Error: Sujet
Df Sum of Sq Mean Sq F Value Pr(F)
Residuals 23 4984396 216712.9
Error: Lien %in% Sujet
Df Sum of Sq Mean Sq F Value Pr(F)
Lien 2 62062.6 31031.30 2.917880 0.0645736
Repetition 2 35207.5 17603.76 1.655285 0.2027117
Residuals 44 467934.7 10634.88
> anova(lme( fixed=Reponse~Lien+Repetition, data=SPenSLR1, random= ~1| Sujet),
> type="marginal")
Problem in .C("mixed_EM",: Singularity in backsolve, while calling subroutine
mixed_EM
Use traceback() to see the call stack
# Simulated and real data:
"esfr5" <-
structure(.Data = list(c("s1", "s1", "s1", "s2", "s2", "s2", "s3", "s3", "s3",
"s4", "s4", "s4", "s5", "s5", "s5", "s6", "s6", "s6", "s7", "s7", "s7",
"s8", "s8", "s8", "s9", "s9", "s9", "s10", "s10", "s10", "s11", "s11",
"s11", "s12", "s12", "s12", "s13", "s13", "s13",
"s14", "s14", "s14", "s15", "s15", "s15", "s16", "s16", "s16", "s17",
"s17", "s17", "s18", "s18", "s18", "s19", "s19",
"s19", "s20", "s20", "s20", "s21", "s21", "s21", "s22", "s22", "s22",
"s23", "s23", "s23", "s24", "s24", "s24", "s25",
"s25", "s25", "s26", "s26", "s26", "s27", "s27", "s27", "s28", "s28",
"s28", "s29", "s29", "s29", "s30", "s30", "s30")
, "value" = c("a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2",
"a3", "a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3",
"a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3",
"a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2",
"a3", "a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2",
"a3", "a1", "a2", "a3", "a1", "a2", "a3", "a1",
"a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3", "a1",
"a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3",
"a1", "a2", "a3", "a1", "a2", "a3", "a1", "a2", "a3")
, "value" = structure(.Data = c(1, 2, 3, 1, 3, 2, 2, 1, 3, 2, 3, 1, 3, 1, 2, 3,
2, 1, 1, 2, 3, 1, 3, 2, 2, 1, 3, 2, 3, 1, 3, 1, 2, 3, 2, 1, 1, 2, 3, 1, 3,
2, 2, 1, 3, 2, 3, 1, 3, 1, 2, 3, 2, 1, 1, 2, 3, 1, 3, 2, 2, 1, 3, 2, 3,
1, 3, 1, 2, 3, 2, 1, 1, 2, 3, 1, 3, 2, 2, 1,
3, 2, 3, 1, 3, 1, 2, 3, 2, 1)
, levels = c("b1", "b2", "b3")
, class = "factor"
)
, c(0.024646461009979248, 0.82596908044070005, 0.39935583062469959,
0.42794170277193189, 0.79330063052475452, 0.41740808542817831,
0.42689276579767466, 0.97277549048885703, 0.86306705698370934,
0.20230660401284695, 0.33371580624952912,
0.61999470554292202, 0.74001671373844147, 0.07387995719909668,
0.30246197292581201, 0.35847964184358716,
0.46099079307168722, 0.16043389076367021, 0.58081981865689158,
0.22667973302304745, 0.16983740823343396,
0.33375463215634227, 0.41230997582897544, 0.35489231906831264,
0.1146096084266901, 0.66199403256177902,
0.086105858441442251, 0.92259258404374123, 0.033672890160232782,
0.97718893829733133, 0.54872243478894234,
0.063163440208882093, 0.18481955723837018, 0.69088406767696142,
0.00048982584848999977, 0.9480914156883955,
0.18231043219566345, 0.006615795660763979, 0.67624678136780858,
0.31703903619199991, 0.38895514793694019,
0.21916917618364096, 0.80299417721107602, 0.61142651597037911,
0.60505218058824539, 0.75811595609411597,
0.59731168439611793, 0.90336279338225722, 0.66943846317008138,
0.86255600396543741, 0.58649370679631829,
0.64013873925432563, 0.35643625492230058, 0.1302819293923676,
0.76756644900888205, 0.39251976786181331,
0.19385471474379301, 0.01424878928810358, 0.21007986413314939,
0.059990861918777227, 0.7896508164703846,
0.25749274902045727, 0.77087316988036036, 0.50042055826634169,
0.96990229655057192, 0.95752570079639554,
0.68353205313906074, 0.7501282668672502, 0.41614685719832778,
0.24187999684363604, 0.68989087454974651,
0.80914638098329306, 0.44907701155170798, 0.83931859489530325,
0.93146510748192668, 0.44827951770275831,
0.66668012505397201, 0.82564923400059342, 0.090691552963107824,
0.61973605211824179, 0.74060451937839389,
0.033790021203458309, 0.62395105883479118, 0.56522526312619448,
0.32841544877737761, 0.67579050175845623,
0.029812175780534744, 0.0094875660724937916, 0.4210558719933033,
0.067989643197506666)
)
, 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", "81", "82",
"83", "84", "85", "86", "87", "88", "89", "90")
, names = c("S", "A", "B", "resp")
, class = "data.frame"
)
"SPenSLR1" <-
structure(.Data = list(c(1164.5, 968., 931.25, 1263.25, 1091., 1189.5, 802.,
655.5, 689.75, 1067.5, 1178.75, 1370.75, 909.5, 878.5, 1070., 649.25,
594.25, 580.75, 769., 776.75, 688., 1489.75, 1361.25, 1428.25, 718.5,
608.25, 612.5, 892.25, 1509.75, 1372.75, 1477.5,
1284., 1425.5, 1014.25, 981.5, 1218.5, 798.75, 840.25, 1022.75, 966.75,
1049., 1218., 1212.5, 1027.5, 1149.5, 798.75,
938.25, 919., 1034.5, 947.5, 1028.5, 997., 900.25, 1120.25, 611.5,
688., 621.75, 680.25, 722., 685.5, 822.75, 914.,
973.75, 912.75, 969., 992.75, 1602.25, 1596.5, 1716.25, 1200.25, 1281.,
1275.)
, structure(.Data = c(7, 7, 7, 8, 8, 8, 9, 9, 9, 10, 10, 10, 11, 11, 11, 12,
12, 12, 13, 13, 13, 14, 14, 14, 1, 1, 1, 15, 15, 15, 16, 16, 16, 17,
17, 17, 18, 18, 18, 19, 19, 19, 20, 20, 20, 21, 21, 21, 22, 22, 22, 23,
23, 23, 2, 2, 2, 24, 24, 24, 25, 25, 25, 26,
26, 26, 27, 27, 27, 28, 28, 28)
, levels = c("2", "3", "4", "6", "7", "9", "10", "11", "12", "13", "15", "17",
"18", "19", "20", "21", "22", "23", "24", "25", "26", "27",
"29", "31", "33", "34", "35", "36", "37", "38", "39", "40", "41", "42",
"43", "44", "45", "46", "47", "49", "50", "51",
"52", "53", "55", "59", "60", "61", "62", "64", "65")
, class = "factor"
)
, structure(.Data = c(1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1,
2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2,
3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2,
3, 1, 2, 3, 1, 2, 3)
, levels = c("1", "2", "3", "4")
, class = "factor"
)
, structure(.Data = c(1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 3, 2, 1, 3, 2, 1,
3, 2, 1, 3, 2, 2, 1, 3, 2, 1, 3, 2, 1, 3, 2, 1, 3, 2, 3, 1, 2, 3,
1, 2, 3, 1, 2, 3, 1, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 1, 2, 3, 2, 1, 3, 2,
1, 3, 2, 1, 3, 2, 1)
, levels = c("1", "2", "3", "4")
, class = "factor"
)
)
, names = c("Reponse", "Sujet", "Lien", "Repetition")
, class = "data.frame"
, row.names = c("1", "6", "11", "17", "22", "27", "33", "38", "43", "49", "54",
"59", "65", "71", "74", "81", "87", "90", "97", "103", "106",
"113", "119", "122", "130", "133", "139", "146", "149", "155", "162",
"165", "171", "178", "181", "187", "194", "199",
"201", "210", "215", "217", "226", "231", "233", "242", "247", "249",
"259", "261", "266", "275", "277", "282", "291",
"293", "298", "307", "309", "314", "323", "326", "329", "339", "342",
"345", "355", "358", "361", "371", "374", "377")
)
|