Bagaimana cara mendapatkan nilai-p "keseluruhan" dan ukuran efek untuk faktor kategorikal dalam model campuran (lme4)?

Saya ingin mendapatkan nilai-p dan ukuran efek dari variabel kategori independen (dengan beberapa level) - yaitu "keseluruhan" dan tidak untuk setiap level secara terpisah, seperti output normal dari lme4dalam R. Ini seperti hal yang dilaporkan orang saat menjalankan ANOVA.

Bagaimana saya bisa mendapatkan ini?

Kedua konsep yang Anda sebutkan (nilai p dan ukuran efek model campuran linier) memiliki masalah yang melekat. Sehubungan dengan ukuran efek , mengutip Doug Bates, penulis asli lme4,

Dengan asumsi bahwa seseorang ingin mendefinisikan ukuran , saya pikir argumen dapat dibuat untuk memperlakukan jumlah residu kuadrat dari model campuran linier dengan cara yang sama seperti kita mempertimbangkan jumlah residu kuadrat dari model linier. Atau seseorang dapat menggunakan hanya jumlah kuadrat residu tanpa penalti atau jumlah kuadrat residu minimum yang dapat diperoleh dari sekumpulan istilah tertentu, yang sesuai dengan matriks presisi tanpa batas. Saya tidak tahu, sungguh. Itu tergantung pada apa yang Anda coba cirikan.$R^2$

Untuk informasi lebih lanjut, Anda dapat melihat thread ini , thread ini , dan pesan ini . Pada dasarnya, masalahnya adalah bahwa tidak ada metode yang disepakati untuk inklusi dan dekomposisi varians dari efek acak dalam model. Namun, ada beberapa standar yang digunakan. Jika Anda melihat Wiki yang diatur untuk / oleh milis r-sig-model campuran , ada beberapa pendekatan yang terdaftar.

Salah satu metode yang disarankan melihat korelasi antara nilai yang dipasang dan yang diamati. Ini dapat diimplementasikan dalam R seperti yang disarankan oleh Jarrett Byrnes di salah satu utas tersebut:

r2.corr.mer <- function(m) {
  lmfit <-  lm(model.response(model.frame(m)) ~ fitted(m))
  summary(lmfit)$r.squared
}

Jadi misalnya, katakanlah kami memperkirakan model campuran linier berikut:

set.seed(1)
d <- data.frame(y = rnorm(250), x = rnorm(250), z = rnorm(250),
                g = sample(letters[1:4], 250, replace=T)       )
library(lme4)
summary(fm1 <- lmer(y ~ x + (z | g), data=d))
# Linear mixed model fit by REML ['lmerMod']
# Formula: y ~ x + (z | g)
#    Data: d
# REML criterion at convergence: 744.4
# 
# Scaled residuals: 
#     Min      1Q  Median      3Q     Max 
# -2.7808 -0.6123 -0.0244  0.6330  3.5374 
# 
# Random effects:
#  Groups   Name        Variance Std.Dev. Corr 
#  g        (Intercept) 0.006218 0.07885       
#           z           0.001318 0.03631  -1.00
#  Residual             1.121439 1.05898       
# Number of obs: 250, groups: g, 4
# 
# Fixed effects:
#             Estimate Std. Error t value
# (Intercept)  0.02180    0.07795   0.280
# x            0.04446    0.06980   0.637
# 
# Correlation of Fixed Effects:
#   (Intr)
# x -0.005

Kita dapat menghitung ukuran efek menggunakan fungsi yang didefinisikan di atas:

r2.corr.mer(fm1)
# [1] 0.0160841

Alternatif serupa direkomendasikan dalam sebuah makalah oleh Ronghui Xu , disebut sebagai $\Omega^{2}_{0}$ , dan dapat dihitung dalam R hanya:

1-var(residuals(fm1))/(var(model.response(model.frame(fm1))))
# [1] 0.01173721  # Usually, it would be even closer to the value above

Sehubungan dengan nilai-p , ini adalah masalah yang jauh lebih kontroversial (setidaknya di R / lme4komunitas). Lihat diskusi dalam pertanyaan di sini , di sini , dan di sini antara banyak lainnya. Merujuk halaman Wiki lagi, ada beberapa pendekatan untuk menguji hipotesis tentang efek dalam model campuran linier. Terdaftar dari "terburuk hingga terbaik" (menurut penulis halaman Wiki yang saya yakini termasuk Doug Bates serta Ben Bolker yang banyak berkontribusi di sini):

• Wald Z-tes
• Untuk LMM seimbang dan bersarang tempat df dapat dihitung: Wald t-tes
• Uji rasio kemungkinan , baik dengan mengatur model sehingga parameter dapat diisolasi / dijatuhkan (melalui anovaatau drop1), atau melalui komputasi profil kemungkinan
• Interval kepercayaan bootstrap MCMC atau parametrik

Mereka merekomendasikan rantai pengambilan sampel rantai Markov, Monte Carlo, dan juga mendaftar sejumlah kemungkinan untuk menerapkan ini dari pendekatan semu dan sepenuhnya Bayesian, yang tercantum di bawah ini.

Pseudo-Bayesian:

• Pengambilan sampel post-hoc, biasanya (1) dengan asumsi nilai dasar dan (2) mulai dari MLE, mungkin menggunakan perkiraan varians-kovarians perkiraan untuk memilih distribusi kandidat
• Via mcmcsamp(jika tersedia untuk masalah Anda: yaitu LMM dengan efek acak sederhana - bukan GLMM atau efek acak kompleks)
  Melalui pvals.fncdalam languageRpaket, pembungkus untukmcmcsamp )
• Dalam Pembuat Model AD, mungkin melalui glmmADMBpaket (gunakan mcmc=TRUEopsi) atau R2admbpaket (tulis definisi model Anda sendiri dalam Pembuat Model AD), atau di luar R
• Melalui simfungsi dari armpaket (mensimulasikan posterior hanya untuk koefisien beta (efek tetap)

Pendekatan sepenuhnya Bayesian:

• Melalui MCMCglmm paket
• Menggunakan glmmBUGS(bungkus WinBUGS / R antarmuka )
• Menggunakan JAGS / WinBUGS / OpenBUGS dll., Melalui rjags/ r2jags/ R2WinBUGS/ BRugspaket

Demi ilustrasi untuk menunjukkan seperti apa ini, di bawah ini adalah MCMCglmmperkiraan menggunakan MCMCglmmpaket yang Anda akan melihat menghasilkan hasil yang sama seperti model di atas dan memiliki semacam nilai p Bayesian:

library(MCMCglmm)
summary(fm2 <- MCMCglmm(y ~ x, random=~us(z):g, data=d))
# Iterations = 3001:12991
# Thinning interval  = 10
#  Sample size  = 1000 
# 
#  DIC: 697.7438 
# 
#  G-structure:  ~us(z):g
# 
#       post.mean  l-95% CI u-95% CI eff.samp
# z:z.g 0.0004363 1.586e-17 0.001268    397.6
# 
#  R-structure:  ~units
# 
#       post.mean l-95% CI u-95% CI eff.samp
# units    0.9466   0.7926    1.123     1000
# 
#  Location effects: y ~ x 
# 
#             post.mean l-95% CI u-95% CI eff.samp pMCMC
# (Intercept)  -0.04936 -0.17176  0.07502     1000 0.424
# x            -0.07955 -0.19648  0.05811     1000 0.214

Saya harap ini sedikit membantu. Saya pikir saran terbaik untuk seseorang yang memulai dengan model campuran linier dan mencoba memperkirakannya dalam R adalah membaca faq Wiki dari mana sebagian besar informasi ini diambil. Ini adalah sumber yang bagus untuk semua jenis tema efek campuran dari dasar hingga lanjutan dan dari pemodelan hingga perencanaan.

Berkenaan dengan menghitung nilai signifikansi ( p ), Luke (2016) Mengevaluasi signifikansi dalam model efek campuran linier dalam R melaporkan bahwa metode optimal adalah pendekatan Kenward-Roger atau Satterthwaite untuk derajat kebebasan (tersedia dalam R dengan paket seperti lmerTestatau afex).

Abstrak

Model efek-campuran semakin sering digunakan dalam analisis data eksperimental. Namun, dalam paket lme4 dalam R standar untuk mengevaluasi signifikansi efek tetap dalam model ini (yaitu, mendapatkan nilai-p) agak kabur. Ada alasan bagus untuk ini, tetapi karena para peneliti yang menggunakan model-model ini diperlukan dalam banyak kasus untuk melaporkan nilai-p, beberapa metode untuk mengevaluasi signifikansi output model diperlukan. Makalah ini melaporkan hasil simulasi yang menunjukkan bahwa dua metode yang paling umum untuk mengevaluasi signifikansi, menggunakan tes rasio kemungkinan dan menerapkan distribusi z ke nilai Wald dari model output (t-as-z), agak anti-konservatif, terutama untuk ukuran sampel yang lebih kecil. Metode lain untuk mengevaluasi signifikansi,Hasil simulasi ini menunjukkan bahwa tingkat kesalahan Tipe 1 paling dekat dengan 0,05 ketika model dipasang menggunakan REML dan nilai-p diturunkan menggunakan pendekatan Kenward-Roger atau Satterthwaite, karena pendekatan ini menghasilkan tingkat kesalahan Tipe 1 yang dapat diterima bahkan untuk yang lebih kecil sampel.

(penekanan ditambahkan)

Saya menggunakan lmerTestpaket. Ini termasuk perkiraan nilai p dalam anova()output untuk analisis MLM saya, tetapi tidak memberikan ukuran efek karena alasan yang diberikan dalam posting lain di sini.