Comparing Multiple Model Types (Poisson, Negative Binomial, ZINB)

Question

I'm working with OTU count data (essentially counts of clustered sequence reads), and trying to determine differential abundance of counts between two groups. As of now, we're using a method defined in Romero, et al's paper titled "The composition and stability of the vaginal microbiota of normal pregnant women is different from that of non-pregnant women", whereby we fit the model using Poisson, NB, and ZINB (applying mixed effects when we have repeated measures), find the model with the best fit (using AIC/BIC), and use those coefficients and p-values to generate q-values and determine significantly-differently abundant OTUs between the two groups that we're comparing.

My questions are as follows: (1) If my model looks like OTU ~ treatment or OTU ~ treatment + (1|studyid) is AIC or BIC better? The same individuals will be in all comparisons, so I don't think that I'm concerned about differences in sample size. Romero uses AIC, but there also seems to be some argument for BIC, and I'm not sure in which instances one is more correct than the other.

(2) Along the same lines, I assume that I can use AIC/BIC for the mixed-effects versions of the Poisson, NB, and ZINB? I'm using the versions implemented in R using the lme4 and glmmADMB packages.

(3) Some of my models won't converge, and to this point, I've taken the models that don't converge and consider them to be non-biologically-relevant anyways. I therefore toss out the resulting coefficients/p-values, and only compare the converging models. Does this make sense? Or do I need to consider that the models that didn't converge may still be useful?

I know that these questions can't be validated on a bench, but I'm not sure how to create a script that will be as true as possible to the bacteria that we find in our samples.

Answer 1

(1) Which one is better is a very theoretical discussion. In practice BIC penalizes the number of parameters more and thus has a stronger preference for parsimonious models (e.g Poisson rather than NB) than AIC. I don't like this approach though, since goodness of fit may be related to the hypothesis tested. The Poisson distribution is definitely too optimistic in terms of variability for microbiome data. I think you should use the same distribution (and the same covariates in the models) for all taxa for comparability.

(2) I guess so, but see (1)

(3) That's a very dangerous approach. Again non-convergence could be related to the hypothesis tested. This could more specifically throw off the multiple testing correction(conversion to q-values) since it relies on the ensemble of P-values