Survival curves after multiple imputation

[Hans Pinnschmidt]

After having done multiple imputation on some covariates before running stcox, I'd like to plot the survival curves for the individual levels of the covariates "averaged" across imputed data sets.
Example:
xi:mim:stcox i.group
How do I now obtain such "averaged" curves for the various levels of covariate "group"? Or, if averaging across imputed data sets is not possible, can I plot curves for each individual imputed data set? I checked all posts that I could find on that issue - to no avail. Command "stcurve" yields the message "last estimates not found"...
Some advice on this would be highly appreciated.

[Steve Samuels]

Sorry, but I don't know mim, as its functionality has been superseded by Stata's mi commands in recent versions. Here's how I would solve your problem in Stata 13. You'll have to translate as best you can. Basically, I extract each imputed data set, estimate the stcox model on that extract, run stcurve, and save the results. At the end I collapse the estimated survival curves over extracts and plot the resulting averages.

clear
save results, emptyok replace // file to hold results
sysuse auto, clear
set seed 43228226
recode rep78 1/2=3
gen id = _n
rename id pid
mi set flong
mi register impute rep78

mi impute ologit rep78  mpg, add(5) double

mi stset mpg, fail(foreign)
mi estimate, post: stcox headroom i.rep78
tempfile d0
save `d0', replace

tempfile d1 d2 d3 d4 d5
forvalues i = 1/5{
  use d0, clear
  mi extract `i' /* get i-th imputed data set */
  qui stcox headroom i.rep78
  stcurve , surv at1(headroom = 5) at2(headroom = 10) ///
  outfile(`d`i'', replace)
  use results, clear
  append using `d`i''
  save, replace
}

use results, clear
collapse (mean) surv2 (mean) surv3, by(_t)
sort _t
twoway scatter surv2 _t, c(stairstep) ms(i) ///
 || scatter surv3 _t, c(stairstep) ms(i) ///
  ti("Averaged Curves")

[Hans Pinnschmidt]

After a few trials & errors, it seems to work just fine. Thanks a lot, Steve!

[Steve Samuels]

You are very welcome, Hans. Note that the code could have been made automatic for any number of imputations by(1) putting the number of imputed data sets into a local macro ("isets", below) and (2). Defining the temporary data sets d1, d2,... within each iteration of the forvalues loop, instead of outside the loop. Here is the improved version.

clear
save results, emptyok replace // file to hold results
sysuse auto, clear
set seed 43228226
recode rep78 1/2=3
gen id = _n
rename id pid
mi set flong
mi register impute rep78


mi impute ologit rep78  mpg, add(5) double
local isets = r(M) /* number of imputed data sets */
mi stset mpg, fail(foreign)
return list
mi estimate: stcox headroom i.rep78
return list


tempfile d0
save `d0', replace
qui{
forvalues i = 1/`isets'{
  tempfile d`i'
  use `d0', clear
  mi extract `i' /* get i-th imputed data set */
  qui stcox headroom
  stcurve , surv at1(headroom = 5) at2(headroom = 10) ///
  outfile(`d`i'', replace)
  use results, clear
  append using `d`i''
  save, replace
}
}

use results, clear
collapse (mean) surv2 (mean) surv3, by(_t)
labe var surv2 "headroom 5" surv3 "headroom 10"
sort _t
twoway scatter surv2 _t, c(stairstep) ms(i) ///
 || scatter surv3 _t, c(stairstep) ms(i) ///
  ti("Averaged Curves")

[Hans Pinnschmidt]

Hi Steven,

your code runs just fine as is. However, when I try to implement my own “needs” (for instance, I want to read already-imputed datasets [10 imputations] from my data file location as indicated in your code below), two things go wrong (I have used “.set trace on” to track this):

1. Invalid syntax error message (r198):

. qui{
- forvalues i = 1/`isets'{
= forvalues i = 1/{

So I changed “forvalues i = 1/`isets'” into “forvalues i = 1/10{“ and this error disappears but now I got another invalid file specification error r(198):
. qui{
- forvalues i = 1/10{
- tempfile d`i'
= tempfile d1
- use `d0', clear
= use , clear
invalid file specification

So I changed “use `d0', clear” into “use d0, clear” and everything seems to work. However, I’m aware of the fact that now I’ve lost some nice features of your code as I now have to type the number of imputations, instead of the program determining that by itself. Can you see what’s wrong here?

Again, many thanks for your help!

Hans

***** Code *****
cd "<my data file path>”
clear
save results, emptyok replace // file to hold results
use "<my data file>", clear

local isets = r(M) /* number of imputed data sets */
mi stset beobachtungszeit2, failure(status2==1)
return list
mi estimate: stcox i.hyperlipidamie
*mi estimate, post: stcox i.hyperlipidamie
return list

tempfile d0
save `d0', replace

*set trace on
qui{
*forvalues i = 1/`isets'{
forvalues i = 1/10{
tempfile d`i'
*use `d0', clear
use d0, clear
mi extract `i' /* get i-th imputed data set */
qui stcox i.hyperlipidamie
stcurve , surv at1(hyperlipidamie = 0) at2(hyperlipidamie = 1) ///
outfile(`d`i'', replace)
use results, clear
append using `d`i''
save, replace
}
}

use results, clear
collapse (mean) surv2 (mean) surv3, by(_t)
label variable surv2 "No hyperlipidemia"
label variable surv3 "Hyperlipidemia"
sort _t
twoway scatter surv2 _t, c(stairstep) ms(i) ///
|| scatter surv3 _t, c(stairstep) ms(i) ///
ti("Survival curves averaged over imputed data sets")

[Steve Samuels]

Patrick,
r(M)is a macro left-behind by mi impute, and so would not be available for defining isets at the start of your code. My suggestion is to define isets directly before the quietly block and the forvalues loop:

local isets = 10
quietly {
forvalues i = 1/`isets' {
...

[Mareika Minder]

Hi Steve, thanks very much for this Syntax! It helps me a lot!!!! one question: if i want to draw more than 2 curves I write

stcurve , surv at1(headroom = 1) at2(headroom = 2) at3(headroom = 3) at4(headroom = 4)

But how do I continue? I dont understand how I have to change the second part of the Syntax:
use results, clear
collapse (mean) surv2 (mean) surv3, by(_t)
labe var surv2 "headroom 5" surv3 "headroom 10"
sort _t
twoway scatter surv2 _t, c(stairstep) ms(i) ///
|| scatter surv3 _t, c(stairstep) ms(i) ///
ti("Averaged Curves")

Thanks a lot!

[Saeid Bitaraf]

hi, how can i have baseline hazard after multiple imputation