Plot hazard ration into Kaplan Meier graph

[Dorothea Ekoka Mbassi]

Dear community,

I am running a survival analysis and for now, I was plotting the p-value of log-rank test into my graph:

sts test rand
di chi2tail(`r(df)',`r(chi2)')
local a=round(chi2tail(`r(df)',`r(chi2)'), 0.001)
sts graph, by(rand) ci ciopts(fcolor(%10)) text(.95 518400000 "Log-rank p=`a'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))

However, I would now like to plot the Hazard Ratio as per Cox regression into my graph:

stcox i.rand

If I understand the -stcox- helpfile correctly, the Hazard Ratio is not automatically stored in e(), but should be in r(table). By the looks of it, the Hazard ratio is 'b' in r(table). I tried to code this, but cannot figure out the syntax to specify 'b' (if this is really the HR).

My try:

stcox i.rand
local b=round(stcox(`r(table)'), 0.001)
sts graph, by(rand) ci ciopts(fcolor(%10)) text(.95 518400000 "Cox HR=`b'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))

How do I get the Hazard Ratio to appear in my graph?

Thank You!

[Andrew Musau]

r(table) is a matrix, so see

help matrix subscripting

to see how to refer to elements of a matrix.

webuse stan3, clear
stset
stcox age posttran surg year
mat list r(table)
display r(table)["b", "age"]

Res.:

. mat list r(table)

r(table)[9,4]
               age    posttran     surgery        year
     b   1.0302239   .97872425   .37382783    .8873107
    se   .01432009   .30325967   .16320401   .05980796
     z   2.1421714  -.06940519  -2.2538148  -1.7737928
pvalue   .03217969    .9446671   .02420782   .07609743
    ll   1.0025359   .53322911    .1588759   .77750223
    ul   1.0586766   1.7964157       .8796   1.0126277
    df           .           .           .           .
  crit    1.959964    1.959964    1.959964    1.959964
 eform           1           1           1           1

. 
. display r(table)["b", "age"]
1.0302239

[Dorothea Ekoka Mbassi]

Dear Andrew,

I start understanding better. However, see what happens:

. stcox i.rand

         Failure _d: censored10TBS
   Analysis time _t: (exitdatetime10TBS-origin)
             Origin: time firstdosedatetime
  Enter on or after: time firstdosedatetime

Iteration 0:   log likelihood = -70.601632
Iteration 1:   log likelihood = -70.122713
Iteration 2:   log likelihood =  -70.12271
Refining estimates:
Iteration 0:   log likelihood =  -70.12271

Cox regression with Breslow method for ties

No. of subjects =         29                            Number of obs =     29
No. of failures =         27
Time at risk    = 5336156000
                                                        LR chi2(1)    =   0.96
Log likelihood = -70.12271                              Prob > chi2   = 0.3277

------------------------------------------------------------------------------
          _t | Haz. ratio   Std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
        rand |
    Placebo  |   .6827018    .267299    -0.97   0.330      .316927    1.470628
------------------------------------------------------------------------------

. 
. mat list r(table)

r(table)[9,2]
                4b.          5.
              rand        rand
     b           1   .68270183
    se           .   .26729896
     z           .  -.97488331
pvalue           .   .32961815
    ll           .     .316927
    ul           .   1.4706282
    df           .           .
  crit    1.959964    1.959964
 eform           1           1

. 
. display r(table) ["b", "rand"]
type mismatch
r(109);

. display r(table) ["b", "5"]
type mismatch
r(109);

. display r(table) ["b", "5. rand"]
type mismatch
r(109);

. display r(table) ["b", "5"]
type mismatch
r(109);

. display r(table) ["b", "Placebo"]
type mismatch
r(109);

. display r(table) ["b", "i.rand"]
type mismatch
r(109);

. display r(table) ["b", "rand Placebo"]
type mismatch
r(109);

Do You have an idea what variable that "type mismatch" might refer to?

But adopting Your example, might my code to display HR in graph be:?

stcox age
local b=round(stcox(`r(table)["b", "age"]'), 0.001)
sts graph, by(rand) ci ciopts(fcolor(%10)) text(.95 518400000 "Cox HR=`b'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))

[Andrew Musau]

No spaces and the variable name is "5.rand"

display r(table)["b", "5.rand"]

[Andrew Musau]

But adopting Your example, might my code to display HR in graph be:?

stcox age
local b=round(stcox(`r(table)["b", "age"]'), 0.001)
sts graph, by(rand) ci ciopts(fcolor(%10)) text(.95 518400000 "Cox HR=`b'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))

If you just want to refer to the coefficient, then a local and r(table) are unnecessary. In Stata, you refer directly to the coefficient on age as "_b[age]", the coefficient on female as "_b[female]" and so on. See in particular the option -coeflegend- within

help estimation options

webuse stan3, clear
stset
stcox age
set scheme s1color
sts graph, ci ciopts(fcolor(%10)) text(.7 1000 "Cox HR=`:di %4.3f _b[age]'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))

[Dorothea Ekoka Mbassi]

Dear Andrew,

sts graph, ci ciopts(fcolor(%10)) text(.7 1000 "Cox HR=`:di %4.3f _b[age]'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))

inserts a different value than my HR into the graph, one that I cannot see in my stcox result table.
However,

stcox i.rand
local b=round(r(table)["b", "5.rand"], 0.001)

works and prints my Hazard ratio into the graph with

sts graph, ci ciopts(fcolor(%10)) text(.7 1000 "Cox HR=`b'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))

Thank You once again!

[Andrew Musau]

sts graph, ci ciopts(fcolor(%10)) text(.7 1000 "Cox HR=`:di %4.3f _b[age]'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))

inserts a different value than my HR into the graph, one that I cannot see in my stcox result table.

The hazard ratio is the exponentiated coefficient, so you want

sts graph, ci ciopts(fcolor(%10)) text(.7 1000 "Cox HR=`:di %4.3f exp(_b[age])'") title("") subtitle("") ytitle("") xtitle("") ylabel(0(0.1)1, angle(horizontal))