Mixed model with zero truncated Poisson distributed data
[I'm taking the liberty of cc'ing this to the r-sig-mixed-models mailing list, where it can be archived]
On 12-01-12 12:52 PM, Helen Ward wrote:
Dear Ben, Thank you for your speedy response. I've been having a go with MCMCglmm today, but I suspect I have over simplified the model (I'm new to Bayesian statistics), as I now get three suspiciously significant results... The model I've tried (included below with its summary) runs, but my reading of the help file and the vignette make me suspect that I need to give some sort of prior of something to make it sensible. Am I on the right track? I have also had a quick go with glmmADMB, trying the following models,
model1<-glmmadmb(RSperYr~Age+I(Age^2),data,family="truncpoiss",~1|DadID)
Error in switch(link, log = log, logit = qlogis, probit = qnorm, inverse
= function(x) { :
EXPR must be a length 1 vector
Here the problem is that you are specifying the random effect without giving its name (i.e. random=~1|DadID), so R is trying to interpret it as your desired link function (which is the fourth argument to the glmmadmb() function: see ?glmmadmb).
model1<-glmmadmb(RSperYr~Age+I(Age^2)+(1|DadID),data,family="truncpoiss")
Error in UseMethod("droplevels") :
no applicable method for 'droplevels' applied to an object of class
"c('integer', 'numeric')",
but as you can see from the error mesddages I am obviously doing
something a bit wrong here as well.
As mentioned in another message on this list this morning, the problem here is that you need to explicitly convert "DadID" to a factor from a numeric value [e.g. data$DadID <- factor(data$DadID)] -- the current error message is not very informative.
I'll keep trying!, but are there any more obvious pointers you can give me please? All the best, Helen This is the model I've tried using MCMCglmm.
model2<-MCMCglmm(RSperYr~Age+I(Age^2),random=~DadID,family="ztpoisson",data=data)
summary(model2)
Iterations = 3001:12991 Thinning interval= 10 Sample size= 1000 DIC: 808.5343 G-structure:~DadID post.meanl-95% CI u-95% CI eff.samp DadID0.07409 0.00057460.182239.75 R-structure:~units post.mean l-95% CI u-95% CI eff.samp units0.06205 0.0014050.176322.55 Location effects: RSperYr ~ Age + I(Age^2) post.meanl-95% CIu-95% CI eff.samp pMCMC (Intercept) -0.694625 -1.321090 -0.16384737.17 0.004 ** Age0.1978390.0667140.33041768.64 0.002 ** I(Age^2)-0.008292 -0.014887 -0.00158971.46 0.010 **
At a *quick* glance this looks reasonable. Try plot(model2$Sol) and
plot(model2$VCV) to see if the trace and density plots look reasonable
(compare them with the results from example("MCMCglmm"). You should take
a quick look at the Overview and CourseNotes vignettes (e.g.
vignette("Overview",package="MCMCglmm")) (the latter in particular is a
bit overwhelming but well worth digging into if you're going to use
these methods)
On 11/01/2012 18:46, Ben Bolker wrote:
Helen Ward<h.l.ward at ...> writes:
I would like to describe the relationship between age and male reproductive success in a population of greater horseshoe bats. My data consists of three columns: MaleID, Age, NumberofPups (at that age). Many of the males appear multiple times in the data set, so I believe I need to derive a mixed model with MaleID as a random variable. The data is Poisson distributed, but zero-truncated. So far I have only succeeded in making a mixed model with a poisson distribution (using glmmPQL in the MASS package), and a zero truncated poisson model (using vglm in the VGAM package), but not a mixed model capable of handling zero truncated Poisson data. It has been suggested that I could just minus 1 from each value in the NumberofPups column to make a more usual Poisson distribution, so I can ignore the zero truncated bit. I have tried this and it changes the results of the model, but is this an acceptable transformation? If not, can anyone advise me on a mixed model that can handle zero truncated Poisson data please?
Thanks for letting us know about cross-posting. You should be able to do this in either the MCMCglmm package or (recent versions of) the glmmADMB package. In MCMCglmm, use family="ztpoisson"; in glmmADMB, use family="truncpoiss" ... Ben Bolker
_______________________________________________ R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models