MCMCglmm poisson / not poisson

Hi,

library(MCMCglmm)
l<-rnorm(100,1,2)
y<-rpois(100, exp(l))

dat<-data.frame(y=y)

m1<-MCMCglmm(y~1, data=dat, family="poisson", pl=TRUE)

pp.poisson(y, colMeans(m1$Liab))

# looks bad

pp.poisson(y, colMeans(m1$Liab)+0.5*apply(m1$Liab, 2, var))

# looks better

for(i in 1:1000){
pp.poisson(y, m1$Liab[i,])
}

# looks good

Cheers,

Jarrod





Quoting Mikhail Matz <matz at utexas.edu> on Tue, 28 Aug 2012 22:06:45 -0500:

Hello -

I am playing with ways to justify that the MCMCglmm model fits my data well, which is quite important for me since I am hoping to be able to suggest MCMCglmm-based modeling as a general solution for a particular type of analysis.

I am running "poisson" family on counts data, with two random effects. Following Elston, D. A., R. Moss, et al. (2001). Parasitology 122: 563-569., I am checking whether my lognormal residuals (latent variable minus predicted value) are normally distributed (check), if my random effects (saved with pr=T) are normally distributed (more or less check), and then I try to see if the observed counts really look like Poisson samples based on the latent variables. Again, following Elston et al, I am making a p-p plot using this script (expert coders, please don't judge):

pp.poisson=function(counts,latents) {
	sim=c()
	for(i in 1:length(counts)){
		if (is.na(counts[i])) next
		data=counts[i]
		low=ppois(data,exp(latents[i]))-dpois(data,exp(latents[i]))
		up=ppois(data,exp(latents[i]))
		ss=seq(low,up,(up-low)/100)
		sim=append(sim,sample(ss,1))
	}
	sims=sort(sim)
	xx=(rank(sims)-0.5)/length(sims)
	plot(sims~xx)
	abline(0,1)
}

? and unfortunately it looks really ugly, like a very strongly bent  ' ~ ' rather than a line.
The little script above seems to work; here is a sanity check:

psim=c()
nnn=rnorm(500,10,10)
for (i in 1:length(nnn)){
	psim=append(psim,rpois(1,exp(nnn[i])))
}
pp.poisson(psim,nnn)

I will be extremely grateful for any comments on this.

cheers

Misha
UT Austin

-- 
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.