Saya adalah pengguna yang lebih terbiasa dengan R, dan telah mencoba memperkirakan lereng acak (koefisien seleksi) untuk sekitar 35 individu selama 5 tahun untuk empat variabel habitat. Variabel respons adalah apakah suatu lokasi merupakan habitat "bekas" (1) atau "tersedia" (0) ("digunakan" di bawah).
Saya menggunakan komputer Windows 64-bit.
Dalam R versi 3.1.0, saya menggunakan data dan ekspresi di bawah ini. PS, TH, RS, dan HW adalah efek tetap (standar, jarak diukur ke tipe habitat). lme4 V 1.1-7.
str(dat)
'data.frame': 359756 obs. of 7 variables:
$ use : num 1 1 1 1 1 1 1 1 1 1 ...
$ Year : Factor w/ 5 levels "1","2","3","4",..: 4 4 4 4 4 4 4 4 3 4 ...
$ ID : num 306 306 306 306 306 306 306 306 162 306 ...
$ PS: num -0.32 -0.317 -0.317 -0.318 -0.317 ...
$ TH: num -0.211 -0.211 -0.211 -0.213 -0.22 ...
$ RS: num -0.337 -0.337 -0.337 -0.337 -0.337 ...
$ HW: num -0.0258 -0.19 -0.19 -0.19 -0.4561 ...
glmer(use ~ PS + TH + RS + HW +
(1 + PS + TH + RS + HW |ID/Year),
family = binomial, data = dat, control=glmerControl(optimizer="bobyqa"))
glmer memberi saya perkiraan parameter untuk efek tetap yang masuk akal bagi saya, dan lereng acak (yang saya interpretasikan sebagai koefisien seleksi untuk setiap tipe habitat) juga masuk akal ketika saya menyelidiki data secara kualitatif. Log-kemungkinan untuk model adalah -3050.8.
Namun, sebagian besar penelitian dalam ekologi hewan tidak menggunakan R karena dengan data lokasi hewan, autokorelasi spasial dapat membuat kesalahan standar rentan terhadap kesalahan tipe I. Sementara R menggunakan kesalahan standar berbasis model, kesalahan standar empiris (juga Huber-putih atau sandwich) lebih disukai.
Sementara R saat ini tidak menawarkan opsi ini (sepengetahuan saya - TOLONG, koreksi saya jika saya salah), SAS tidak - meskipun saya tidak memiliki akses ke SAS, seorang rekan setuju untuk membiarkan saya meminjam komputernya untuk menentukan apakah kesalahan standar berubah secara signifikan ketika metode empiris digunakan.
Pertama, kami ingin memastikan bahwa ketika menggunakan kesalahan standar berbasis model, SAS akan menghasilkan perkiraan yang serupa dengan R - untuk memastikan bahwa model tersebut ditentukan dengan cara yang sama di kedua program. Saya tidak peduli apakah mereka persis sama - sama saja. Saya mencoba (SAS V 9.2):
proc glimmix data=dat method=laplace;
class year id;
model use = PS TH RS HW / dist=bin solution ddfm=betwithin;
random intercept PS TH RS HW / subject = year(id) solution type=UN;
run;title;
Saya juga mencoba berbagai bentuk lain, seperti menambahkan garis
random intercept / subject = year(id) solution type=UN;
random intercept PS TH RS HW / subject = id solution type=UN;
Saya mencoba tanpa menentukan
solution type = UN,
atau berkomentar
ddfm=betwithin;
Tidak peduli bagaimana kami menentukan model (dan kami telah mencoba banyak cara), saya tidak bisa mendapatkan kemiringan acak di SAS untuk menyerupai keluaran dari R - meskipun efek tetapnya cukup mirip. Dan ketika saya maksud berbeda, maksud saya bahwa bahkan tanda-tandanya tidak sama. -2 Log Likelihood dalam SAS adalah 71344,94.
Saya tidak dapat mengunggah set data lengkap saya; jadi saya membuat dataset mainan hanya dengan catatan dari tiga orang. SAS memberi saya hasil dalam beberapa menit; di R dibutuhkan lebih dari satu jam. Aneh. Dengan dataset mainan ini, saya sekarang mendapatkan perkiraan berbeda untuk efek tetap.
Pertanyaan saya: Adakah yang bisa menjelaskan mengapa estimasi lereng acak mungkin sangat berbeda antara R dan SAS? Apakah ada yang bisa saya lakukan di R, atau SAS, untuk memodifikasi kode saya sehingga panggilan menghasilkan hasil yang sama? Saya lebih suka mengubah kode dalam SAS, karena saya "percaya" perkiraan R saya lebih banyak.
Saya benar-benar peduli dengan perbedaan-perbedaan ini dan ingin menyelesaikan masalah ini!
Output saya dari dataset mainan yang hanya menggunakan tiga dari 35 individu dalam dataset lengkap untuk R dan SAS dimasukkan sebagai jpegs.
EDIT DAN PEMBARUAN:
Ketika @JakeWestfall membantu menemukan, lereng di SAS tidak termasuk efek tetap. Ketika saya menambahkan efek tetap, inilah hasilnya - membandingkan lereng R dengan lereng SAS untuk satu efek tetap, "PS", di antara program: (Koefisien pemilihan = kemiringan acak). Perhatikan peningkatan variasi dalam SAS.
ID
itu bukan faktor dalam R; periksa dan lihat apakah itu mengubah sesuatu.0
s dan1
s,R
akan memodelkan probabilitas respons "1" sementara SAS akan memodelkan probabilitas respons "0". Untuk membuat model SAS probabilitas "1" Anda harus menuliskan variabel respons Anda sebagaiuse(event='1')
. Tentu saja, bahkan tanpa melakukan ini saya percaya kita masih harus mengharapkan perkiraan yang sama dari varian efek acak, serta perkiraan efek tetap yang sama meskipun dengan tanda-tanda mereka terbalik.ranef()
fungsi daripadacoef()
. Yang pertama memberikan efek acak aktual, sedangkan yang kedua memberikan efek acak ditambah vektor efek tetap. Jadi ini menyumbang banyak alasan mengapa angka-angka yang diilustrasikan dalam posting Anda berbeda, tetapi masih ada perbedaan besar yang belum bisa saya jelaskan sepenuhnya.Jawaban:
Tampaknya saya seharusnya tidak mengharapkan kemiringan acak untuk menjadi serupa antara paket, menurut Zhang et al 2011. Dalam makalah mereka Pada Pemasangan Generalized Linear Efek Campuran Model untuk Respons Biner menggunakan Paket Statistik Berbeda , mereka menggambarkan:
Abstrak:
Saya harap @BenBolker dan tim akan mempertimbangkan pertanyaan saya sebagai suara untuk R untuk menggabungkan kesalahan standar empiris dan kemampuan Gauss-Hermite Quadrature untuk model dengan beberapa istilah kemiringan acak yang akan dimiringkan, karena saya lebih suka antarmuka R dan ingin dapat menerapkan beberapa analisis lebih lanjut dalam program itu. Untungnya, walaupun R dan SAS tidak memiliki nilai yang sebanding untuk lereng acak, tren keseluruhannya sama. Terima kasih atas masukan Anda, saya sangat menghargai waktu dan pertimbangan yang Anda masukkan ke dalam ini!
sumber
Campuran jawaban dan komentar / pertanyaan lainnya:
Saya melengkapi kumpulan data 'mainan' dengan tiga pilihan pengoptimal yang berbeda. (* Catatan 1: mungkin akan lebih berguna untuk tujuan komparatif untuk membuat set data kecil dengan subsampling dari dalam setiap tahun dan id, daripada dengan subsampling variabel pengelompokan. Seperti itu, kita tahu bahwa GLMM tidak akan melakukan sangat baik dengan sejumlah kecil tingkat variabel pengelompokan. Anda dapat melakukan ini melalui sesuatu seperti:
Kode pas batch:
Kemudian saya membaca hasilnya di sesi baru:
Waktu dan penyimpangan yang telah berlalu:
Penyimpangan ini jauh di bawah penyimpangan yang dilaporkan oleh OP dari R (6101.7), dan sedikit di bawah yang dilaporkan oleh OP dari SAS (6078.9), meskipun membandingkan penyimpangan antar paket tidak selalu masuk akal.
Saya benar-benar terkejut bahwa SAS menyatu hanya dalam sekitar 100 evaluasi fungsi!
Rentang waktu dari 17 menit (
nloptwrap
) hingga 80 menit (bobyqa
) pada Macbook Pro, konsisten dengan pengalaman OP. Penyimpangan sedikit lebih baik untuknloptwrap
.Jawaban muncul sangat berbeda dengan
nloptwrap
- meskipun kesalahan standar cukup besar ...(kode di sini memberikan beberapa peringatan tentang hal
year:id
yang harus saya lacak)Bersambung ... ?
sumber