[R-sig-ME] MCMCglmm: how are the predictions computed (conditional and marginal)?

[Joelle Mbatchou]

Hi everyone,

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!

Cheers,
Joelle

	[[alternative HTML version deleted]]