[R-sig-ME] MCMCglmm for binomial models?

[Bob Farmer]

Thanks for the reply.  Unfortunately, I'm not sure that that's the
whole story.  Upping the iterations to achieve an effective sample
size of 1000 (output pasted below) still leads to some differences
between the two models (total runtime on a quad-core Windows box is
under 2 minutes, in case you want to try).  In fact, running it a
second time shooting for an n.eff of 2000 (thin = 325) leads to even
more divergent estimates (e.g. mean.int == -72.96).  Again, assuming
that the glmer() output in this case is a `gold standard', maybe
there's another puzzle piece missing?

Thanks again.
--Bob


#################now binomial from a reference dataset
gm3 <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd),
              family = binomial, data = cbpp)
mc3<-MCMCglmm(cbind(incidence, size - incidence) ~ period,
  random = ~ herd,
  family = "multinomial2",
  data = cbpp, verbose = FALSE,
  nitt = 750E3, thin = 650, burnin = 100E3
)
summary(gm3)@coefs
summary(mc3)

Leads to:

> summary(gm3)@coefs
              Estimate Std. Error   z value     Pr(>|z|)
(Intercept) -1.3985351  0.2278906 -6.136871 8.416284e-10
period2     -0.9923347  0.3053852 -3.249452 1.156274e-03
period3     -1.1286754  0.3260491 -3.461673 5.368286e-04
period4     -1.5803739  0.4288037 -3.685542 2.282169e-04

Versus:

> summary(mc3)

 Iterations = 749351
 Thinning interval  = 100001
 Sample size  = 1000
....
            post.mean l-95% CI u-95% CI eff.samp pMCMC
(Intercept)   -1.6102  -2.2092  -0.8357     1000 0.002 **
period2       -1.2476  -2.2146  -0.1284     1000 0.020 *
period3       -1.3743  -2.3387  -0.2746     1000 0.008 **
period4       -1.9677  -3.2814  -0.8065     1000 0.008 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1