Mengapa nilai taksiran dari Prediktor Linier Linier Terbaik (BLUP) berbeda dari Estimator Linier Tidak Cocok Terbaik (BIRU)?

20

Saya mengerti bahwa perbedaan di antara mereka terkait dengan apakah variabel pengelompokan dalam model diperkirakan sebagai efek tetap atau acak, tetapi tidak jelas bagi saya mengapa mereka tidak sama (jika mereka tidak sama).

Saya secara khusus tertarik pada bagaimana ini bekerja ketika menggunakan estimasi area kecil, jika itu relevan, tetapi saya menduga pertanyaan itu relevan untuk setiap aplikasi efek tetap dan acak.

Jeremy Miles
sumber

Jawaban:

26

Nilai yang Anda dapatkan dari BLUP tidak diestimasi dengan cara yang sama seperti estimasi BLUE untuk efek tetap; menurut konvensi, BLUP disebut sebagai prediksi . Ketika Anda cocok dengan model efek campuran, yang diperkirakan awalnya adalah mean dan varians (dan mungkin kovarians) dari efek acak. Efek acak untuk unit studi yang diberikan (katakanlah seorang siswa) kemudian dihitung dari estimasi rata-rata dan varians, dan data. Dalam model linier sederhana, rata-rata diperkirakan (seperti juga varian residual), tetapi skor yang diamati dianggap terdiri dari keduanya dan kesalahan, yang merupakan variabel acak. Dalam model efek campuran, efek untuk unit yang diberikan juga merupakan variabel acak (meskipun dalam beberapa hal telah terealisasi).

Anda juga dapat memperlakukan unit tersebut sebagai efek tetap, jika diinginkan. Dalam hal ini, parameter untuk unit tersebut diperkirakan seperti biasa. Namun dalam kasus seperti itu, rata-rata (misalnya) dari populasi dari mana unit diambil tidak diperkirakan.

Selain itu, asumsi di balik efek acak adalah bahwa mereka disampel secara acak dari beberapa populasi, dan itu adalah populasi yang Anda sayangi. Asumsi yang mendasari efek tetap adalah bahwa Anda memilih unit-unit itu dengan sengaja karena itu adalah satu-satunya unit yang Anda pedulikan.

Jika Anda berbalik dan menyesuaikan model efek campuran dan memprediksi efek yang sama, mereka cenderung 'menyusut' terhadap rata-rata populasi relatif terhadap perkiraan efek tetap mereka. Anda dapat menganggap ini sebagai analog dengan analisis Bayesian di mana estimasi rata-rata dan varians menentukan sebelum normal dan BLUP seperti rata-rata posterior yang berasal dari penggabungan data secara optimal dengan sebelumnya.

Jumlah susut bervariasi berdasarkan beberapa faktor. Penentuan penting dari seberapa jauh prediksi efek acak akan berasal dari perkiraan efek tetap adalah rasio varians dari efek acak ke varians kesalahan. Berikut ini adalah Rdemo cepat untuk kasus paling sederhana dengan unit 5 'level 2' dengan hanya sarana (penyadapan) yang pas. (Anda dapat menganggap ini sebagai nilai ujian untuk siswa di dalam kelas.)

library(lme4)   # we'll need to use this package
set.seed(1673)  # this makes the example exactly reproducible
nj = 5;    ni = 5;    g = as.factor(rep(c(1:nj), each=ni))

##### model 1
pop.mean = 16;    sigma.g = 1;    sigma.e = 5
r.eff1   = rnorm(nj,    mean=0, sd=sigma.g)
error    = rnorm(nj*ni, mean=0, sd=sigma.e)
y        = pop.mean + rep(r.eff1, each=ni) + error

re.mod1  = lmer(y~(1|g))
fe.mod1  = lm(y~0+g)
df1      = data.frame(fe1=coef(fe.mod1), re1=coef(re.mod1)$g)

##### model 2
pop.mean = 16;    sigma.g = 5;    sigma.e = 5
r.eff2   = rnorm(nj,    mean=0, sd=sigma.g)
error    = rnorm(nj*ni, mean=0, sd=sigma.e)
y        = pop.mean + rep(r.eff2, each=ni) + error

re.mod2  = lmer(y~(1|g))
fe.mod2  = lm(y~0+g)
df2      = data.frame(fe2=coef(fe.mod2), re2=coef(re.mod2)$g)

##### model 3
pop.mean = 16;    sigma.g = 5;    sigma.e = 1
r.eff3   = rnorm(nj,    mean=0, sd=sigma.g)
error    = rnorm(nj*ni, mean=0, sd=sigma.e)
y        = pop.mean + rep(r.eff3, each=ni) + error

re.mod3  = lmer(y~(1|g))
fe.mod3  = lm(y~0+g)
df3      = data.frame(fe3=coef(fe.mod3), re3=coef(re.mod3)$g)

Jadi rasio varians dari efek acak ke varians kesalahan adalah 1/5 untuk model 1, 5/5 untuk model 2, dan 5/1 untuk model 3. Perhatikan bahwa saya menggunakan level means coding untuk model efek tetap. Kita sekarang dapat memeriksa bagaimana perkiraan efek tetap dan efek acak yang diprediksi dibandingkan untuk ketiga skenario ini.

df1
#         fe1     re1
# g1 17.88528 15.9897
# g2 18.38737 15.9897
# g3 14.85108 15.9897
# g4 14.92801 15.9897
# g5 13.89675 15.9897

df2
#          fe2      re2
# g1 10.979130 11.32997
# g2 13.002723 13.14321
# g3 26.118189 24.89537
# g4 12.109896 12.34319
# g5  9.561495 10.05969

df3
#         fe3      re3
# g1 13.08629 13.19965
# g2 16.36932 16.31164
# g3 17.60149 17.47962
# g4 15.51098 15.49802
# g5 13.74309 13.82224

Cara lain untuk berakhir dengan prediksi efek acak yang lebih dekat dengan perkiraan efek tetap adalah ketika Anda memiliki lebih banyak data. Kita dapat membandingkan model 1dari atas, dengan rasio rendah dari efek acak varians ke varians kesalahan, ke versi ( model 1b) dengan rasio yang sama, tetapi lebih banyak data (perhatikan bahwa ni = 500alih - alih ni = 5).

##### model 1b
nj = 5;    ni = 500;    g = as.factor(rep(c(1:nj), each=ni))
pop.mean = 16;    sigma.g = 1;    sigma.e = 5
r.eff1b  = rnorm(nj,    mean=0, sd=sigma.g)
error    = rnorm(nj*ni, mean=0, sd=sigma.e)
y        = pop.mean + rep(r.eff1b, each=ni) + error

re.mod1b = lmer(y~(1|g))
fe.mod1b = lm(y~0+g)
df1b     = data.frame(fe1b=coef(fe.mod1b), re1b=coef(re.mod1b)$g)

Berikut efeknya:

df1
#         fe1     re1
# g1 17.88528 15.9897
# g2 18.38737 15.9897
# g3 14.85108 15.9897
# g4 14.92801 15.9897
# g5 13.89675 15.9897

df1b
#        fe1b     re1b
# g1 15.29064 15.29543
# g2 14.05557 14.08403
# g3 13.97053 14.00061
# g4 16.94697 16.92004
# g5 17.44085 17.40445

Pada catatan yang agak terkait, Doug Bates (penulis paket R lme4) tidak menyukai istilah "BLUP" dan menggunakan "mode kondisional" sebagai gantinya (lihat hlm. 22-23 dari draft lme4 book pdf-nya ). Secara khusus, ia menunjukkan dalam bagian 1.6 bahwa "BLUP" hanya dapat secara bermakna digunakan untuk model efek campuran linier .

gung - Reinstate Monica
sumber
3
+1. Tapi saya tidak yakin saya sepenuhnya menghargai perbedaan terminologis antara "memprediksi" dan "memperkirakan". Jadi parameter distribusi "diperkirakan", tetapi variabel laten hanya dapat "diprediksi"? Apakah saya kemudian mengerti dengan benar bahwa misalnya pemuatan faktor dalam analisis faktor "diperkirakan", tetapi skor faktor "diprediksi"? Terlepas dari itu, saya merasa sangat membingungkan bahwa sesuatu yang disebut "prediktor linier tidak bias terbaik" sebenarnya adalah penaksir yang bias (karena menerapkan penyusutan dan karenanya harus menjadi bias) jika seseorang memperlakukannya sebagai "penaksir" dari efek tetap. ..
amuba kata Reinstate Monica
@amoeba, apa artinya "yang terbaik"? Terbaik apa? Apakah ini merupakan estimasi terbaik dari rata-rata data, atau kombinasi terbaik dari informasi yang terkandung dalam data & sebelumnya? Apakah analogi Bayes membantu Anda?
gung - Reinstate Monica
2
Setidaknya jelas apa artinya "linear" :-) Serius, saya telah menemukan jawaban yang sangat membantu ini oleh @whuber pada perbedaan terminologis antara "prediksi" dan "estimasi". Saya pikir itu menjelaskan terminologi kepada saya, tetapi bahkan memperkuat perasaan saya bahwa BLUP adalah penduga, terlepas dari namanya. [lanjutan]
amuba mengatakan Reinstate Monica
2
@amoeba, ya itu semua masuk akal. Tetapi saya tidak ingin menggunakan nama yang sama untuk keduanya, karena Anda melakukan sesuatu yang berbeda (yaitu persamaan berbeda) & itu berguna untuk nama-nama yang berbeda.
gung - Reinstate Monica
1
@amoeba, saya mengubah frasa pada paragraf pertama untuk menghilangkan tekanan istilah-istilah tersebut, agar tidak mengulangi "prediksi", tetapi untuk mempertahankan perbedaan. Lihat apakah Anda pikir saya memasang jarum atau apakah itu harus diklarifikasi lebih lanjut.
gung - Reinstate Monica