Apakah struktur-R-struktur G dalam sebuah glmm?

Saya telah menggunakan MCMCglmmpaket baru-baru ini. Saya bingung dengan apa yang disebut dalam dokumentasi sebagai struktur-R dan struktur-G. Ini tampaknya berhubungan dengan efek acak - khususnya menentukan parameter untuk distribusi sebelumnya, tetapi diskusi dalam dokumentasi tampaknya mengasumsikan bahwa pembaca tahu apa istilah-istilah ini. Sebagai contoh:

daftar opsional spesifikasi sebelumnya yang memiliki 3 elemen yang mungkin: R (R-struktur) G (G-struktur) dan B (efek tetap) ............ Prior untuk struktur varians (R dan G ) adalah daftar dengan varian (co) (V) yang diharapkan dan parameter derajat kepercayaan (nu) untuk invers-Wishart

... diambil dari sini .

EDIT: Harap dicatat bahwa saya telah menulis ulang sisa pertanyaan setelah komentar dari Stephane.

Adakah yang bisa menjelaskan apa itu R-struktur dan G-struktur, dalam konteks model komponen varians sederhana di mana prediktor linier adalah dengan dan

β_{0} + e_{0 i j} + u_{0 j}

$\beta_0 + e_{0ij} + u_{0j}$

e_{0 i j} \sim N (0, σ_{0 e}^{2})

$e_{0ij} \sim N(0,\sigma_{0e}^2)$

u_{0 j} \sim N (0, σ_{0 u}^{2})

$u_{0j} \sim N(0,\sigma_{0u}^2)$

Saya membuat contoh berikut dengan beberapa data yang disertakan MCMCglmm

> require(MCMCglmm)
> require(lme4)
> data(PlodiaRB)
> prior1 = list(R = list(V = 1, fix=1), G = list(G1 = list(V = 1, nu = 0.002)))
> m1 <- MCMCglmm(Pupated ~1, random = ~FSfamily, family = "categorical", 
+ data = PlodiaRB, prior = prior1, verbose = FALSE)
> summary(m1)


 G-structure:  ~FSfamily

         post.mean l-95% CI u-95% CI eff.samp
FSfamily    0.8529   0.2951    1.455      160

 R-structure:  ~units

      post.mean l-95% CI u-95% CI eff.samp
units         1        1        1        0

 Location effects: Pupated ~ 1 

            post.mean l-95% CI u-95% CI eff.samp  pMCMC    
(Intercept)   -1.1630  -1.4558  -0.8119    463.1 <0.001 ***
---

> prior2 = list(R = list(V = 1, nu = 0), G = list(G1 = list(V = 1, nu = 0.002)))
> m2 <- MCMCglmm(Pupated ~1, random = ~FSfamily, family = "categorical", 
+ data = PlodiaRB, prior = prior2, verbose = FALSE)
> summary(m2)


 G-structure:  ~FSfamily

         post.mean l-95% CI u-95% CI eff.samp
FSfamily    0.8325   0.3101    1.438    79.25

 R-structure:  ~units

      post.mean l-95% CI u-95% CI eff.samp
units    0.7212  0.04808    2.427    3.125

 Location effects: Pupated ~ 1 

            post.mean l-95% CI u-95% CI eff.samp  pMCMC    
(Intercept)   -1.1042  -1.5191  -0.7078    20.99 <0.001 ***
---

> m2 <- glmer(Pupated ~ 1+ (1|FSfamily), family="binomial",data=PlodiaRB)
> summary(m2)
Generalized linear mixed model fit by the Laplace approximation 
Formula: Pupated ~ 1 + (1 | FSfamily) 
   Data: PlodiaRB 
  AIC  BIC logLik deviance
 1020 1029   -508     1016
Random effects:
 Groups   Name        Variance Std.Dev.
 FSfamily (Intercept) 0.56023  0.74849 
Number of obs: 874, groups: FSfamily, 49

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  -0.9861     0.1344  -7.336  2.2e-13 ***

Jadi berdasarkan komentar dari Stephane saya pikir struktur G adalah untuk . Tetapi komentar juga mengatakan bahwa struktur R adalah untuk namun ini tampaknya tidak muncul dalam output. $\sigma_{0u}^2$ $\sigma_{0e}^2$ lme4

Perhatikan bahwa hasil dari lme4/glmer()konsisten dengan kedua contoh dari MCMC MCMCglmm.

Jadi, apakah struktur R untuk dan mengapa ini tidak muncul pada output ? $\sigma_{0e}^2$ lme4/glmer()

r bayesian mixed-model lme4-nlme Joe King
sumber

Dengan terminologi SAS (tetapi mungkin terminologi yang lebih umum), matriks G adalah matriks varians dari efek acak dan matriks R adalah matriks varians dari "istilah kesalahan" (dalam kasus Anda, mungkin itu adalah perkiraan residual yang diperkirakan). variance ?)

σ_{0 e}^{2}

$\sigma_{0e}^2$

Stéphane Laurent

@ StéphaneLaurent terima kasih. Saya bertanya-tanya apakah mungkin diperkirakan tetapi ketika saya pertama kali belajar tentang model linier umum saya ingat bahwa tidak diperkirakan - hanya "penyimpangan" dihitung (seperti halnya ). Mungkin saya kehilangan sesuatu?

σ_{0 e}^{2}

$\sigma_{0e}^2$

σ_{0 e}^{2}

$\sigma_{0e}^2$ lme4

Joe King

mungkin rasa varian residual tidak jelas ketika keluarga distribusi bukan yang Gaussian

Stéphane Laurent

@ Stéphane Laurent Ya! Silakan lihat komentar saya untuk jawaban Michael semenit yang lalu - untuk hasil biner, ini harus diperbaiki (seperti pada model saya di OP)

Joe King

Saat Anda memiliki model ME / Multilevel, ada beberapa varian. Bayangkan kasus paling sederhana: . Ada variasi dalam intersepsi , dan dalam istilah kesalahan . sering digunakan untuk matriks var-covar dari efek acak (dalam hal ini skalar, ) & adalah untuk matriks var-covar dari varian residual setelah memperhitungkan fixed & that cluster's random efek. Biasanya dipahami sebagai matriks diagonal dari . Juga, kedua dist dianggap sebagai multivariat normal dengan rata-rata = 0.

Y_{i} = β_{0} + β_{1} X + b_{i} + ε_{i}

$Y_i=\beta_0+\beta_1X+b_i+\varepsilon_i$

b_{i}

$b_i$

ε_{i}

$\varepsilon_i$

G

$G$

σ_{b}^{2}

$\sigma^2_b$

R_{i}

$R_i$

ε_{i}

$\varepsilon_i$

σ^{2}

$\sigma^2$

gung - Reinstate Monica

Jawaban:

Saya lebih suka memposting komentar saya di bawah ini sebagai komentar tetapi ini tidak cukup. Ini adalah pertanyaan alih-alih jawaban (sama dengan @gung saya tidak merasa cukup kuat pada topik).

Saya mendapat kesan bahwa MCMCglmm tidak mengimplementasikan Bayesian glmm yang "benar". Model Bayesian yang sebenarnya dijelaskan dalam bagian 2 tulisan ini . Demikian pula dengan model frequentist, kita memiliki dan ada yang diperlukan sebelumnya pada parameter dispersi selain parameter tetap dan varian "G" dari efek acak . $g(E(y \mid u)) = X\beta + Zu$ $\phi_1$ $\beta$ $u$

Tetapi menurut sketsa MCMCglmm ini , model yang diimplementasikan dalam MCMCglmm diberikan oleh , dan tidak melibatkan parameter dispersi . Ini tidak mirip dengan model frequentist klasik. $g(E(y \mid u,e)) = X\beta + Zu + e$ $\phi_1$

Karena itu saya tidak akan terkejut bahwa tidak ada analog dengan glmer. $\sigma_e$

Mohon maaf atas komentar kasar ini, saya hanya melihat sekilas tentang itu.

Stéphane Laurent
sumber

Terima kasih. Apakah topik ini seharusnya sulit, karena saya merasa cukup sulit? Saya pikir saya puas dengan makna struktur R dan G sekarang. Saya masih bingung tentang kurangnya dengan dan saya sangat ingin tahu tentang komentar Anda yang tidak benar-benar Bayesian. Saya tidak dapat dengan jujur mengatakan bahwa saya memahami semua makalah yang Anda tautkan dan saya juga berjuang dengan bagian sketsa, tetapi hanya murni dari perspektif contoh saya, saya percaya parameter dispersi harus konstan (karena contohnya adalah binomial). Apa yang saya lewatkan?

σ_{e}

$\sigma_e$ glmerMCMCglmmMCMCglmm

ϕ_{1}

$\phi_1$

Joe King

Maaf, kata-kata saya tidak sepenuhnya tepat. MCMCglmm benar-benar Bayesian, tetapi tidak persis menerapkan glmm klasik (saya pikir). Selain itu Anda harus menyadari bahwa sulit untuk menetapkan prior menghasilkan inferensi pada komponen varians dekat dengan inferensi frequentist.

Stéphane Laurent

Terima kasih lagi. Dalam penelitian saya, saya telah menemukan bahwa saya dapat menggunakan distribusi invers-wishart default untuk komponen varians dalam MCMCglmmmenggunakan berbagai parameter, dan interval yang kredibel 95% selalu berisi nilai varians untuk perkiraan efek acak glmersehingga saya merasa ini masuk akal , tetapi bagaimana saya harus menafsirkan kasus ini, yang mungkin tidak khas, di mana hasilnya bahwa MCMCglmmintervalnya tidak terlalu sensitif terhadap pilihan sebelumnya? Mungkin saya harus mengajukan pertanyaan baru tentang ini?

Joe King

Mungkin Anda memiliki ukuran sampel yang besar? Berkenaan dengan pertanyaan awal Anda, saya mendapat kesan bahwa, setidaknya untuk kasus binomial, model glmer setara dengan model MCMCglmm dengan

. Apa yang terjadi jika Anda menetapkan sebelumnya pada sangat terkonsentrasi pada ?

σ_{e} = 0

$\sigma_e=0$

σ_{e}

$\sigma_e$

0

$0$

Stéphane Laurent

Ya, saya memiliki ukuran sampel yang cukup besar: 50.000 pengamatan di 225 cluster (data saya sendiri, bukan contoh dalam pertanyaan saya). Ketika saya menetapkan sebelumnya yang sangat terkonsentrasi di dekat nol pada , dengan menetapkan V = 0,01 dan nu = 100 maka saya mendapatkan 0,25 (CI: 0,16, 0,29) untuk dan 0,53 (0,38, 0,73) untuk . Ketika saya menetapkan sebelumnya kurang informatif, dengan V = 10 dan nu = 0,01 maka saya mendapatkan 0,18 (0,12, 0,23) dan 0,49 (0,34, 0,63) masing-masing. Ini dibandingkan dengan 0,51 dari . Saya bahkan mencoba flat yang tidak tepat sebelumnya, yang memberi 0,10 (0,08, 0,13) dan 0,47 (0,25, 0,68).

σ_{e}

$\sigma_e$

σ_{e}

$\sigma_e$

σ_{u}

$\sigma_u$ glmer

Joe King

Saya terlambat ke permainan, tetapi beberapa catatan. The struktur adalah struktur residual. Dalam kasus Anda, "struktur" hanya memiliki satu elemen (tetapi ini tidak harus terjadi). Untuk variabel respon Gaussian, varians residual, biasanya diperkirakan. Untuk hasil biner, ia dipertahankan konstan. Karena cara pengaturan MCMCglmm , Anda tidak dapat memperbaikinya pada nol, tetapi relatif standar untuk memperbaikinya pada (juga berlaku untuk model probit). Untuk menghitung data (misalnya, dengan distribusi poisson), Anda tidak memperbaikinya dan ini secara otomatis memperkirakan parameter penyebaran berlebihan pada dasarnya. $\mathbf{R}$ $\sigma^{2}_{e}$ $1$

The struktur adalah struktur efek acak. Sekali lagi dalam kasus Anda, hanya mencegat acak, tetapi jika Anda memiliki beberapa efek acak, mereka akan membentuk matriks varians-kovarians, . $\mathbf{G}$ $\mathbf{G}$

Catatan akhir, karena varians residual tidak tetap pada nol, estimasi tidak akan cocok dengan yang dari glmer. Anda perlu mengubah skala mereka. Ini adalah sedikit contoh (tidak menggunakan efek acak, tetapi menggeneralisasi). Perhatikan bagaimana varians struktur R tetap pada 1.

# example showing how close the match is to ML without separation
m2 <- MCMCglmm(vs ~ mpg, data = mtcars, family = "categorical",
  prior = list(
    B = list(mu = c(0, 0), V = diag(2) * 1e10),
    R = list(V = 1, fix = 1)),
  nitt = 1e6, thin = 500, burnin = 10000)
summary(m2)

Berikut adalah konstanta penskal untuk keluarga binomial:

k <- ((16*sqrt(3))/(15*pi))^2

Sekarang bagi solusinya, dan dapatkan mode posterior

posterior.mode(m2$Sol/(sqrt(1 + k)))

Yang harusnya cukup dekat dengan apa yang kita dapatkan glm

summary(glm(vs ~mpg, data = mtcars, family = binomial))

Joshua
sumber

apakah Anda tahu cara menentukan heteroskedastisitas pada level satu di MCMCglmm? Apakah itu struktur R? Apa sintaksisnya?

Maxim.K

@ Yosua, bisakah Anda menjelaskan "konstanta penskal untuk keluarga binomial"? PS: Untuk seed 123, saya dapatkan (dengan koreksi) dari m2nilai -8.164dan 0.421; dan dari glmnilai -8.833dan 0.430.

Qaswed

Konstanta penskalaan dapat ditemukan di Diggle et. Al. ( amazon.de/Analysis-Longitudinal-Oxford-Statistics-Science/dp/… ) - menurut cran.r-project.org/web/packages/MCMCglmm/vignettes/… eq. 2.14 di halaman 47.

Qaswed