Menambahkan efek acak memengaruhi estimasi koefisien

Saya selalu diajarkan bahwa efek acak hanya memengaruhi varians (kesalahan), dan bahwa efek tetap hanya memengaruhi rata-rata. Tapi saya telah menemukan contoh di mana pengaruh acak juga mempengaruhi rata-rata - estimasi koefisien:

require(nlme)
set.seed(128)
n <- 100
k <- 5
cat <- as.factor(rep(1:k, each = n))
cat_i <- 1:k # intercept per kategorie
x <- rep(1:n, k)
sigma <- 0.2
alpha <- 0.001
y <- cat_i[cat] + alpha * x + rnorm(n*k, 0, sigma)
plot(x, y)

# simulate missing data
y[c(1:(n/2), (n*k-n/2):(n*k))] <- NA

m1 <- lm(y ~ x)
summary(m1)

m2 <- lm(y ~ cat + x)
summary(m2)

m3 <- lme(y ~ x, random = ~ 1|cat, na.action = na.omit)
summary(m3)

Anda dapat melihat bahwa koefisien estimasi untuk xdari model m1adalah -0.013780, sedangkan dari model m3adalah 0,0011713 - keduanya berbeda secara signifikan dari nol.

Perhatikan bahwa ketika saya menghapus garis yang mensimulasikan data yang hilang, hasilnya sama (ini adalah matriks penuh).

Mengapa demikian?

PS: harap dicatat saya bukan ahli statistik profesional, jadi jika Anda akan merespons dengan banyak matematika maka harap buat juga beberapa ringkasan sederhana untuk boneka :-)

"Saya selalu diajari bahwa efek acak hanya memengaruhi varians (kesalahan), dan bahwa efek tetap hanya memengaruhi rata-rata."

Seperti yang Anda temukan, ini hanya berlaku untuk set data yang seimbang, lengkap (yaitu, tidak ada data yang hilang) tanpa prediktor berkelanjutan. Dengan kata lain, untuk jenis data / model yang dibahas dalam teks ANOVA klasik. Dalam keadaan ideal ini, efek tetap dan efek acak dapat diperkirakan tidak tergantung satu sama lain.

Ketika kondisi ini tidak berlaku (karena sangat sering tidak di "dunia nyata"), efek tetap dan acak tidak independen. Sebagai tambahan yang menarik, inilah mengapa model campuran "modern" diperkirakan menggunakan metode optimasi berulang, daripada diselesaikan dengan sedikit aljabar matriks seperti dalam kasus ANOVA campuran klasik: untuk memperkirakan efek tetap, kita harus mengetahui efek acak, tetapi untuk memperkirakan efek acak, kita harus mengetahui efek tetap! Lebih relevan dengan pertanyaan ini, ini juga berarti bahwa ketika data tidak seimbang / tidak lengkap dan / atau ada prediktor kontinu dalam model, kemudian menyesuaikan struktur efek-acak dari model campuran dapat mengubah perkiraan bagian tetap dari model , dan sebaliknya.

Edit 2016-07-05. Dari komentar: " Bisakah Anda menguraikan atau memberikan kutipan mengapa prediktor terus-menerus akan memengaruhi estimasi bagian tetap dari model? "

Perkiraan untuk bagian yang tetap dari model akan tergantung pada perkiraan untuk bagian acak dari model - yaitu, komponen varians yang diestimasi - jika (tetapi tidak hanya jika) varians dari prediktor berbeda di seluruh cluster. Yang hampir pasti benar jika salah satu prediktor bersifat kontinu (setidaknya dalam data "dunia nyata" - dalam teori akan mungkin hal ini tidak benar, misalnya dalam dataset yang dibangun).

Pada tingkat pertama, saya pikir Anda semua mengabaikan penyusutan terhadap nilai populasi; " lereng dan subjek per subjek dari model efek campuran lebih dekat dengan perkiraan populasi daripada perkiraan kuadrat terkecil dalam subjek. " [ref. 1]. Tautan berikut mungkin juga akan membantu ( Apa deskripsi yang tepat untuk melihat model campuran saya? ), Lihat jawaban Mike Lawrence).

Selain itu, saya pikir Anda sedikit tidak beruntung dalam contoh mainan Anda karena Anda memiliki desain yang seimbang sempurna yang menyebabkan Anda memiliki perkiraan yang sama persis jika tidak ada nilai yang hilang.

Coba kode berikut yang memiliki proses yang sama tanpa nilai yang hilang sekarang:

 cat <- as.factor(sample(1:5, n*k, replace=T) ) #This should be a bit unbalanced.
 cat_i <- 1:k # intercept per kategorie
 x <- rep(1:n, k)
 sigma <- 0.2
 alpha <- 0.001
 y <- cat_i[cat] + alpha * x + rnorm(n*k, 0, sigma) 

 m1 <- lm(y ~ x)  
 m3 <- lme(y ~ x, random = ~ 1|cat, na.action = na.omit) 

 round(digits= 7,fixef(m3)) ==  round(digits=7, coef(m1)) #Not this time lad.
 #(Intercept)           x 
 #      FALSE       FALSE 

Di mana sekarang, karena desain Anda tidak seimbang sempurna, Anda tidak memiliki estimasi koefisien yang sama.

Sebenarnya jika Anda bermain dengan pola nilai Anda yang hilang dengan cara yang konyol (jadi misalnya:) y[ c(1:10, 100 + 1:10, 200 + 1:10, 300 + 1:10, 400 +1:10)] <- NAsehingga desain Anda masih sangat seimbang Anda akan mendapatkan koefisien yang sama lagi.

 require(nlme)
 set.seed(128)
 n <- 100
 k <- 5
 cat <- as.factor(rep(1:k, each = n))
 cat_i <- 1:k # intercept per kategorie
 x <- rep(1:n, k)
 sigma <- 0.2
 alpha <- 0.001
 y <- cat_i[cat] + alpha * x + rnorm(n*k, 0, sigma)
 plot(x, y)

 # simulate missing data in a perfectly balanced way
 y[ c(1:10, 100 + 1:10, 200 + 1:10, 300 + 1:10, 400 +1:10)] <- NA

 m1 <- lm(y ~ x)  
 m3 <- lme(y ~ x, random = ~ 1|cat, na.action = na.omit) 

 round(digits=7,fixef(m3)) ==  round(digits=7, coef(m1)) #Look what happend now...
 #(Intercept)           x 
 #       TRUE        TRUE 

Anda sedikit tersesat oleh desain sempurna percobaan asli Anda. Ketika Anda memasukkan NA dalam non-seimbang Anda mengubah pola berapa banyak "kekuatan" yang bisa dipinjam oleh masing-masing subjek.

Singkatnya, perbedaan yang Anda lihat disebabkan oleh efek penyusutan dan lebih khusus lagi karena Anda mendistorsi desain asli Anda yang seimbang dengan nilai hilang yang tidak seimbang.

Ref 1: Douglas Bates lme4: Pemodelan efek campuran dengan R , halaman 71-72