Membandingkan model efek campuran dan efek tetap (menguji signifikansi efek acak)

10

Diberikan tiga variabel, ydan x, yang positif kontinu, dan z, yang kategorikal, saya memiliki dua model kandidat yang diberikan oleh:

fit.me <- lmer( y ~ 1 + x + ( 1 + x | factor(z) ) )

dan

fit.fe <- lm( y ~ 1 + x )

Saya berharap untuk membandingkan model ini untuk menentukan model mana yang lebih tepat. Tampaknya bagi saya bahwa dalam beberapa hal fit.febersarang di dalam fit.me. Biasanya, ketika skenario umum ini berlaku, tes chi-squared dapat dilakukan. Di R, kita bisa melakukan tes ini dengan perintah berikut,

anova(fit.fe,fit.me)

Ketika kedua model berisi efek acak (dihasilkan oleh lmerdari lme4paket), anova()perintah tersebut berfungsi dengan baik. Karena parameter batas, biasanya disarankan untuk menguji statistik Chi-Square yang dihasilkan melalui simulasi, namun demikian, kita masih dapat menggunakan statistik dalam prosedur simulasi.

Ketika kedua model hanya berisi efek tetap, pendekatan ini --- dan, anova()perintah terkait --- berfungsi dengan baik.

Namun, ketika satu model berisi efek acak dan model yang dikurangi hanya berisi efek tetap, seperti dalam skenario di atas, anova()perintah tidak berfungsi.

Lebih khusus lagi, saya mendapatkan kesalahan berikut:

 > anova(fit.fe, fit.me)
 Error: $ operator not defined for this S4 class

Apakah ada yang salah dengan menggunakan pendekatan Chi-Square dari atas (dengan simulasi)? Atau apakah ini hanya masalah karena anova()tidak tahu bagaimana berurusan dengan model linier yang dihasilkan oleh fungsi yang berbeda?

Dengan kata lain, apakah akan sesuai untuk secara manual menghasilkan statistik Chi-Square yang berasal dari model? Jika demikian, apa derajat kebebasan yang tepat untuk membandingkan model-model ini? Menurut saya:

F=((SSEreducedSSEfull)/(pk))((SSEfull)/(np1))Fpk,np1

Kami memperkirakan dua parameter dalam model efek tetap (kemiringan dan intersepsi) dan dua parameter lebih (parameter varians untuk kemiringan acak dan intersep acak) dalam model efek campuran. Biasanya, parameter intersep tidak dihitung dalam perhitungan derajat kebebasan, sehingga menyiratkan bahwa dan ; karena mengatakan bahwa saya tidak yakin apakah parameter varians untuk parameter efek-acak harus dimasukkan dalam derajat kebebasan komputasi; estimasi varians untuk parameter efek tetap tidak dipertimbangkan , tetapi saya percaya hal itu karena estimasi parameter untuk efek tetap diasumsikan sebagai konstanta yang tidak diketahui sementara mereka dianggap sebagai variabel acak yang tidak diketahuik=1p=k+2=3untuk efek campuran. Saya akan sangat menghargai bantuan untuk masalah ini.

Akhirnya, apakah ada yang punya solusi yang lebih tepat ( Rberbasis) untuk membandingkan model ini?

pengguna9171
sumber
4
Jika Anda mengganti lm()dengan gls()dari nlmepaket, dan lmer()dengan lme()(lagi dari nlmepaket), semuanya akan berfungsi dengan baik. Tetapi perhatikan bahwa Anda akan mendapatkan tes konservatif (nilai p terlalu besar ), karena parameter untuk model yang lebih sederhana adalah pada batas ruang parameter. Dan benar-benar pilihan apakah akan memasukkan efek acak harus didasarkan pada teori (misalnya, rencana pengambilan sampel), bukan pada uji statistik.
Karl Ove Hufthammer
1
Apa yang ingin Anda lakukan dengan model? Satu model mungkin lebih baik untuk beberapa tujuan dan model lainnya lebih baik untuk tujuan lain. Semua model salah, jadi pertanyaannya bukan model mana yang benar, tetapi mana yang lebih berguna untuk masalah khusus Anda.
Kodiologist
1
@Kodiologist Pada dasarnya, saya ingin memastikan bahwa estimasi parameter untuk efek tetap dapat diandalkan. Kesalahan standar yang mungkin tidak dapat diandalkan jika pengamatan diasumsikan independen. Selain itu, akan lebih baik untuk membuat beberapa pernyataan tentang bagaimana variabel efek acak, tapi saya kira ini tidak begitu penting.
user9171
2
@ user9171 Cara yang baik untuk memeriksa stabilitas (keandalan) dalam estimasi parameter model adalah dengan menggunakan bootstrap. Distribusi grafik bootstrap untuk setiap parameter yang dibagi dua model, dengan satu grafik per parameter dan model. Distribusi yang lebih ketat menyiratkan stabilitas yang lebih tinggi. Anda mungkin akan menemukan bahwa model yang lebih sederhana menghasilkan estimasi yang lebih stabil, karena lebih sedikit parameter memungkinkan untuk estimasi yang lebih tepat dari setiap parameter.
Kodiologist

Jawaban:

6

Secara teknis, Anda bisa membuatnya bekerja hanya dengan mengganti urutan parameter:

> anova(fit.me, fit.fe) 

Akan bekerja dengan baik. Jika Anda melewatkan objek yang dihasilkan lmerterlebih dahulu, anova.merModakan dipanggil bukan anova.lm(yang tidak tahu cara menangani lmerobjek). Lihat:

?anova.merMod

Meskipun, memilih model campuran atau model tetap adalah pilihan pemodelan yang perlu mempertimbangkan desain eksperimental, bukan masalah pemilihan model. Lihat @ BenBolker https://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#testing-significance-of-random-effects untuk detail lebih lanjut:

Pertimbangkan untuk tidak menguji signifikansi efek acak.

witek
sumber
+1. Saya mengambil kebebasan untuk memasukkan tautan ke FAQ @ BenBolker yang berisi beberapa diskusi dan referensi lebih lanjut.
amoeba