How to compare melogit coefficients in different models?

[Man Yang]

Hello all, I am running two models using melogit. Both models contain the same variables and a couple of interaction terms as well. I want to test to see whether there are differences of coefficients for a certain variable from one model as opposed to another. Below is the result of model 1:

melogit gedtimehi c.grade##c.grade ib6.RACE ib6.RACE#c.grade || kindergarten:, || studentnum2:, cov(un)

Fitting fixed-effects model:

Iteration 0:   log likelihood = -6232.0048  
Iteration 1:   log likelihood =  -6221.352  
Iteration 2:   log likelihood = -6221.3449  
Iteration 3:   log likelihood = -6221.3449  

Refining starting values:

Grid node 0:   log likelihood = -5295.7391

Fitting full model:

Iteration 0:   log likelihood = -5295.7391  (not concave)
Iteration 1:   log likelihood = -5274.2167  
Iteration 2:   log likelihood = -4906.4134  
Iteration 3:   log likelihood =  -4848.274  
Iteration 4:   log likelihood = -4825.8941  
Iteration 5:   log likelihood = -4824.2776  
Iteration 6:   log likelihood = -4823.8966  
Iteration 7:   log likelihood =  -4823.824  
Iteration 8:   log likelihood = -4823.8111  
Iteration 9:   log likelihood = -4823.8083  
Iteration 10:  log likelihood = -4823.8077  
Iteration 11:  log likelihood = -4823.8077  

Mixed-effects logistic regression               Number of obs      =     11088

-----------------------------------------------------------
                |   No. of       Observations per Group
 Group Variable |   Groups    Minimum    Average    Maximum
----------------+------------------------------------------
   kindergarten |       41          2      270.4        588
    studentnum2 |     4261          1        2.6         10
-----------------------------------------------------------

Integration method: mvaghermite                 Integration points =         7

                                                Wald chi2(8)       =    147.24
Log likelihood = -4823.8077                     Prob > chi2        =    0.0000
------------------------------------------------------------------------------------------
               gedtimehi |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------------------+----------------------------------------------------------------
                gradelvl |   .2788558   .0577974     4.82   0.000      .165575    .3921366
                         |
   c.gradelvl#c.gradelvl |  -.0458036   .0062849    -7.29   0.000    -.0581218   -.0334854
                         |
                    RACE |
       African American  |   -.943746   .2228096    -4.24   0.000    -1.380445   -.5070473
                  Asian  |  -.4564835   .2429516    -1.88   0.060    -.9326599    .0196928
               Hispanic  |   .0272462   .2178825     0.13   0.900    -.3997956     .454288
                         |
         RACE#c.gradelvl |
       African American  |  -.0585262   .0475355    -1.23   0.218     -.151694    .0346415
                  Asian  |   .1767307    .058828     3.00   0.003       .06143    .2920313
               Hispanic  |  -.1326549   .0477209    -2.78   0.005    -.2261861   -.0391237
                         |
                   _cons |    3.26059   .2499694    13.04   0.000     2.770659    3.750522
-------------------------+----------------------------------------------------------------
kindergarten             |
               var(_cons)|   1.161754   .2758725                      .7294345    1.850301
-------------------------+----------------------------------------------------------------
kindergarten>studentnum2 |
               var(_cons)|   13.27539   1.272639                      11.00138    16.01944
------------------------------------------------------------------------------------------
LR test vs. logistic regression:     chi2(2) =  2795.07   Prob > chi2 = 0.0000

Note: LR test is conservative and provided only for reference.

Below is the output of model 2:

. melogit gedtimehi c.grade##c.grade ib6.RACE ib6.RACE#c.grade || studentnum2:, vce(cluster kindergarten)

Fitting fixed-effects model:

Iteration 0:   log likelihood = -6232.0048  
Iteration 1:   log likelihood =  -6221.352  
Iteration 2:   log likelihood = -6221.3449  
Iteration 3:   log likelihood = -6221.3449  

Refining starting values:

Grid node 0:   log likelihood = -5496.5717

Fitting full model:

Iteration 0:   log pseudolikelihood = -5496.5717  
Iteration 1:   log pseudolikelihood = -5073.0641  
Iteration 2:   log pseudolikelihood = -4933.6525  
Iteration 3:   log pseudolikelihood =  -4906.986  
Iteration 4:   log pseudolikelihood = -4905.6049  
Iteration 5:   log pseudolikelihood = -4905.6955  
Iteration 6:   log pseudolikelihood = -4905.7086  
Iteration 7:   log pseudolikelihood = -4905.7108  
Iteration 8:   log pseudolikelihood = -4905.7114  
Iteration 9:   log pseudolikelihood = -4905.7116  

Mixed-effects logistic regression               Number of obs      =     11088
Group variable:     studentnum2                 Number of groups   =      4261

                                                Obs per group: min =         1
                                                               avg =       2.6
                                                               max =        10

Integration method: mvaghermite                 Integration points =         7

                                                Wald chi2(8)       =     71.68
Log pseudolikelihood = -4905.7116               Prob > chi2        =    0.0000
                                   (Std. Err. adjusted for 41 clusters in kindergarten)
---------------------------------------------------------------------------------------
                      |               Robust
            gedtimehi |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
----------------------+----------------------------------------------------------------
             gradelvl |   .2876227   .0737832     3.90   0.000     .1430102    .4322352
                      |
c.gradelvl#c.gradelvl |  -.0457968    .009275    -4.94   0.000    -.0639754   -.0276181
                      |
                 RACE |
    African American  |  -.8329582   .2451435    -3.40   0.001    -1.313431   -.3524858
               Asian  |  -.5298281   .2384022    -2.22   0.026    -.9970879   -.0625684
            Hispanic  |   .2371843   .2034578     1.17   0.244    -.1615857    .6359544
                      |
      RACE#c.gradelvl |
    African American  |  -.0596044   .0575145    -1.04   0.300    -.1723307     .053122
               Asian  |   .1852394   .0773328     2.40   0.017     .0336699    .3368089
            Hispanic  |  -.1416672   .0580928    -2.44   0.015     -.255527   -.0278073
                      |
                _cons |   3.325987   .2647707    12.56   0.000     2.807046    3.844928
----------------------+----------------------------------------------------------------
studentnum2           |
            var(_cons)|   16.41135   2.886854                      11.62553    23.16732
---------------------------------------------------------------------------------------

My question is how should I compare to see whether the racial disparities in model 1 is sig different from that in model 2? Many thanks!

[Clyde Schechter]

Both models have a mis-specification, so you should not compare them. It is incorrect to include RACE#c.gradelvl without also including RACE#c.gradelvl#c.gradelvl. So the bottom-level of both models should look like:

melogit gedtimehi ib6.RACe##c.grade##c.grade || ...whatever

General principle: When you incorporate polynomial representation in a model, the various terms of the polynomial must be treated as a single unit. Whatever you do to any of them you must do to all of them. It makes no sense to interact race with the linear component of the gradelvl effect but not the quadratic component. Similarly any tests of parameters that make reference to gradelvl, should also include the corresponding gradelvl#gradelvl term.

I also don't understand the rationale for comparing these two models, even after the quadratic interaction is put in. If the students are, in fact, nested in kindergartens, then omitting the kindergarten level from the model isn't justified unless the results show that there is no consequential kindergarten effect. The output from the first model shows that the kindergarten-level effect, while not enormous, is at least comparable to that of the race effect at low levels of gradelvl; so they don't seem negligible. The use of vce(cluster kindergarten) is not a substitute way of dealing with this. You might want to use vce(cluster kindergarten) in addition to having a kindergarten: level in the model, but vce(cluster kindergarten) alone does not overcome possible omitted variable bias from leaving out kindergarten-level effefcts. All it does is adjust the standard errors. [You could run the model with and without the kindergarten-level (and no vce(cluster...)) and compare those with a likelihood ratio test, to see if that says that the kindergarten-level effects are not statistically significant. (Though, in truth, I am at best lukewarm about allowing hypothesis tests to determine what goes into a model.)]

Also, to get a helpful response I think you need to define what you mean by "racial disparity." Almost certainly you will get at it by running some -margins- command, but without knowing what kind of disparity you're interested in, it's impossible to say anything more specific.

[Man Yang]

Hi Clyde, thanks for pointing out the quadratic question. I didn't add them in because I've done tests to see whether I should keep the interaction terms that have quadratic grade and the results are not significant (that's why I dropped them in both models). To answer your last question about what do I mean by saying racial disparities, I want to know whether the racial disparities regarding the outcome gedtimehi (whether spending the majority of the time in general education classroom per day) differs from the 1st model as opposed to the 2nd model.

[Clyde Schechter]

I didn't add them in because I've done tests to see whether I should keep the interaction terms that have quadratic grade and the results are not significant (that's why I dropped them in both models).

No, that's not legitimate. The statistical significance of the race#quadratic terms by themselves are completely meaningless for any purpose, least of all for deciding what to include in the model. The implications of doing that are complicated and too long for me to discuss here but you are creating a Frankenstein model by doing that.

To answer your last question about what do I mean by saying racial disparities, I want to know whether the racial disparities regarding the outcome gedtimehi

You have simply restated the question without answering it. What do you mean by "disparities regarding the outcome?" Does disparity refer to racial differences in the levels of the outcome? Or to marginal effects? Or to the location of the vertex of the parabola? Or to something else?

In any case, my advice for understanding any model this complicated is to run -margins- and -marginsplot- afterward. In particular, I would do:

margins race, at(gradelvl = 0(1)12) // OR WHATEVER THE RANGE OF gradelvl IS
marginsplot

and study the graph. If you want to compare the graphs of the two models, you may well get insight into how and how much they make different predictions.

That said, I stand by my comment in #2 that the model without the kindergarten: level effects looks to me like it is, based on what I have seen so far of your results, illegitimate and should not be considered; although it may turn out that the two model predictions are quite similar, in which case I wouldn't make much of it.

[Man Yang]

Hi Clyde, I see your point. Following your suggestion, I ran two models and then ran lrtest to see the model fit but STATA gives me an error saying that it is not valid to conduct lrtest since I am using vce option. It seems I should use Wald test instead. But could you let me know how I should do it in STATA? Below are the output I just got.

 melogit gedtimehi c.grade##c.grade ib6.RACE ib6.RACE#c.grade ib6.RACE#c.grade#c.grade || studentnum2:,
>  cov(un)

Fitting fixed-effects model:

Iteration 0:   log likelihood = -6228.8803  
Iteration 1:   log likelihood = -6218.0137  
Iteration 2:   log likelihood = -6218.0061  
Iteration 3:   log likelihood = -6218.0061  

Refining starting values:

Grid node 0:   log likelihood = -5494.3156

Fitting full model:

Iteration 0:   log likelihood = -5494.3156  
Iteration 1:   log likelihood = -5071.0547  
Iteration 2:   log likelihood = -4931.4398  
Iteration 3:   log likelihood = -4904.8067  
Iteration 4:   log likelihood =  -4903.391  
Iteration 5:   log likelihood = -4903.4923  
Iteration 6:   log likelihood = -4903.5071  
Iteration 7:   log likelihood = -4903.5097  
Iteration 8:   log likelihood = -4903.5103  
Iteration 9:   log likelihood = -4903.5105  

Mixed-effects logistic regression               Number of obs      =     11088
Group variable:     studentnum2                 Number of groups   =      4261

                                                Obs per group: min =         1
                                                               avg =       2.6
                                                               max =        10

Integration method: mvaghermite                 Integration points =         7

                                                Wald chi2(11)      =    145.36
Log likelihood = -4903.5105                     Prob > chi2        =    0.0000
--------------------------------------------------------------------------------------------
                 gedtimehi |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
---------------------------+----------------------------------------------------------------
                  gradelvl |   .4107856   .0985781     4.17   0.000     .2175761     .603995
                           |
     c.gradelvl#c.gradelvl |  -.0625514   .0125194    -5.00   0.000    -.0870889   -.0380139
                           |
                      RACE |
         African American  |  -.6441697   .2395709    -2.69   0.007     -1.11372   -.1746194
                    Asian  |  -.5015813   .2627944    -1.91   0.056    -1.016649    .0134863
                 Hispanic  |   .3611794   .2346533     1.54   0.124    -.0987326    .8210914
                           |
           RACE#c.gradelvl |
         African American  |  -.3007209   .1353421    -2.22   0.026    -.5659866   -.0354552
                    Asian  |   .1517206   .1601457     0.95   0.343    -.1621592    .4656004
                 Hispanic  |  -.3020238   .1343413    -2.25   0.025    -.5653278   -.0387197
                           |
RACE#c.gradelvl#c.gradelvl |
         African American  |   .0327346   .0171376     1.91   0.056    -.0008544    .0663236
                    Asian  |   .0042971    .020607     0.21   0.835    -.0360918     .044686
                 Hispanic  |   .0216956   .0168871     1.28   0.199    -.0114025    .0547937
                           |
                     _cons |   3.242404   .2153224    15.06   0.000      2.82038    3.664429
---------------------------+----------------------------------------------------------------
studentnum2                |
                 var(_cons)|   16.49526   1.662907                      13.53782    20.09879
--------------------------------------------------------------------------------------------
LR test vs. logistic regression: chibar2(01) =  2628.99 Prob>=chibar2 = 0.0000

. estimates store m1

. 
. melogit gedtimehi c.grade##c.grade ib6.RACE ib6.RACE#c.grade ib6.RACE#c.grade#c.grade || kindergarten:
> , || studentnum2:, vce(cluster kindergarten)

Fitting fixed-effects model:

Iteration 0:   log likelihood = -6228.8803  
Iteration 1:   log likelihood = -6218.0137  
Iteration 2:   log likelihood = -6218.0061  
Iteration 3:   log likelihood = -6218.0061  

Refining starting values:

Grid node 0:   log likelihood = -5292.8935

Fitting full model:

Iteration 0:   log pseudolikelihood = -5292.8935  (not concave)
Iteration 1:   log pseudolikelihood = -5270.4464  
Iteration 2:   log pseudolikelihood = -4904.9404  
Iteration 3:   log pseudolikelihood = -4849.1997  
Iteration 4:   log pseudolikelihood = -4824.4203  
Iteration 5:   log pseudolikelihood = -4821.7993  
Iteration 6:   log pseudolikelihood = -4821.2431  
Iteration 7:   log pseudolikelihood = -4821.1977  
Iteration 8:   log pseudolikelihood = -4821.1847  
Iteration 9:   log pseudolikelihood = -4821.1825  
Iteration 10:  log pseudolikelihood = -4821.1821  

Mixed-effects logistic regression               Number of obs      =     11088

-----------------------------------------------------------
                |   No. of       Observations per Group
 Group Variable |   Groups    Minimum    Average    Maximum
----------------+------------------------------------------
   kindergarten |       41          2      270.4        588
    studentnum2 |     4261          1        2.6         10
-----------------------------------------------------------

Integration method: mvaghermite                 Integration points =         7

                                                Wald chi2(11)      =     80.25
Log pseudolikelihood = -4821.1821               Prob > chi2        =    0.0000
                                        (Std. Err. adjusted for 41 clusters in kindergarten)
--------------------------------------------------------------------------------------------
                           |               Robust
                 gedtimehi |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
---------------------------+----------------------------------------------------------------
                  gradelvl |   .4047653   .1362668     2.97   0.003     .1376873    .6718433
                           |
     c.gradelvl#c.gradelvl |  -.0629216   .0202618    -3.11   0.002     -.102634   -.0232092
                           |
                      RACE |
         African American  |  -.7413208   .2868672    -2.58   0.010     -1.30357   -.1790714
                    Asian  |  -.4393565   .2556202    -1.72   0.086    -.9403628    .0616499
                 Hispanic  |   .1475097   .2929123     0.50   0.615    -.4265878    .7216072
                           |
           RACE#c.gradelvl |
         African American  |  -.3189533   .1830427    -1.74   0.081    -.6777105    .0398039
                    Asian  |   .1585661   .1710993     0.93   0.354    -.1767824    .4939146
                 Hispanic  |  -.2916834   .1937665    -1.51   0.132    -.6714587    .0880918
                           |
RACE#c.gradelvl#c.gradelvl |
         African American  |   .0353231   .0234798     1.50   0.132    -.0106964    .0813425
                    Asian  |   .0021023   .0251702     0.08   0.933    -.0472305     .051435
                 Hispanic  |   .0215085   .0254938     0.84   0.399    -.0284585    .0714755
                           |
                     _cons |   3.174641   .2332627    13.61   0.000     2.717454    3.631827
---------------------------+----------------------------------------------------------------
kindergarten               |
                 var(_cons)|   1.170544   .4149911                       .584269    2.345108
---------------------------+----------------------------------------------------------------
kindergarten>studentnum2   |
                 var(_cons)|   13.31292   1.670895                      10.40973    17.02578
--------------------------------------------------------------------------------------------

. estimates store m2

. 
. lrtest m1 m2
LR test likely invalid for models with robust vce
r(498);

Also to respond to your second suggestion in terms of marginsplot, I actually did that but I was wondering whether there are other ways besides visually detecting the differences from the marginsplot regarding racial disparities in terms of the predicted probabilities of spending most of their time in general education classroom? I apologize for not stating my question clear enough.

[Clyde Schechter]

With regard to the likelihood ratio test, if you re-read what I wrote in #2, you will see that it said that you needed to leave out the -vce(cluster...)- in order to do that. There is no Wald test suitable for this particular task of determining whether you need the kindergarten: level in your model.

For my part, with complicated models like this, I am much more comfortable with looking at the graphs or tables of predicted probabilities and making a judgment call as to whether the model predictions are materially different. Statistical tests are extremely sensitive to sample size, and that makes them only minimally useful for selecting what to include in a model. But if you really want to do a formal statistical test, ditch the -vce(cluster...)- part and do the likelihood ratio test.

Well, in addition to doing it graphically, you can present tables of the output from -margins- showing the expected outcomes and their confidence intervals, with the different racial group results side by side for easy comparison. But tables with that many numbers tend to be hard to read, and trends and patterns are difficult to perceive: if ever a picture is worth 10,000 words (or numbers) it is in models that have quadratic dependencies and interaction terms!

If at some point you feel pressed to do a "significance test," to test the null hypothesis that the outcome-gradelvl relationship is the same among all racial groups, that can be done with a Wald test:

testparm i.RACE i.RACE#c.gradelvl i.RACE#c.gradelvl#c.gradelvl

But with such complex relationships and four racial groups, boiling it all down to a single test statistic and p-value really does your subject matter and analysis a disservice. Your audience will deserve better.

[Man Yang]

Thanks so much Clyde!! Your explanations help so much.