How to interpret squared term in GLMER model

I'm trying to explain Y (a count variable), in terms of X, Xsquared, and Size, with random effects for Industry/Firm, by year. My hypothesis is that there is a positive relationship between Y and X, and that this relationship has a inverted-U shape (thus the Xsquared term). I normally do this with normal regression models, but I want to know if it would mean the same in a GLMER one (I'm new with these models)

Do my results, here below, support my inverted-U hypothesis, or something else? Any tips on how to show this inverted-U effect graphically in r? Lastly, how can I interpret the random effects part?

My GLMER model looks like this:

model <- glmer(Y ~ Year + X + Xsquared + Size + (1 + Year|Industry/Firm), data = mydata, family = poisson)


Results:

Random effects:
Groups      Name        Variance    Std.Dev.  Corr
Firm:Industry (Intercept) 1.787e+04 1.337e+02
Year          4.434e-03 6.659e-02 -1.00
Industry      (Intercept) 7.749e-01 8.803e-01
Year          1.923e-07 4.385e-04 -1.00
Number of obs: 436, groups:  Firm:Industry, 109; Industry, 37

Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 51.697639   9.178129   5.633 1.77e-08 ***
Year        -0.027132   0.004535  -5.983 2.19e-09 ***
X            1.322702   0.334092   3.959 7.52e-05 ***
Xsquared    -0.277335   0.141129  -1.965   0.0494 *
Size         0.321026   0.043825   7.325 2.39e-13 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
(Intr) Year   DOI    DOI2
Year     -0.999
DOI       0.004 -0.019
DOI2     -0.009  0.020 -0.953
LNAssets -0.162  0.119  0.003 -0.007


Updated results after changing Year to cYear

Random effects:
Groups        Name        Variance  Std.Dev. Corr
Firm:Industry (Intercept) 5.840e-01 0.764205
cYear       5.546e-03 0.074474 0.38
Industry      (Intercept) 8.243e-05 0.009079
cYear       2.011e-05 0.004485 0.54
Number of obs: 436, groups:  Firm:Industry, 109; Industry, 37

Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.87753    0.47970  -5.999 1.99e-09 ***
cYear       -0.02303    0.02778  -0.829    0.407
X            1.32358    0.33632   3.935 8.30e-05 ***
Xsquared    -0.27628    0.14217  -1.943    0.052 .
FirmSize     0.31789    0.04989   6.371 1.87e-10 ***


Going through your results a bit at a time:

Random effects:
Groups      Name        Variance    Std.Dev.  Corr
Firm:Industry (Intercept) 1.787e+04 1.337e+02
Year          4.434e-03 6.659e-02 -1.00
Industry      (Intercept) 7.749e-01 8.803e-01
Year          1.923e-07 4.385e-04 -1.00


The very small standard deviations of Year, and the perfect negative correlation between year and intercept within each group, suggest that your model is overfitted. I would recommend dropping the Year term and reverting to an intercept-only model (although it probably won't hurt anything if you insist on leaving it there). The among-firm-within-industry standard deviation is very large (approx. 130); if this is really a Poisson regression with a log link, it suggests that your response variable consists of very large and variable numbers of counts.

I would definitely check for outliers in your data, and strongly consider checking for overdispersion/fitting a negative-binomial model/if the counts are large enough, consider just a linear model on log-transformed data.

The variation among Industries is much smaller than the variability at the level of Firm-within-Industry -- not negligible (it's larger than some of your effect sizes).

Number of obs: 436, groups:  Firm:Industry, 109; Industry, 37


This is a reasonable-sized (I would call it 'small to medium') data set for mixed modeling.

Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 51.697639   9.178129   5.633 1.77e-08 ***
Year        -0.027132   0.004535  -5.983 2.19e-09 ***
X            1.322702   0.334092   3.959 7.52e-05 ***
Xsquared    -0.277335   0.141129  -1.965   0.0494 *
Size         0.321026   0.043825   7.325 2.39e-13 ***


Everything is highly significant (except for your X-squared term, which is below the magic $p<0.05$ but just barely!). The intercept is huge, but this is an artifact of the fact that you have Year in your model (probably) as years-since-AD 0; see below. When X=0, the effect of X on the log-counts is positive (1.32 log-counts per unit increase in X), but the curve decelerates/effect of X decreases with increasing X (the effect of X is zero, or the predicted curve of log-counts vs. X peaks, at $X=1.32/0.28 \approx 5$).

Correlation of Fixed Effects:
(Intr) Year   DOI    DOI2
Year     -0.999
DOI       0.004 -0.019
DOI2     -0.009  0.020 -0.953
LNAssets -0.162  0.119  0.003 -0.007


You're probably putting year in as a Gregorian-calendar year (e.g. 2015); this puts the intercept of the model at AD 0, which results in a very strong correlation between the intercept and the Year slope. This probably isn't messing anything up, but it would be clearer and would reduce the chances of numerical problems to set cYear=Year-mean(Year) (or subtract the minimum value, or the median, or any other reasonable value that means the zero value of Year is somewhere near your data).

For graphical solutions, try using predict() to plot predicted values from the model vs. X.

• Hi Ben, I changed the Year variable to your proposed cYear and some things appear better now (except the significance level of Xsquared that it is no longer above the 0.05 threshold...). The among-firm-within-industry standard deviation is now 0.76, which is better right? Any ideas why the Xsquared is less significant now? Is the inverted-U shape still there? (I updated the question with the new results) Commented Aug 15, 2016 at 21:09
• Well, the value has barely changed. The fact that p has changed from 0.049 to 0.052 only means it has crossed the magic "p=0.05" line, not that the results are actually very different. Commented Aug 15, 2016 at 21:16
• That's right. What about the other values, does everything is more or less in reasonable range now (e.g., the standard deviation of among-firm-within-in‌​dustry)? Commented Aug 15, 2016 at 21:21
• Yes. Among-Industry variation is now very small (1-3 orders of magnitude smaller than all other terms in the model). You really need to check for overdispersion whenever you're using a Poisson model, though (see e.g. tinyurl.com/glmmFAQ ) Commented Aug 15, 2016 at 21:24

The interpretation of a squared term is a difference in differences comparing groups differing in 1 unit of the outcome. For instance if the model is linear, and we fit $E[Y|X] = \beta_0 + \beta_1 X$ then $\beta_1 = E[Y|X=9] - E[Y|X=8]$ as well as $\beta_1 = E[Y|X=1] - E[Y|X=0]$. But if the true model is quadratic, then that difference is not constant among all values.

Case 1: model misspecification, if the true model is $E[Y|X] = \beta_0 + \beta_1 X + \beta_2 X^2$ then $E[Y|X=x+1] - E[Y|X=x] = \beta_1 + \beta_2 (2X+1)$. Kind of complicated, but you see the "linear slope" is now a linear function of the $X$. But if you do a difference in differences for $(E[Y|X=x+2] - E[Y|X=x+1]) - (E[Y|X=x+1] - E[Y|X=x]) = 2\beta_2$. So basically the $\beta_1$ is the tangent slope of the quadratic curve at the origin, and $\beta_2$ is a quadratic slope.

With Poisson GLMs, the coefficients are exponentiated and interpreted as relative rates. So $e^{\beta_2}$ would be called a ratio of ratio of rates comparing groups differing by 1 unit differing by 1 unit of $X$. It's basically an interaction term.

The only distinction in mixed models is the interpretation of individual level effects.