I haven't used priors specifically for count data, but the conjugate priors for the Poisson and negative-binomial distributions are the gamma and beta distributions, respectively.

See: https://en.wikipedia.org/wiki/Conjugate_prior#When_likelihood_function_is_a_discrete_distribution

David Kaplan's book on Bayesian Statistics also has a brief chapter on multilevel models, and there is this:

Article Bayesian estimation in multilevel modeling.

Regression coefficients usually take Normal priors.

Also, these papers (which I only skimmed) seem to talk about models in a maximum likelihood framework that might be relevant:

Article Multilevel zero-inflated negative binomial regression modeli...

Article Multi-Level Zero-Inflated Poisson Regression Modelling of Co...

Hi Heiko,

thanks a lot. My question was more related to the overall specification of the whole R syntax in the MCMCglmm function in R. This thing is a beast, and the multi-level structure makes it even harder as the fixed and random effects have a covariance matrix for which you have to specify the prior. But even beyond the prior thing, I have absolutely no clue how to specify the syntax for my purposes.

After sending my post, I found perhaps a better solution. The glmmPQL package allows to run GLMM with a either negative binomial and quasipoisson error distribution which would in my case almost do the job. I have no idea how to do a ZINP model in a multilevel context, though but my motivation really sinks every minute to do it perfectly. The two papers you posted are for normal generalized linear models - my problem is the mixed effects /multilevel case... :(

Apart from that I know of no single applied paper that reports a GLMM. I simply have no clue what to put in the table (the output of glmmPQL does not report many things, hence, the usual things that, e.g., Aguinis recommens are not possible.

If you or someone has some example, I would be pleased.

Best,

Holger

packages lme4 , nlme , MCMCglmm deal with poisson (well... algorithms are different, syntax is different )

http://scs.math.yorku.ca/index.php/Mixed_Models_with_R

It has been about three years, but i remember running into the same problem. I attached an overview/tutorial i had saved in the project. But i'm not sure if it's that helpful. If it's not useful enough, I might be able to look through the code and figure out what i actually ended up doing.

Best

Yodi

Dear Fjodor Melnikov

thanks a lot. I will take a look into your files. But actually, I have almost given up. I found a further function (glmmPQL) in the MASS package which look actually promising for estimating neg bin models but the models don't converge.

Finally, I am now up to using a final package before I tend to give up (glmmADMB). This was mentioned by someone who has too convergence problems.

But as mentioned, I will take look

Best,

Hogler

Hi Holger,

Instead of using glmmADMB I recommend using glmmTMB:

https://github.com/glmmTMB/glmmTMB

I think it's written by the same authors as glmmADMB (at least Ben Bolker is involved in both), and I found it relatively easy to use. As far as I can remember, it uses the same syntax as lme4. And the last time I used the two packages I ran into a problem with glmmADMB, which I did not have with glmmTMB.

Good luck,

Urs

Wow Urs Kalbitzer thanks so much!

Best,

Holger

Dear Urs Kalbitzer

I cannot thank you enough. I tried glmmTMB and it worked. It is a tremendous package and offers a simple code that can be applied to all forms of GLMs or GLMMs.

I have, however still a few questions from a rather pragmatical perspective. If you have any references, that would be fine.

1) I started estimating the whole series of possibilities (Poisson, NegBin, ZIP, ZINB and finally a hurdle model). I compared the AIC, LL and deviances which all showed that the hurdle is the best. Any comments?

2) When I then moved from a "random intercept + fixed effects" model (which showed very plausible effects) to a "random slope" model, the effects became either n.s. or were still significant but changed their sign. More importantly, the standard errors blew up. Looks like a typical multicolinarity problem. My guess would be that the combination of a L1-binary predictor and rather low number of clusters (51) leads to a hugh slope-intercept correlation. Comment?

To get the job done, I decided to stick to the former model and don't use this model. Likewise, I avoid estimating the interaction effect.

However, when I tried to estimate the interaction (while the random slope is still fixed) I got a significant interaction. Does that make sense?

3) I have no experiences what to put in a table.

4) ICCs and R-squares are not possible with zero-inflated models. I used the ones from a NegBin model. I guess these estimates are too high due to the wrong zeros cann not be identified. Correct?

I would be very thankful for any ideas.

Best,

Holger

Hi Holger,

I'm glad to hear that the recommendation of the package was helpful.

With regard to your question 2), I recommend to have a look at the help of "isSingular" of the lme4 package. The important part is:

"There is not yet consensus about how to deal with singularity, or more generally to choose which random-effects specification (from a range of choices of varying complexity) to use. Some proposals include:

• avoid fitting overly complex models in the first place, i.e. design experiments/restrict models a priori such that the variance-covariance matrices can be estimated precisely enough to avoid singularity (Matuschek et al 2017)
• use some form of model selection to choose a model that balances predictive accuracy and overfitting/type I error (Bates et al 2015, Matuschek et al 2017)
• “keep it maximal”, i.e. fit the most complex model consistent with the experimental design, removing only terms required to allow a non-singular fit (Barr et al. 2013), or removing further terms based on p-values or AIC
• use a partially Bayesian method that produces maximum a posteriori(MAP) estimates using regularizing priors to force the estimated random-effects variance-covariance matrices away from singularity (Chung et al 2013, blme package)
• use a fully Bayesian method that both regularizes the model via informative priors and gives estimates and credible intervals for all parameters that average over the uncertainty in the random effects parameters (Gelman and Hill 2006, McElreath 2015; MCMCglmm, rstanarm and brms packages)"

(from lme4_1.1-21)

Thus, it's debated which way is the best and it may depend on the dataset. In my opinion, it illustrates one of the issues with making dichotomous decision about being significant or not.

For your other questions, I recommend asking on Stackoverflow/Stackexchange, where I think you are more likely to get answers than here.

All the best,

Urs