Menggunakan lmer untuk model efek-campuran linier yang diulang-ulang

41

EDIT 2: Awalnya saya pikir saya perlu menjalankan ANOVA dua faktor dengan tindakan berulang pada satu faktor, tapi saya sekarang berpikir model efek campuran linear akan bekerja lebih baik untuk data saya. Saya pikir saya hampir tahu apa yang perlu terjadi, tetapi saya masih bingung dengan beberapa poin.

Eksperimen yang saya perlu analisis terlihat seperti ini:

  • Subjek ditugaskan ke salah satu dari beberapa kelompok perlakuan
  • Pengukuran setiap mata pelajaran dilakukan dalam beberapa hari
  • Begitu:
    • Subjek bersarang dalam perawatan
    • Perawatan disilangkan dengan hari

(masing-masing subjek hanya ditugaskan untuk satu perawatan, dan pengukuran dilakukan pada setiap subjek setiap hari)

Dataset saya berisi informasi berikut:

  • Subjek = faktor pemblokiran (faktor acak)
  • Hari = dalam subjek atau faktor ukuran berulang (faktor tetap)
  • Pengobatan = antara faktor subjek (faktor tetap)
  • Obs = variabel yang diukur (tergantung)

UPDATE OK, jadi saya pergi dan berbicara dengan ahli statistik, tapi dia pengguna SAS. Menurutnya model itu harus:

Perawatan + Hari + Subjek (Perawatan) + Hari * Subjek (Pengobatan)

Jelas notasinya berbeda dengan sintaks R, tetapi model ini seharusnya menjelaskan:

  • Pengobatan (tetap)
  • Hari (tetap)
  • interaksi pengobatan * hari
  • Subjek bersarang dalam Pengobatan (acak)
  • Hari dilintasi dengan "Subjek dalam Pengobatan" (acak)

Jadi, apakah ini sintaks yang benar untuk digunakan?

m4 <- lmer(Obs~Treatment*Day + (1+Treatment/Subject) + (1+Day*Treatment/Subject), mydata)

Saya khususnya prihatin tentang apakah Hari dilintasi dengan bagian "Subjek dalam Pengobatan" benar. Adakah yang akrab dengan SAS, atau yakin bahwa mereka memahami apa yang terjadi dalam modelnya, dapat mengomentari apakah upaya sedih saya pada sintaksis R cocok?

Berikut ini adalah upaya saya sebelumnya dalam membangun model dan menulis sintaksis (dibahas dalam jawaban & komentar):

m1 <- lmer(Obs ~ Treatment * Day + (1 | Subject), mydata)

Bagaimana cara saya menangani fakta bahwa subjek bersarang dalam perawatan? Apa bedanya dengan m1:

m2 <- lmer(Obs ~ Treatment * Day + (Treatment|Subject), mydata)
m3 <- lmer(Obs ~ Treatment * Day + (Treatment:Subject), mydata)

dan apakah m2dan m3setara (dan jika tidak, mengapa)?

Juga, apakah saya harus menggunakan nlme daripada lme4 jika saya ingin menentukan struktur korelasinya (seperti correlation = corAR1)? Menurut Tindakan Berulang , untuk analisis tindakan berulang dengan tindakan berulang pada satu faktor, struktur kovarians (sifat korelasi antara pengukuran subjek yang sama) adalah penting.

Ketika saya mencoba melakukan ANOVA tindakan berulang, saya memutuskan untuk menggunakan Tipe II SS; apakah ini masih relevan, dan jika demikian, bagaimana saya menentukannya?

Berikut ini contoh tampilan data:

mydata <- data.frame(
  Subject  = c(13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 29, 30, 31, 32, 33, 
               34, 35, 36, 37, 38, 39, 40, 62, 63, 64, 65, 13, 14, 15, 16, 17, 18, 
               19, 20, 21, 22, 23, 24, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 
               40, 62, 63, 64, 65, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 
               29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 62, 63, 64, 65), 
  Day       = c(rep(c("Day1", "Day3", "Day6"), each=28)), 
  Treatment = c(rep(c("B", "A", "C", "B", "C", "A", "A", "B", "A", "C", "B", "C", 
                      "A", "A", "B", "A", "C", "B", "C", "A", "A"), each = 4)), 
  Obs       = c(6.472687, 7.017110, 6.200715, 6.613928, 6.829968, 7.387583, 7.367293, 
                8.018853, 7.527408, 6.746739, 7.296910, 6.983360, 6.816621, 6.571689, 
                5.911261, 6.954988, 7.624122, 7.669865, 7.676225, 7.263593, 7.704737, 
                7.328716, 7.295610, 5.964180, 6.880814, 6.926342, 6.926342, 7.562293, 
                6.677607, 7.023526, 6.441864, 7.020875, 7.478931, 7.495336, 7.427709, 
                7.633020, 7.382091, 7.359731, 7.285889, 7.496863, 6.632403, 6.171196, 
                6.306012, 7.253833, 7.594852, 6.915225, 7.220147, 7.298227, 7.573612, 
                7.366550, 7.560513, 7.289078, 7.287802, 7.155336, 7.394452, 7.465383, 
                6.976048, 7.222966, 6.584153, 7.013223, 7.569905, 7.459185, 7.504068, 
                7.801867, 7.598728, 7.475841, 7.511873, 7.518384, 6.618589, 5.854754, 
                6.125749, 6.962720, 7.540600, 7.379861, 7.344189, 7.362815, 7.805802, 
                7.764172, 7.789844, 7.616437, NA, NA, NA, NA))
terfosforilasi
sumber

Jawaban:

18

Saya pikir pendekatan Anda benar. Model m1menentukan intersep yang terpisah untuk setiap subjek. Model m2menambahkan kemiringan terpisah untuk setiap subjek. Kemiringan Anda melewati beberapa hari karena subjek hanya berpartisipasi dalam satu kelompok perlakuan. Jika Anda menulis model m2sebagai berikut, lebih jelas bahwa Anda memodelkan intersep dan kemiringan yang terpisah untuk setiap subjek

m2 <- lmer(Obs ~ Treatment * Day + (1+Day|Subject), mydata)

Ini setara dengan:

m2 <- lmer(Obs ~ Treatment + Day + Treatment:Day + (1+Day|Subject), mydata)

Yakni efek utama dari perawatan, hari dan interaksi antara keduanya.

Saya pikir Anda tidak perlu khawatir tentang bersarang selama Anda tidak mengulangi ID subjek dalam kelompok perlakuan. Model mana yang benar, sangat tergantung pada pertanyaan riset Anda. Adakah alasan untuk meyakini bahwa kemiringan subjek berbeda di samping efek perawatan? Anda dapat menjalankan kedua model dan membandingkannya dengan anova(m1,m2)untuk melihat apakah data mendukung salah satu.

Saya tidak yakin apa yang ingin Anda ungkapkan dengan model m3? Sintaks bersarang menggunakan a /, mis (1|group/subgroup).

Saya tidak berpikir bahwa Anda perlu khawatir tentang autokorelasi dengan sejumlah kecil titik waktu.

pengguna12719
sumber
Ini tidak benar. Pengobatan adalah variabel level-2, itu tidak dapat disarangkan dalam Subjek.
Patrick Coulombe
Tentang autokorelasi dan jumlah titik waktu: Saya hanya menunjukkan tiga dalam contoh data ini, tetapi data asli saya berisi pengamatan pada 8 hari yang berbeda, jadi saya pikir itu mungkin akan menjadi masalah. Ada ide bagaimana memasukkannya?
Difosforilasi
1
Juga, saya sekarang cukup bingung tentang bersarang; Apakah (1 + Perawatan | Subjek) berbeda dengan (1 + Perawatan / Subjek)? Apa artinya "|" maksud, dalam bahasa Inggris biasa? Saya tidak mengerti penjelasan yang saya baca.
Difosforilasi
Hai. Apa yang ada di sini "kemiringan terpisah untuk setiap subjek"? karena subjek adalah variabel faktor, bukan variabel kontinu.
skan
12

Saya tidak merasa cukup nyaman untuk mengomentari masalah kesalahan autokorelasi Anda (atau tentang implementasi yang berbeda di lme4 vs nlme), tetapi saya dapat berbicara dengan yang lain.

Model Anda m1adalah model intersepsi acak, di mana Anda telah menyertakan interaksi lintas level antara Pengobatan dan Hari (efek Hari diizinkan berbeda di antara kelompok-kelompok Perawatan). Untuk memungkinkan perubahan dari waktu ke waktu berbeda di antara peserta (yaitu untuk secara eksplisit memodelkan perbedaan individu dalam perubahan dari waktu ke waktu), Anda juga perlu membiarkan efek Hari menjadi acak . Untuk melakukan ini, Anda akan menentukan:

m2 <- lmer(Obs ~ Day + Treatment + Day:Treatment + (Day | Subject), mydata)

Dalam model ini:

  • Mencegat jika skor yang diprediksi untuk kategori referensi pengobatan pada Hari = 0
  • Koefisien untuk Hari adalah perubahan yang diperkirakan dari waktu ke waktu untuk setiap peningkatan 1 unit dalam hari untuk kategori referensi pengobatan
  • Koefisien untuk dua kode boneka untuk kelompok perlakuan (secara otomatis dibuat oleh R) adalah perbedaan yang diperkirakan antara masing-masing kelompok perlakuan yang tersisa dan kategori referensi pada Hari = 0
  • Koefisien untuk dua istilah interaksi adalah perbedaan dalam pengaruh waktu (Hari) pada skor yang diprediksi antara kategori referensi dan kelompok perlakuan yang tersisa.

Baik intersep dan efek Day pada skor adalah acak (setiap subjek diperbolehkan memiliki skor prediksi yang berbeda pada Hari = 0 dan perubahan linear yang berbeda dari waktu ke waktu). Kovarians antara intersep dan lereng juga dimodelkan (mereka diizinkan untuk kovari).

Seperti yang Anda lihat, interpretasi koefisien untuk dua variabel dummy tergantung pada Hari = 0. Mereka akan memberi tahu Anda jika skor prediksi pada Hari = 0 untuk kategori referensi secara signifikan berbeda dari dua kelompok perlakuan yang tersisa. Karena itu, di mana Anda memutuskan untuk memusatkan variabel Hari Anda adalah penting. Jika Anda berpusat pada Hari 1, maka koefisien memberi tahu Anda apakah skor prediksi untuk kategori referensi pada Hari 1 berbeda secara signifikan dari skor prediksi kedua kelompok yang tersisa. Dengan cara ini, Anda bisa melihat apakah ada perbedaan yang sudah ada sebelumnya antara kelompok . Jika Anda memusatkan pada Hari 3, maka koefisien memberi tahu Anda apakah skor yang diprediksi untuk kategori referensi pada Hari 3secara signifikan berbeda dari skor prediksi dari dua kelompok yang tersisa. Dengan cara ini, Anda bisa melihat apakah ada perbedaan antara kelompok-kelompok di akhir intervensi .

Akhirnya, perhatikan bahwa Subjek tidak bersarang dalam Pengobatan. Tiga perawatan Anda bukan level acak dari populasi level yang ingin Anda generalisasi hasil Anda - melainkan, seperti yang Anda sebutkan, level Anda tetap, dan Anda ingin menggeneralisasi hasil Anda ke level ini saja. (Belum lagi, Anda tidak boleh menggunakan pemodelan bertingkat jika Anda hanya memiliki 3 unit tingkat atas; lihat Maas & Hox, 2005.) Sebagai gantinya, pengobatan adalah pemrediksi tingkat-2, yaitu pemrediksi yang mengambil nilai tunggal melintasi Hari. (level-1 unit) untuk setiap mata pelajaran. Karena itu, ini hanya dimasukkan sebagai prediktor dalam model Anda.

Referensi:
Maas, CJM, & Hox, JJ (2005). Ukuran sampel yang memadai untuk pemodelan multilevel. Metodologi: Jurnal Eropa Metode Penelitian untuk Ilmu Perilaku dan Sosial , 1 , 86-92.

Patrick Coulombe
sumber
1
Ini tidak dapat diperkirakan oleh LER karena jumlah efek acak <= angka dan varians residual mungkin tidak dapat diidentifikasi.
Shuguang
Struktur rumus dalam jawaban sudah benar. Untuk mengganti kesalahan yang disebutkan oleh @Shuguang, Anda harus menambahkan ...,control=lmerControl(check.nobs.vs.nRE="ignore"). lihat tautan ini untuk penjelasan lebih lanjut oleh Ben Bolker.
NiuBiBang
Penjelasan yang bagus. Bisakah Anda jelaskan sedikit lebih banyak mengapa "Subjek tidak bersarang dalam Pengobatan" dan mengapa Anda tidak membuat istilah kesalahan + (Pengobatan | Subjek) dan mengapa tidak hanya (1 | Subjek) atau bahkan (1 | Pengobatan * Hari )?
skan
Secara teknis Anda bisa membuat sarang subjek dalam perawatan, tetapi jika prediktornya adalah yang akan sama tidak peduli berapa kali Anda menjalankan percobaan, itu harus efek tetap (bukan acak). Faktor-faktor yang akan berbeda setiap kali Anda menjalankan percobaan, seperti karakteristik individu dari subjek - misalnya nilai awal mereka atau respons khas mereka terhadap perubahan dalam pengobatan dari waktu ke waktu - adalah efek acak. (1 + Day|Subject)berarti model lereng acak, yang memungkinkan nilai awal masing-masing subjek (Intercept) dan tingkat perubahan pada hasil menjadi berbeda.
llewmills