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