Mengapa 2 pendekatan ini untuk menerapkan model campuran menghasilkan hasil yang berbeda?

Saya menganalisis ulang data kolega. Data dan kode R ada di sini .

Ini adalah desain 2x2x2x2x3 sepenuhnya dalam-Ss. Salah satu variabel prediktor cue,, adalah variabel dua tingkat yang ketika diciutkan ke skor perbedaan mencerminkan nilai yang berhubungan dengan teori. Dia sebelumnya jatuh cueke skor perbedaan dalam setiap subjek dan kondisi, kemudian menghitung ANOVA, menghasilkan MSE yang kemudian dapat ia gunakan untuk perbandingan yang direncanakan dari skor perbedaan rata-rata setiap kondisi terhadap nol. Anda harus percaya kepada saya bahwa dia tidak memancing dan memang memiliki dasar teori yang baik untuk melakukan semua 24 tes.

Saya pikir saya akan melihat apakah ada perbedaan ketika menggunakan model efek campuran untuk mewakili data. Seperti yang ditunjukkan dalam kode, saya mengambil dua pendekatan:

Metode 1 - Memodelkan data sebagai desain 2x2x2x2x3, memperoleh sampel posteriori dari model ini, menghitung cueskor perbedaan untuk setiap kondisi dalam setiap sampel, menghitung interval prediksi 95% untuk skor perbedaan isyarat dalam setiap kondisi.

Metode 2 - Perkecil cueskor selisih dalam setiap subjek dan kondisi, modelkan data sebagai desain 2x2x2x3, dapatkan sampel posteriori dari model ini, hitung interval prediksi 95% untuk skor selisih isyarat dalam setiap kondisi.

Tampaknya metode 1 menghasilkan interval prediksi yang lebih luas daripada metode 2 dengan konsekuensi bahwa jika seseorang menggunakan tumpang tindih dengan nol sebagai kriteria untuk "signifikansi", hanya 25% dari skor cuing yang "signifikan" di bawah metode 1 sementara 75% dari skor cuing adalah "signifikan" dengan metode 2. Secara khusus, pola signifikansi yang diperoleh dengan metode 2 lebih mirip dengan hasil berbasis ANOVA asli daripada pola yang diperoleh dengan metode 1.

Adakah yang tahu apa yang terjadi di sini?

Tidaklah mengejutkan melihat perbedaan dengan lmer atau lme. Model sederhana dengan intersep acak (misalnya, (1 | id) dalam kasus Anda) kadang-kadang mungkin gagal untuk sepenuhnya menangkap efek acak. Untuk melihat mengapa hal ini terjadi, izinkan saya menggunakan dataset yang jauh lebih sederhana daripada milik Anda untuk menunjukkan perbedaan yang halus. Dengan data 'dat' dari utas yang saya salin ke sini:

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")

uji-t berpasangan (atau kasus khusus ANOVA satu arah dalam subjek / tindakan berulang) akan seperti Metode 2 Anda:

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

Versi saya yang sesuai dengan Metode 1 Anda adalah:

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

Hal yang sama untuk mitra lmer:

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

Meskipun perbedaan dengan contoh sederhana ini kecil, tetapi ini menunjukkan bahwa uji-t berpasangan memiliki asumsi yang jauh lebih kuat tentang dua tingkat ("f" dan "m") dari faktor ("jenis kelamin"), bahwa kedua tingkat tersebut berkorelasi, dan asumsi tersebut tidak ada dalam model lme / lmer diatas. Perbedaan asumsi semacam itu juga ada di antara kedua metode dalam kasus Anda.

Untuk merekonsiliasi perbedaan, kita dapat melanjutkan pemodelan 'dat' dengan kemiringan acak (atau matriks simetris atau bahkan gabungan simetri) di 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

Namun, dengan beberapa faktor dalam kasus Anda, beberapa lereng acak (atau spesifikasi struktur efek-acak lainnya) dapat menjadi sulit dengan lme / lmer jika bukan tidak mungkin.