Interpreting fixed and random effects from a negative binomial mixed model

[Dougie Jones]

Hi everyone.

For a problem I am trying to tackle, I have created a negative binomial mixed model. The model contains several fixed effects, and only contains a random intercept for each cluster within the data.

use http://www.stata-press.com/data/r13/drvisits

meglm numvisit reform age educ married badh loginc || id:, family(nbinomial) link(log) intpoints(30) irr
predict remeans, remeans

In order to ease the interpretation of the regression coefficients, I'm creating a table (see example below) that describes multiple scenarios in which the fixed effects take different values, and therefore the total estimate of the dependent variable changes. The random intercept is of course a cluster-level effect, so that it does not change at the individual level, and shall be reported as a range.

	Scenario 1	Scenario 2	Scenario 3
Fixed effect 1	Value = 1	Value = 2	Value = 3
Fixed effect 2	Value = 3	Value = 2	Value = 1
Fixed effect 3	Value = 2	Value = 3	Value = 2
Random intercept	[-X, +Y]	[-X, +Y]	[-X, +Y]
Total estimate	sum of above rows

Considering that my model reports incidence rate ratios, I can manually calculate the fixed effects part, but I am unsure how to interpret the reported values of the random intercept.
Can I just add them to the estimate from the fixed effects, or do they need to be rescaled?

[Clyde Schechter]

I don't understand what the [-X, +Y] notation in your table means. In any case, if you have a putative value for the random effect, you can just add it to the fixed-effects xb and then exponentiate to get the expected IRR for that scenario.

That said, this is a very unusual thing to do. It is of no prospective value because if you are dealing with a person outside your estimation sample you do not know what the value of his or her random intercept is: you can only find that out by actually observing her outcomes and then "back-fitting" the model. In other words, it obliterates the predictive value of the model. I suppose you could use it to predict future outcomes for people already in your sample in the event that their covariate values were to change.

The other thing you might want to do, at least for the fixed-effects calculations, is use the -margins- command. You will have to first revise our original -meglm- command to use factor notation so that -margins- can identify which variables are discrete and which are continuous. Then you can specify all the scenarios you are interested in using -at()- options in the -margins- command. It will save you a lot of time, and is also recommended because getting the -lincom- expressions right for all these scenarios is rather error-prone.

Read -help fvvarlist- to learn about factor variable notation. Then, for an exceptionally lucid introduction to the -margins- command, read the excellent Richard Williams' https://www3.nd.edu/~rwilliam/stats/Margins01.pdf.

[Dougie Jones]

Thanks, Clyde.

I don't understand what the [-X, +Y] notation in your table means. In any case, if you have a putative value for the random effect, you can just add it to the fixed-effects xb and then exponentiate to get the expected IRR for that scenario.

The [-X, +Y] notation was my way of indicating the range of values for the random intercept. For my understanding, the values retrieved by --predict randomintercept, remeans-- are on the log scale and not on the linear scale? So they need to be exponentiated first?

That said, this is a very unusual thing to do. It is of no prospective value because if you are dealing with a person outside your estimation sample you do not know what the value of his or her random intercept is: you can only find that out by actually observing her outcomes and then "back-fitting" the model. In other words, it obliterates the predictive value of the model. I suppose you could use it to predict future outcomes for people already in your sample in the event that their covariate values were to change.

My goal was to somehow describe what the effect of the random intercept on the dependent variable at the individual level was. What would be a good way to report the effect of the random effects on individuals?

The other thing you might want to do, at least for the fixed-effects calculations, is use the -margins- command. You will have to first revise our original -meglm- command to use factor notation so that -margins- can identify which variables are discrete and which are continuous. Then you can specify all the scenarios you are interested in using -at()- options in the -margins- command. It will save you a lot of time, and is also recommended because getting the -lincom- expressions right for all these scenarios is rather error-prone.

Read -help fvvarlist- to learn about factor variable notation. Then, for an exceptionally lucid introduction to the -margins- command, read the excellent Richard Williams' https://www3.nd.edu/~rwilliam/stats/Margins01.pdf.

Thanks for the heads-up regarding the factor notation. So far I had, manually, calculated y_hat for each scenario using the IRR (y_hat = B0*B1*...*B99), but perhaps using margins is the more reliable and efficient way.

[Clyde Schechter]

For my understanding, the values retrieved by --predict randomintercept, remeans-- are on the log scale and not on the linear scale? So they need to be exponentiated first?

The values get from -predict remeans- are predictions of the random effects themselves, not their logarithms. Now, when you want to get predicted values of the IRR taking the random effect into account, you add the random effect to the fixed effects and then exponentiate that total. You do not separately exponentiate the predicted remeans.

There is no need to do any of this just to get predicted outcomes, anyhow. You can just -predict, mu- and you will get predicted values that include the random effects. If, however, you want to calculate predicted outcomes at stipulated values of the random effect, then, you will have to go through the process of adding the stipulated random effect value to the predicted xb, adding those, and exponentiating the sum.

[Dougie Jones]

So I am not sure what the best approach is to calculate the how the dependent variable was impacted by variables at the individual level (at specified values), while also showing the effect of the random intercept on the dependent variable. Is that only possible if I would do that for an individual at each cluster?

[Clyde Schechter]

If you have a complete specification of values of all the predictor variables that you are interested in, along with a value of the random effect that you are interested in, you can calculate exp(b0 + b1x1 + b2x2 + ... + bnxn + u). [Here x1 is the value of X1 you are interested in, etc., and u is the value of the random effect you are interested in.] If you have several sets of such specifications, you can just do that in a loop.. Or, better still, create "synthetic" observations having those values and then use -predict- for the -xb- part to simplify the code.