GLMMs for count data with glmmTMB: random slopes specification, cross-level-interaction and strange results

Question

folks,

I recently found the great glmmTMB package which I hoped would help me with my models. My data are 60,000 facebook posts that are nested in 51 companies (i.e., the posts by these companies). My questions are 1) how to incorporate a random slope component exactly in a hurdle model and 2) what's the best overall procedure would look like. 3) I have found some strange effects and wonder if these are trustworthy

The background of the model is:

On the post level, I have a binary predictor "CSRcommunication" (yes vs. no) and the DV is the number of likes (NoL) that the post receives. On the organizational level, have the "reputation of the company as irresponsible" (RR) (wrt. to diverse dimensions, working place, environment etc.). I would like to model a cross-level interaction (posts with a CSR content should backfire when they stem from a bad reputation firm).

The NoL has zero-inflation as well as overdispersion. Hence, I tested a series of models -- but focused on the random intercept model containing my predictors (industry, CSRCom, and RR) - that is Poisson, NegBin, ZIP, ZINB and hurdle model. The winner of this series was the hurdle model. But when I now want to proceed with the random slope model and adding cross-level interaction, I do not know whether I have to incorporate the random slope component in the conditional part only or also in the ziformula.

My code looks like this:

glmmTMB(NoL ~ Industry +  RR_gmcnt*CSRCom_cnt +  (CSRCom_cnt | orgID), 
        ziformula=~RR_gmcnt*CSRCom_cnt +  (CSRCom_cnt | orgID),                         
        data=mldata, 
        family=list(family="truncated_nbinom2",link="log"))

Notes: RR was grand mean centered, and CSRCom was within-firm centered. I had to leave out industry from the ziformula due to convergence problems.

Beyond this question: If I run this code, I get strange results - namely a negative effect of CSRCom in the conditional part indicating a lower NoL when posting about CSR (contrary to expectation). Likewise, the effect of RR is negative in the conditional part but positive in the zero-inflation part. This is especially notable as the effects in the random intercept model look completely different (as expected). In addition, when I add the random slope component, the standard errors explode. There is something problematic here.

As I haven't seen random slopes and cross-level interactions with glmmTMB, I would like to receive some advice. I don't want to draw the wrong conclusions.

Thanks Holger

Answer 1

I do not know whether I have to incorporate the random slope component in the conditional part only or also in the ziformula.

This is also a "theoretical" question. If it makes more sense that Industry predicts whether your outcome is zero or not, maybe it's better to use a "parsimonious" zi-part, namely w/o random effects.

But your model formula in general (cross-level interaction) looks correct.

If I run this code, I get strange results - namely a negative effect of CSRCom in the conditional part indicating a lower NoL when posting about CSR (contrary to expectation).

Have you just looked at the coefficient of the main effect? In presence of interaction, you can't interprete main effects w/o the effect of the interaction. An easy way to see the "real" effect are marginal effects plots. You can use the ggeffects-package, especially for glmmTMB- and GLMMadaptive-models with zero-inflation (because here the confidence intervals for marginal effects need some special attention, see this vignette).

In addition, when I add the random slope component, the standard errors explode. There is something problematic here.

I'm not sure, but this might be because you have a factor (and with only two levels) as random slope. Maybe, but this is something I'm just guessing, not knowing, it's because in such cases more variance-covariance need to be estimated compared to scalar random effect vectors (see http://www.stat.wisc.edu/~bates/UseR2008/WorkshopD.pdf, slide 95).

Answer 2

It is a bit difficult to answer your questions without more information, such as output from the fitted model(s). Nonetheless, a couple of points:

• Your syntax seems to be correct, namely on how to include random effects in the zero-part.
• Package glmmTMB fits the model using the Laplace approximation to calculate the integrals involved in the specification of the likelihood function of the model. However, the Laplace approximation is known not to work optimally, especially as the distribution of data deviate from normality. An alternative procedure that allows you to control the approximation of the integrals is the adaptive Gaussian quadrature. This is implemented in the package GLMMadaptive. Examples of mixed models for count data with over-dispersion and extra zeros can be found here, here, and here.