Por que essas duas abordagens para aplicar modelos mistos produzem resultados diferentes?

Estou analisando novamente os dados de um colega. Os dados e o código R estão aqui .

É um projeto 2x2x2x2x3 completamente dentro do Ss. Uma das variáveis ​​preditoras,, cueé uma variável de dois níveis que, quando diminuída para uma pontuação de diferença, reflete um valor pertinente à teoria. Anteriormente, ela reduziu cuea uma pontuação de diferença dentro de cada sujeito e condição, depois calculou uma ANOVA, produzindo um MSE que ela poderia usar para comparações planejadas da pontuação de diferença média de cada condição contra zero. Você terá que confiar em mim que ela não estava pescando e realmente tinha uma boa base teórica para realizar todos os 24 testes.

Eu pensei em ver se havia alguma diferença ao usar modelos de efeitos mistos para representar os dados. Como mostrado no código, tomei duas abordagens:

Método 1 - Modele os dados como um projeto 2x2x2x2x3, obtenha amostras a posteriori deste modelo, calcule a cuepontuação da diferença para cada condição em cada amostra, calcule o intervalo de previsão de 95% para a pontuação da diferença da sugestão em cada condição.

Método 2 - Reduza cuepara uma pontuação de diferença dentro de cada sujeito e condição, modele os dados como um projeto 2x2x2x3, obtenha amostras a posteriori desse modelo, calcule o intervalo de previsão de 95% para a pontuação de diferença de sugestão em cada condição.

Parece que o método 1 gera intervalos de previsão mais amplos do que o método 2, com a conseqüência de que, se se usar sobreposição com zero como critério de "significância", apenas 25% das pontuações no cuing serão "significativas" no método 1, enquanto 75% das pontuações no cuing são "significativos" no método 2. Sem dúvida, os padrões de significância obtidos pelo método 2 são mais semelhantes aos resultados originais baseados em ANOVA do que os padrões obtidos pelo método 1.

Alguma idéia do que está acontecendo aqui?

Não é surpreendente ver essa diferença com lmer ou lme. Um modelo simples com uma interceptação aleatória (por exemplo, (1 | id) no seu caso) às vezes pode não conseguir capturar completamente os efeitos aleatórios. Para ver por que isso acontece, deixe-me usar um conjunto de dados muito mais simples que o seu para demonstrar a diferença sutil. Com os dados 'dat' do thread que copio para aqui:

dat <- structure(list(sex = structure(c(1L, 2L, 1L, 2L, 1L, 2L, 1L,
2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L,
2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L), .Label = c("f",
"m"), class = "factor"), prevalence = c(0, 0.375, 0.133333333333333,
0.176470588235294, 0.1875, 0, 0, 1, 1, 0.5, 0.6, 0.333333333333333,
0.5, 0, 0.333333333333333, 0, 0.5, 0, 0.625, 0.333333333333333,
0.5, 0, 0.333333333333333, 0.153846153846154, 0.222222222222222,
0.5, 1, 0.5, 0, 0.277777777777778, 0.125, 0, 0, 0.428571428571429,
0.451612903225806, 0.362068965517241), tripsite = structure(c(1L,
1L, 4L, 4L, 14L, 14L, 5L, 5L, 8L, 8L, 15L, 15L, 6L, 6L, 9L, 9L,
11L, 11L, 16L, 16L, 2L, 2L, 7L, 7L, 10L, 10L, 13L, 13L, 17L,
17L, 3L, 3L, 12L, 12L, 18L, 18L), .Label = c("1.2", "4.2", "5.2",
"1.3", "2.3", "3.3", "4.3", "2.4", "3.4", "4.4", "3.5", "5.5",
"4.6", "1.9", "2.9", "3.9", "4.9", "5.9"), class = "factor")), .Names =
c("sex","prevalence", "tripsite"), row.names = c(1L, 2L, 3L, 4L, 9L,
10L, 11L, 12L, 13L, 14L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L,
27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L, 38L, 39L, 40L,
41L, 42L, 43L, 45L, 46L), class = "data.frame")

um teste t emparelhado (ou um caso especial de ANOVA unidirecional / medidas repetidas) seria como o método 2:

t0 <- with(dat,t.test(prevalence[sex=="f"],prevalence[sex=="m"],paired=TRUE,var.equal=TRUE))
(fstat0 <- t0$statistic^2)         #0.789627

Sua versão lme correspondente ao seu método 1 seria:

a1 <- anova(lme(prevalence~sex,random=~1|tripsite,data=dat,method="REML"))
(fstat1 <- a1[["F-value"]][2])   # 0.8056624

A mesma coisa para a contraparte mais antiga:

a2 <- anova(lmer(prevalence~sex+(1|tripsite), data=dat))
(fstat2 <- a2[["F value"]][2])  # 0.8056624

Embora a diferença neste exemplo simples seja pequena, mas mostra que o teste t emparelhado tem uma suposição muito mais forte sobre os dois níveis ("f" e "m") do fator ("sexo"), que os dois níveis estão correlacionados e essa suposição está ausente no modelo lme / lmer acima. Essa diferença de suposição também existe entre os dois métodos no seu caso.

Para reconciliar a diferença, podemos continuar modelando 'dat' com uma inclinação aleatória (ou matriz simétrica ou mesmo simetria composta) em lme / lmer:

a3 <- anova(lme(prevalence~sex,random=~sex-1|tripsite,data=dat,method="REML"))
(fstat3 <- a3[["F-value"]][2]) # 0.789627

a31 <- anova(lme(prevalence~sex,random=list(tripsite=pdCompSymm(~sex-1)),data=dat,method="REML")))
(fstat31 <- a31[["F-value"]][2]) # 0.789627

a4 <- anova(lmer(prevalence~sex+(sex-1|tripsite), data=dat))
(fstat4 <- a4[["F value"]][2]) # 0.789627

No entanto, com vários fatores no seu caso, várias pistas aleatórias (ou outras especificações da estrutura de efeitos aleatórios) podem ficar difíceis de manejar com lme / lmer, se não impossível.