[R-sig-ME] MCMCglmm: how are the predictions computed (conditional and marginal)?
jmbatchou at uchicago.edu
Fri Jun 10 23:13:53 CEST 2016
I have familial binary data that I fit using a logistic mixed model with MCMCglmm:
prior1 <- list(R = list(V = 1, fix = 1),
G = list(G1 = list(V = 1, nu = 1, alpha.mu=0, alpha.V=10^2)))
mMCMC<-MCMCglmm(Y ~ 1 + age, random =~idU, family = "categorical",
ginverse = list(idU = Prec), prior=prior1,
data = fam_data, nitt=500e3, thin=500, burnin=5e3, pr = T)
To get the predictions integrating out the random effects, I use:
mar_est_pi <- predict(mMCMC, marginal = ~idU, type = "response")
To get the predictions conditional on the posterior means of the parameters, I use:
cond_est_pi <- predict(mMCMC, posterior = "mean", type = "response")
I am trying to understand how I can compute these by hand. For the conditional predictions, I used:
u <- mMCMC$Sol[,-(1:2)]
b <- as.matrix(mMCMC$Sol[,1:2])
cond_pred <- colMeans(plogis(tcrossprod(b, as.matrix(mMCMC$X)) + u))
This doesn't produce the same predictions as in cond_est_pi. I think it's because I am not also incorporating the residual error that was fixed at 1 but I don't see where to incorporate it.
Also, for the marginal predictions, I tried to use the approximation (2.14) in MCMCglmm course notes but I don't get the same predictions as in mar_est_pi...
k <- ((16*sqrt(3))/(15*pi))^2
sig <- rowSums(mMCMC$VCV)
mar_pred <- colMeans(plogis(tcrossprod(b, as.matrix(mMCMC$X))/ sqrt(1 + k * sig)))
Any guidance on how these are computed would be highly appreciated. Thanks!
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models