Hi Users,

I have a panel data set covering the years 2011-18 but mainly focusing on 2015-16. I am attempting to use the xtologit command in Stata 17.0.

I have an ordinal dependent variable which is the number of visits to the GP (general practitioner doctor surgery in the UK) in the last 12 months. Categorised as follows: 1 (Very poor health (10+ visits)); 2 (Poor health (4-6 visits)); 3 (Average health (2-3 visits); 4 (Good health (1-2 visits); 5 (Perfect health (0 visits). Whereby this variable is a proxy for health - those in a low category (1, or 2) are visiting the GP more times and are believed to be in worse health.

In my dataset I have a group of individuals who are 'treated' (by receiving an exogenous increase in income) and a group who are 'control' (who are the untreated / counterfactual). I have control variables for female (1= female; 0=male); furthereduc (1= achieved higher education; 0= did not) and age (continuous). I am using a two period difference-in-differences approach to determine if the increase in income has a causal effect on decreasing the number of GP visits - which is why I am most interested in the interpretation of the estimate on the interaction term.

Below is my code and I am struggling to understand how to get the marginal effect of the interaction between post2016 and treated2016_basecase - as I understand it there is little value in reporting the odds ratios or standalone coefficients.

My intention is to be able to make the following comments from the results "The ATET (being treated and in the post period) means an individual is X% less likely to be in category 1 of visit_gp (i.e. very poor health) and the individual is Y% less likely to be in the category 2 of visit_gp (i.e. poor health) and the individual is Z% more likely to be in the category 5 (i.e. perfect health).

Code:
clear 
use Dataset_Regressions.dta
xtset pidp year
keep if year==2015 | year==2016 
*Remove observations which have missing values
drop if treated2016_basecase==. | visit_gp==. | age==. | agesq==. | furthereduc==. | female==. | married==. | homeowner==. | child==. | extraminch==. | permanent==. | commute==. | smoker==. | numcigs==.
*Any individual that only appears once must be dropped because otherwise regression may capture effects of leaving the labour market
bysort pidp: drop if _N==1
*Check to ensure that the size of treatment and control groups is the same in both the pre and post intervention periods
tab year treated2016_basecase

*Ordered logistic regression for GP visits
eststo drop *
*Regression
eststo: xtologit visit_gp i.post2016##i.treated2016_basecase c.age i.furthereduc i.female
*Marginal effects
eststo: margins, dydx(treated2016_basecase) at(treated2016_basecase=1) at(post2016=1)


The output from the above regression is posted below. And I believe after the regression the coefficients can be interpreted as "being treated in the post intervention period there is a higher change of being in one of the lower categories".

The subsequent margins analysis I am trying to isolate the marginal effects of treatment when treatment =1 and post2016 =1 but this output does not seem correct as I would expect it to be:

1 1 1
2 1 1
3 1 1
4 1 1
5 1 1

Please can someone help clarify the approach in particular to correctly specifying the margins command for my purpose.

Thank you!

Code:
Fitting comparison model:

Iteration 0:   log likelihood = -2355.6421  
Iteration 1:   log likelihood = -2337.7227  
Iteration 2:   log likelihood = -2337.6938  
Iteration 3:   log likelihood = -2337.6938  

Refining starting values:

Grid node 0:   log likelihood = -2286.9491

Fitting full model:

Iteration 0:   log likelihood = -2286.9491  
Iteration 1:   log likelihood = -2209.1284  
Iteration 2:   log likelihood = -2204.8774  
Iteration 3:   log likelihood = -2204.8211  
Iteration 4:   log likelihood = -2204.8211  

Random-effects ordered logistic regression           Number of obs    =  1,690
Group variable: pidp                                 Number of groups =    845

Random effects u_i ~ Gaussian                        Obs per group:
                                                                  min =      2
                                                                  avg =    2.0
                                                                  max =      2

Integration method: mvaghermite                      Integration pts. =     12

                                                     Wald chi2(6)     =  28.67
Log likelihood = -2204.8211                          Prob > chi2      = 0.0001

-----------------------------------------------------------------------------------------------
                     visit_gp | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
------------------------------+----------------------------------------------------------------
                              |
                     post2016 |
           Post Intervention  |  -.1433743   .1195572    -1.20   0.230    -.3777021    .0909535
                              |
         treated2016_basecase |
                     Treated  |   .3348412   .2164677     1.55   0.122    -.0894277    .7591101
                              |
post2016#treated2016_basecase |
   Post Intervention#Treated  |  -.3885202   .2085959    -1.86   0.063    -.7973605    .0203202
                              |
                          age |   -.007706   .0080193    -0.96   0.337    -.0234236    .0080116
                              |
                  furthereduc |
  Achieved Further Education  |    .244885   .1809999     1.35   0.176    -.1098683    .5996383
                              |
                       female |
                      female  |  -.6936186   .1887752    -3.67   0.000    -1.063611   -.3236261
------------------------------+----------------------------------------------------------------
                        /cut1 |  -5.400528   .4615378                     -6.305126   -4.495931
                        /cut2 |   -3.87221   .4405636                     -4.735699   -3.008721
                        /cut3 |  -1.822361     .42432                     -2.654013   -.9907095
                        /cut4 |   1.126816   .4201844                      .3032702    1.950363
------------------------------+----------------------------------------------------------------
                    /sigma2_u |    4.43263   .5114511                      3.535469    5.557453
-----------------------------------------------------------------------------------------------
LR test vs. ologit model: chibar2(01) = 265.75        Prob >= chibar2 = 0.0000
(est1 stored)

. *Marginal effects
. eststo: margins, dydx(treated2016_basecase) at(treated2016_basecase=1) at(post2016=1)

Average marginal effects                                 Number of obs = 1,690
Model VCE: OIM

dy/dx wrt: 1.treated2016_basecase

1._predict: Pr(1.visit_gp), predict(pr outcome(1))
2._predict: Pr(2.visit_gp), predict(pr outcome(2))
3._predict: Pr(3.visit_gp), predict(pr outcome(3))
4._predict: Pr(4.visit_gp), predict(pr outcome(4))
5._predict: Pr(5.visit_gp), predict(pr outcome(5))

1._at: treated2016_basecase = 1
2._at: post2016             = 1

-----------------------------------------------------------------------------------------
                        |            Delta-method
                        |      dy/dx   std. err.      z    P>|z|     [95% conf. interval]
------------------------+----------------------------------------------------------------
0.treated2016_basecase  |  (base outcome)
------------------------+----------------------------------------------------------------
1.treated2016_basecase  |
           _predict#_at |
                   1 1  |  -.0040163   .0062129    -0.65   0.518    -.0161933    .0081607
                   1 2  |   .0019696   .0078703     0.25   0.802    -.0134559    .0173951
                   2 1  |  -.0052691   .0074747    -0.70   0.481    -.0199193    .0093811
                   2 2  |    .002248   .0089558     0.25   0.802    -.0153049     .019801
                   3 1  |  -.0087369   .0112432    -0.78   0.437    -.0307732    .0132995
                   3 2  |   .0031155   .0123593     0.25   0.801    -.0211083    .0273393
                   4 1  |   .0007934    .003308     0.24   0.810    -.0056903     .007277
                   4 2  |  -.0013864   .0056086    -0.25   0.805     -.012379    .0096062
                   5 1  |   .0172289   .0218619     0.79   0.431    -.0256197    .0600775
                   5 2  |  -.0059467   .0235787    -0.25   0.801      -.05216    .0402666
-----------------------------------------------------------------------------------------
Note: dy/dx for factor levels is the discrete change from the base level.
(est2 stored)

A snippet of my code is:

Code:
clear
input long pidp int year float(post2016 treated2016_basecase female furthereduc age)
68017687 2015 0 1 1 0 36
68017687 2016 1 1 1 1 37
68116287 2015 0 0 1 0 33
68116287 2016 1 0 1 0 34
68150967 2015 0 0 0 0 54
68150967 2016 1 0 0 0 55
68150975 2015 0 0 0 0 27
68150975 2016 1 0 0 0 29
68155047 2015 0 0 1 0 59
68155047 2016 1 0 1 0 60
68166607 2015 0 1 1 1 60
68166607 2016 1 1 1 1 61
68197895 2015 0 1 0 0 24
68197895 2016 1 1 0 0 25
68244807 2015 0 0 1 0 38
68244807 2016 1 0 1 0 38
68273371 2015 0 1 1 0 49
68273371 2016 1 1 1 0 50
68293099 2015 0 0 0 1 25
68293099 2016 1 0 0 1 26
68317567 2015 0 1 0 0 36
68317567 2016 1 1 0 0 37
68347487 2015 0 1 0 0 58
68347487 2016 1 1 0 0 59
68351571 2015 0 0 0 1 41
68351571 2016 1 0 0 1 42
68356331 2015 0 0 0 0 32
68356331 2016 1 0 0 0 33
68439289 2015 0 0 1 0 53
68439289 2016 1 0 1 0 54
68469887 2015 0 0 1 1 48
68469887 2016 1 0 1 1 49
68474647 2015 0 1 1 0 41
68474647 2016 1 1 1 0 42
68569173 2015 0 0 1 1 24
68569173 2016 1 0 1 1 25
68573253 2015 0 0 1 1 25
68573253 2016 1 0 1 1 26
68622207 2015 0 0 1 0 57
68622207 2016 1 0 1 0 58
68636487 2015 0 0 1 0 43
68636487 2016 1 0 1 0 44
68757527 2015 0 0 0 1 53
68757527 2016 1 0 0 1 55
69007091 2015 0 0 0 1 44
69007091 2016 1 0 0 1 45
69078487 2015 0 0 1 1 36
69078487 2016 1 0 1 1 37
69224007 2015 0 0 0 1 37
69224007 2016 1 0 0 1 39
69264127 2015 0 0 1 0 39
69264127 2016 1 0 1 0 40
69269571 2015 0 0 0 0 55
69269571 2016 1 0 0 0 56
69280447 2015 0 1 0 1 62
69280447 2016 1 1 0 1 63
69625205 2015 0 0 1 1 45
69625205 2016 1 0 1 1 46
69831249 2015 0 1 1 0 47
69831249 2016 1 1 1 0 48
70069925 2015 0 0 1 0 37
70069925 2016 1 0 1 0 38
70625485 2015 0 0 1 0 45
70625485 2016 1 0 1 0 46
70786645 2015 0 1 1 0 54
70786645 2016 1 1 1 0 55
71148405 2015 0 1 0 0 50
71148405 2016 1 1 0 0 51
71148409 2015 0 0 1 0 48
71148409 2016 1 0 1 0 49
71434005 2015 0 1 1 1 39
71434005 2016 1 1 1 1 40
71461209 2015 0 1 1 1 26
71461209 2016 1 1 1 1 27
77031164 2015 0 1 0 0 29
77031164 2016 1 1 0 0 30
81959729 2015 0 0 0 0 48
81959729 2016 1 0 0 0 49
82171885 2015 0 1 1 0 48
82171885 2016 1 1 1 0 49
88785569 2015 0 0 0 1 57
88785569 2016 1 0 0 1 58
88916809 2015 0 1 1 1 38
88916809 2016 1 1 1 1 39
89099729 2015 0 1 1 0 33
89099729 2016 1 1 1 0 34
89300329 2015 0 0 1 1 51
89300329 2016 1 0 1 1 52
89358125 2015 0 0 0 0 52
89358125 2016 1 0 0 0 53
95278889 2015 0 0 0 1 39
95278889 2016 1 0 0 1 40
95416925 2015 0 0 1 0 41
95416925 2016 1 0 1 0 42
95575369 2015 0 0 0 0 43
95575369 2016 1 0 0 0 44
96029605 2015 0 0 0 1 46
96029605 2016 1 0 0 1 47
96029609 2015 0 0 1 0 47
96029609 2016 1 0 1 0 48
end
label values year b_pregno
label values post2016 post2016
label def post2016 0 "Pre Intervention", modify
label def post2016 1 "Post Intervention", modify
label values treated2016_basecase treated2016_basecase
label def treated2016_basecase 0 "Control", modify
label def treated2016_basecase 1 "Treated", modify
label values female female
label def female 0 "male", modify
label def female 1 "female", modify
label values furthereduc furthereduc
label def furthereduc 0 "No Further Education", modify
label def furthereduc 1 "Achieved Further Education", modify