How to use margins to calculate difference in differences (i.e. interaction) probability

[Paul Sanfilippo]

I have a mixed-effects logistic model with two primary factors (two levels each) of interest; period and DMT_eff. I can calculate margins and pairwise comparisons simple enough:

. margins period#DMT_eff

Predictive margins                                       Number of obs = 2,410
Model VCE: OIM

Expression: Marginal predicted mean, predict()

--------------------------------------------------------------------------------
               |            Delta-method
               |     Margin   std. err.      z    P>|z|     [95% conf. interval]
---------------+----------------------------------------------------------------
period#DMT_eff |
      pre#low  |   .5094355   .0216589    23.52   0.000     .4669848    .5518862
     pre#high  |   .2034331   .0158642    12.82   0.000     .1723399    .2345263
     post#low  |   .4436992   .0230592    19.24   0.000      .398504    .4888943
    post#high  |   .1077383   .0114164     9.44   0.000     .0853626    .1301141
--------------------------------------------------------------------------------

. margins period#DMT_eff, pwcompare(effects)

Pairwise comparisons of predictive margins               Number of obs = 2,410
Model VCE: OIM

Expression: Marginal predicted mean, predict()

--------------------------------------------------------------------------------------------
                           |            Delta-method    Unadjusted           Unadjusted
                           |   Contrast   std. err.      z    P>|z|     [95% conf. interval]
---------------------------+----------------------------------------------------------------
            period#DMT_eff |
  (pre#high) vs (pre#low)  |  -.3060024    .027144   -11.27   0.000    -.3592036   -.2528012
  (post#low) vs (pre#low)  |  -.0657363   .0293404    -2.24   0.025    -.1232425   -.0082302
 (post#high) vs (pre#low)  |  -.4016972   .0245983   -16.33   0.000     -.449909   -.3534854
 (post#low) vs (pre#high)  |   .2402661   .0292099     8.23   0.000     .1830156    .2975165
(post#high) vs (pre#high)  |  -.0956948   .0191122    -5.01   0.000     -.133154   -.0582355
(post#high) vs (post#low)  |  -.3359608   .0259999   -12.92   0.000    -.3869198   -.2850019
--------------------------------------------------------------------------------------------

But I also want to be able to calculate the marginal (probability) difference of differences (e.g. -.3060024 - -.3359608 = 0.0299584). I have tried the following:

. contrast period#DMT_eff, effects

Contrasts of marginal linear predictions

Margins: asbalanced

--------------------------------------------------
               |         df        chi2     P>chi2
---------------+----------------------------------
initiation_bin |
period#DMT_eff |          1        5.23     0.0222
--------------------------------------------------

------------------------------------------------------------------------------------------------
                               |   Contrast   Std. err.      z    P>|z|     [95% conf. interval]
-------------------------------+----------------------------------------------------------------
initiation_bin                 |
                period#DMT_eff |
(post vs base) (high vs base)  |  -.5313962   .2324251    -2.29   0.022    -.9869412   -.0758513
------------------------------------------------------------------------------------------------

. margins r.DMT_eff#r.period, contrast(effects)

Contrasts of predictive margins                          Number of obs = 2,410
Model VCE: OIM

Expression: Marginal predicted mean, predict()

--------------------------------------------------
               |         df        chi2     P>chi2
---------------+----------------------------------
DMT_eff#period |          1        0.76     0.3833
--------------------------------------------------

----------------------------------------------------------------------------------------------
                             |            Delta-method
                             |   Contrast   std. err.      z    P>|z|     [95% conf. interval]
-----------------------------+----------------------------------------------------------------
              DMT_eff#period |
(high vs low) (post vs pre)  |  -.0299584   .0343603    -0.87   0.383    -.0973034    .0373865
----------------------------------------------------------------------------------------------

The contrast command above reproduces the output from the regression model for the interaction term, but displays the estimate as the linear predictor. The margins command below it seems to display the correct probability estimate but the SE, p value, CI don't appear to be correct.

How do I do this?

[Jeff Wooldridge]

Paul: I've been working on nonlinear DiD lately and I know full well how confusing this can be. I'm assuming the treated group/period is represented by post#high. It's easiest to define the treatment variable explicitly. If you have panel data, you should add vce(cluster id) to the probit. The vce(robust) option means you'd be admitting the probit model is incorrect.

gen treat = post*high
probit y high post i.treat, vce(robust)
margins, dydx(treat) subpop(if treat == 1)

What you're estimating is the average treatment effect on the treated -- the parameter of most interest in these applications. In the linear case, there's no distinction between that and the overall average treatment effect if you don't have controls.

You can see more examples and my shared Dropbox, as well as my forthcoming paper in the Econometrics Journal:

https://www.dropbox.com/sh/zj91darud...bgsnxS6Za?dl=0

[Paul Sanfilippo]

Thanks very much Jeff. Unfortunately it doesn't work for me. I manually created the interaction variable and tried to follow your margins example, but it doesn't give the DiD probability estimate I was expecting. Not sure why...

[Jeff Wooldridge]

Can you show your whole set of commands? I'm wondering if your calculation is just reproducing the linear model estimates. This happens if you are using the Ai-Norton definition of the effect of the interaction. You'll always get the same as the linear model.

[Paul Sanfilippo]

Thanks Jeff, TBH I could very easily be doing something inconsistent with your advice. The output is below:

. gen period2 = period - 1         // pre = 0, post = 1

. gen DMT_eff2 = DMT_eff - 1       // low = 0, high = 1

. gen interact = period2*DMT_eff2  // manually create interaction variable

. melogit initiation_bin i.period2 i.DMT_eff2 i.interact || PATIENT_ID:

Fitting fixed-effects model:

Iteration 0:  Log likelihood = -1319.1925  
Iteration 1:  Log likelihood = -1317.4895  
Iteration 2:  Log likelihood = -1317.4831  
Iteration 3:  Log likelihood = -1317.4831  

Refining starting values:

Grid node 0:  Log likelihood = -1319.0203

Fitting full model:

Iteration 0:  Log likelihood = -1319.0203  
Iteration 1:  Log likelihood = -1311.9241  
Iteration 2:  Log likelihood = -1308.5544  
Iteration 3:  Log likelihood = -1308.4624  
Iteration 4:  Log likelihood = -1308.4623  

Mixed-effects logistic regression               Number of obs     =      2,525
Group variable: PATIENT_ID                      Number of groups  =      1,628

                                                Obs per group:
                                                              min =          1
                                                              avg =        1.6
                                                              max =         36

Integration method: mvaghermite                 Integration pts.  =          7

                                                Wald chi2(3)      =     269.71
Log likelihood = -1308.4623                     Prob > chi2       =     0.0000
--------------------------------------------------------------------------------
initiation_bin | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
---------------+----------------------------------------------------------------
     1.period2 |  -.1669533   .1453507    -1.15   0.251    -.4518355    .1179288
    1.DMT_eff2 |  -1.896001   .1591477   -11.91   0.000    -2.207925   -1.584077
    1.interact |  -.4988261   .2125128    -2.35   0.019    -.9153436   -.0823087
         _cons |   .3379402   .1060661     3.19   0.001     .1300544     .545826
---------------+----------------------------------------------------------------
PATIENT_ID     |
     var(_cons)|   .4738878   .1765465                       .228329    .9835353
--------------------------------------------------------------------------------
LR test vs. logistic model: chibar2(01) = 18.04       Prob >= chibar2 = 0.0000

. margins period2#DMT_eff2

Predictive margins                                       Number of obs = 2,525
Model VCE: OIM

Expression: Marginal predicted mean, predict()

----------------------------------------------------------------------------------
                 |            Delta-method
                 |     Margin   std. err.      z    P>|z|     [95% conf. interval]
-----------------+----------------------------------------------------------------
period2#DMT_eff2 |
            0 0  |   .5409645   .0322548    16.77   0.000     .4777462    .6041828
            0 1  |   .1746616   .0123434    14.15   0.000     .1504689    .1988542
            1 0  |   .5037402   .0215782    23.34   0.000     .4614478    .5460327
            1 1  |   .1533483   .0165154     9.29   0.000     .1209787    .1857178
----------------------------------------------------------------------------------

. margins, dydx(interact) subpop(if interact == 1)
warning: prediction constant over observations.

Average marginal effects                               Number of obs   = 2,525
Model VCE: OIM                                         Subpop. no. obs =   784

Expression: Marginal predicted mean, predict()
dy/dx wrt:  1.interact

------------------------------------------------------------------------------
             |            Delta-method
             |      dy/dx   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
    interact |
          0  |          0  (empty)
          1  |   -.056891   .0260061    -2.19   0.029    -.1078621     -.00592
------------------------------------------------------------------------------
Note: dy/dx for factor levels is the discrete change from the base level.

.
end of do-file

[Paul Sanfilippo]

Interestingly I just tested some of these commands in the linear model case, and:

contrast period#DMT_eff, effects                           // from model with interaction syntax
margins r.DMT_eff#r.period, contrast(effects)       // from model with interaction syntax
margins, dydx(interact) subpop(if interact == 1)   // from model with manually created interaction

all give the exact same contrast estimates, SE's, p values and CI's (as the interaction term output from the regression model). So clearly in the non-linear case they are estimating different things and I am just not quite sure what...

[Jeff Wooldridge]

Sorry, I missed it was melogit. I have to figure out what Stata is doing here. I think the same issues of properly defining the partial effect arise ....