[R-sig-ME] trouble getting data scale predictions from a binomial MCMCglmm model

[Kari Ruohonen]

Hi,
I have fitted a binomial model (logit link) with MCMCglmm by fixing the
unit variance to 1 as discussed in the CourseNotes of MCMCglmm.

Now, I would like to get predictions on the data scale. I know there is
a predict function that does this but I need to be able to predict
outside R. With the usual X%*%b where X is the model matrix and b the
vector of coefficients I get predictions on the link scale. The
predictions I get match those from the predict.MCMCglmm method with the
type="terms" option. However, I cannot figure out how to get the
predictions on the data scale. CourseNotes discuss and give an example
of scaling the coefficients suitably with c2 <- (16*sqrt(3)/(15*pi))^2
and then dividing with (sqrt(1+c2*residual var)) to obtain estimates
under residual variance =0 but I cannot figure out how to apply this to
the case of predictions.

Here is a small example based on one in the CourseNotes:

treatment<-gl(2,25)
y<-rbinom(50,1,c(.5,.001)[treatment])
data.bin<-data.frame(treatment=treatment,y=y)
prior.m2c.3<-list(R=list(V=1,fix=1),G=list(G1=list(V=1,nu=.002)))
m2c.3<-MCMCglmm(y~treatment,data=data.bin,family="categorical",prior=prior.m2c.3,verbose=FALSE)

# the next predictions correspond to those from predict.MCMCglmm with "
# type="terms"
preds.link<-m2c.3$X%*%summary(m2c.3$Sol)$statistics[,1]

# but these do not correspond those with type="response"
preds.data<-plogis(preds.link[,1]/sqrt(1+c2*1))

I am obviously missing something and would be grateful to someone
helping me out.

Many thanks, Kari