Margins after stpm2

[Jonviea Chamberlain]

Good afternoon,

I would like to have point estimates of 5-, 10-, and 15-year survival based on the model estimated from the stpm2 command.

Example for 5 year survival:

xi: stpm2 i.sex i.age i.decade i.var_x , scale(hazard) df(3) eform

gen t5=5
predict surv_t5 , surv timevar(t5) at(variable_x 1) ci
slist surv_* in 1/1

The issue is that I would like to hold all the other variables at a population average -- i.e., not at baseline values, which is what would happen when using the zeros option (if I understand this correctly). Ideally, I would like to use the margins command, implementing something like:

margins var_x, predict(surv)

When I try this though, I get the error message "prediction is a function of possibly stochastic quantities other than e(b)"

Any help with this would be greatly appreciated!

Post script: I know that my username does not reflect my real name. I have sent a message to get this fixed.

Jonviea Camberlain

[Paul Lambert]

Hi,

stpm2 has a meansurv option. This will calculate marginal survival curves, i.e. predict a survival curve for each subject in you cohort and then average these curves. You can fix selected covariates to take certain values using that at option, e.g. first everyone to be exposed and then everyone to be unexposed.

Note that a marginal survival curve is different to a predicted survival curve for someone who "happens" to have the average value of each covariate (e.g.. when you use stcurve after stcox or streg).

Example code below

webuse brcancer, clear
stset rectime, f(censrec==1) scale(365.24)

range tt 0 6 100
stpm2 hormon x1 x5e x6, df(4) scale(h)
predict ms0, meansurv ci timevar(tt) at(hormon 0)
predict ms1, meansurv ci timevar(tt) at(hormon 1)

line ms0 ms1 tt

If you want contrasts between these standardized (marginal) survival curves then you can use predictnl or use the stpm2_standsurv command.

[Jonviea Chamberlain]

Dear Prof. Lambert,

Thank you for very much for your quick response -- I did eventually realize that the meansurv must be what I was looking for. So, to now get at the 5-, 10-, 15, and 20-year "marginal" survival estimates, I am running the loop below after my model. I hope that this is now correct!

gen t5=5 
gen t10=10
gen t15=15
gen t20=20

forvalues j=5(5)20 {
disp "`j'-year survival"
predict surv_`j' , meansurv timevar(t`j') ci
summ surv_`j'
summ surv_`j'_lci
summ surv_`j'_uci
drop surv_*
 }

[Paul Lambert]

Hi,

That is a bit inefficient. When you use the timevar() option it evaluates the marginal survival (i.e. averaged over all subjects in your dataset) at each value specified in timevar(), so in your approach you are doing the same calculations N times each time you use predict.

So if you just want the marginal survival at 4 time points create a variable of 4 observations as below,

webuse brcancer, clear
stset rectime, f(censrec==1) scale(365.24)

range tt 1 4 4
stpm2 hormon x1 x5e x6, df(4) scale(h)

predict ms0, meansurv ci timevar(tt) at(hormon 0)
predict ms1, meansurv ci timevar(tt) at(hormon 1)

list tt ms0 ms1 in 1/4

Paul

[Jonviea Chamberlain]

Hi Paul,

Thank you for the tip, I just tried it and it works perfectly. I see now how inefficient my original code was!

Best,
Jonviea