HR using margins command after Cox regression (interaction terms)

[Megan Moreton]

I am comparing mortality between two treatment groups using a Cox model. I want to assess for subgroup differences using an interaction term (between the treatment group and the subgroup) as a heterogeneity test.

Eg,:

. stcox randgroup##age_binary, allbase

         Failure _d: death==1
   Analysis time _t: day_eos
  Exit on or before: time 168

Iteration 0:   log likelihood = -743.21888
Iteration 1:   log likelihood = -738.52467
Iteration 2:   log likelihood = -736.31973
Iteration 3:   log likelihood = -736.29266
Iteration 4:   log likelihood = -736.29266
Refining estimates:
Iteration 0:   log likelihood = -736.29266

Cox regression with Breslow method for ties

No. of subjects =    555                                Number of obs =    555
No. of failures =    120
Time at risk    = 76,776
                                                        LR chi2(3)    =  13.85
Log likelihood = -736.29266                             Prob > chi2   = 0.0031

--------------------------------------------------------------------------------------
                  _t | Haz. ratio   Std. err.      z    P>|z|     [95% conf. interval]
---------------------+----------------------------------------------------------------
           randgroup |
              Arm A  |          1  (base)
              Arm B  |   2.820427   1.076045     2.72   0.007      1.33526      5.9575
                     |
          age_binary |
          <20 years  |          1  (base)
         >=20 years  |   1.028425   .3618726     0.08   0.937     .5160121    2.049677
                     |
randgroup#age_binary |
    Arm A#<20 years  |          1  (base)
   Arm A#>=20 years  |          1  (base)
    Arm B#<20 years  |          1  (base)
   Arm B#>=20 years  |   .3911957   .1710112    -2.15   0.032       .16607    .9215034
--------------------------------------------------------------------------------------

When there appears to be a subgroup effect - like in the above example - I want to report the HR for the treatment group effect (ie, Arm A vs Arm B) at each of the subgroup levels (ie, at age <20 and age >=20).

I thought I was getting there, but I'm getting very muddled with the margins command that I need to run to give these HRs. Do you have any advice?

Is it the following:

. margins randgroup#age_binary

Adjusted predictions                                       Number of obs = 555
Model VCE: OIM

Expression: Predicted hazard ratio, predict()

--------------------------------------------------------------------------------------
                     |            Delta-method
                     |     Margin   std. err.      z    P>|z|     [95% conf. interval]
---------------------+----------------------------------------------------------------
randgroup#age_binary |
    Arm A#<20 years  |          1          .        .       .            .           .
   Arm A#>=20 years  |   1.028425   .3618726     2.84   0.004     .3191679    1.737682
    Arm B#<20 years  |   2.820427   1.076045     2.62   0.009     .7114179    4.929437
   Arm B#>=20 years  |   1.134702   .3959151     2.87   0.004     .3587222    1.910681
--------------------------------------------------------------------------------------

Thank you in advance for any help!

[Clyde Schechter]

The -margins- output you show gives you the ratio of hazard ratios (relative hazard ratio, if you prefer) in each of the four subgroups defined by age and randomization group, relative to the reference group Arm A < 20 years of age. But that isn't what you say you want. To get the Arm A vs Arm B comparison in each age group, here your best bet is to use -lincom-, because you want the comparisons as relative hazard ratios, rather than difference of hazard ratios. So, assuming that randgroup and age_binary are both coded 0/1:

lincom 1.randgroup, hr // TREATMENT EFFECT WHEN age_binary == 0
lincom 1.randgroup + 1.age_binary#1.randgroup, hr // TREATMENT EFFECT WHEN age_binary == 1

Note that the first result will be exactly the same as the hr for randgroup shown in the -stcox- output table. The other result will be the product of the hr's for randgroup and the age_binary#randgroup interaction term, but you will get the associated standard error, test statistics, and confidence interval as well.

[Megan Moreton]

That's brilliantly helpful - thank you.

One of the subgroups that I want to estimate the treatment group HRs for has 3 levels (coded 1, 2, and 3). Can I check I have the lincom command correct for those 3 analyses:

lincom 1.randgroup, hr
lincom 1.randgroup + 2.site#1.randgroup, hr
lincom 1.randgroup + 3.site#1.randgroup, hr

The full output is:

. stcox randgroup##site, allbase

         Failure _d: death==1
   Analysis time _t: day_eos
  Exit on or before: time 84

Iteration 0:   log likelihood =  -603.6713
Iteration 1:   log likelihood =  -600.6684
Iteration 2:   log likelihood = -600.63867
Iteration 3:   log likelihood = -600.63867
Refining estimates:
Iteration 0:   log likelihood = -600.63867

Cox regression with Breslow method for ties

No. of subjects =    555                                Number of obs =    555
No. of failures =     97
Time at risk    = 40,934
                                                        LR chi2(5)    =   6.07
Log likelihood = -600.63867                             Prob > chi2   = 0.2999

--------------------------------------------------------------------------------
            _t | Haz. ratio   Std. err.      z    P>|z|     [95% conf. interval]
---------------+----------------------------------------------------------------
     randgroup |
        Arm A  |          1  (base)
        Arm B  |   1.324861   .4881331     0.76   0.445      .643506    2.727645
               |
          site |
            1  |          1  (base)
            2  |   .9865466   .3634831    -0.04   0.971     .4791825    2.031114
            3  |   1.894118   .7430822     1.63   0.103     .8779543    4.086413
               |
randgroup#site |
      Arm A#1  |          1  (base)
      Arm A#2  |          1  (base)
      Arm A#3  |          1  (base)
      Arm B#1  |          1  (base)
      Arm B#2  |   1.180754   .5718767     0.34   0.732     .4569778    3.050868
      Arm B#3  |   .7454271   .4055815    -0.54   0.589     .2566094    2.165398
--------------------------------------------------------------------------------

. lincom 1.randgroup, hr

 ( 1)  1.randgroup = 0

------------------------------------------------------------------------------
          _t | Haz. ratio   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         (1) |   1.324861   .4881331     0.76   0.445      .643506    2.727645
------------------------------------------------------------------------------

. lincom 1.randgroup + 2.site#1.randgroup, hr

 ( 1)  1.randgroup + 1.randgroup#2.site = 0

------------------------------------------------------------------------------
          _t | Haz. ratio   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         (1) |   1.564334   .4917911     1.42   0.155     .8447555    2.896863
------------------------------------------------------------------------------

. lincom 1.randgroup + 3.site#1.randgroup, hr

 ( 1)  1.randgroup + 1.randgroup#3.site = 0

------------------------------------------------------------------------------
          _t | Haz. ratio   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
         (1) |   .9875871   .3953919    -0.03   0.975     .4505957    2.164531
------------------------------------------------------------------------------

Thank you very much.