How to get c-index and calibration plot for survival models on imputed data?

[Markos Valsamis]

I have been developing survival models using stcox and stpm2 commands and have been using estat concordance and stcstat2 commands to obtain the c-index of apparent model performance, respectively, as post-estimation commands.

However, I am unable to obtain the c-index or calibration plot if I develop the models while also doing multiple imputation, as the above commands do not work.

For instance:

stset spell, failure(event==1) scale(365)
mi set flong
mi register imputed missing_variable
mi impute chained (mlogit) missing_variable = v1 v2 v3 v4, add(10) noisily rseed(123)
mi estimate, saving(my_estimates, replace): stcox v1 v2 v3 v4
estat concordance

or with stepwise selection:

stset spell, failure(event==1) scale(365)
mi set flong
mi register imputed missing_variable
mi impute chained (mlogit) missing_variable = v1 v2 v3 v4, add(10) noisily rseed(123)
xi: mim: stepwise, pr(0.1): stcox v1 v2 v3 v4, nohr
estat concordance

I presume this may need to be done manually if there is no package that can achieve, and would be grateful for your help in how I would code this?


Many thanks in advance

[Rich Goldstein]

see #6 in https://www.statalist.org/forums/for...ng-mi-estimate for a related situation for the c-statistic;

I am a little nervous about your saying, "the c-index or calibration plot" as these are not at all the same things; further, a calibration plot can be done after mi by using the appropriate version on the predict command and thus uses a different logic

[Markos Valsamis]

Many thanks for your reply Rich. I tried running the below code on my survival data but Stata doesn't seem to be registering the list of covariates under local rhs so does not recognise them in the second line (error: "varlist specification required"). Even when I add them manually, the same problem arises for r(M). Would you be able to advise?

local rhs "ageatprimary i.patientgender"
noi mi estimate, saving(my_estimates, replace); stcox `rhs'
qui mi query
local M= r(M)
scalar cstat=0
qui mi xeq 1/`M': stcox `rhs'; estat concordance; scalar cstat=r(C)
scalar cstat=cstat/`M'
noi di "C statistic over imputed data=" cstat

Indeed, a calibration plot is very different to obtaining a c-statistic- what I meant to say is that I cannot get the c-statistic and also I cannot create a calibration plot either.

[Markos Valsamis]

Please ignore my comment above- I realised my error. The below code works, thanks!

local model = "ageatprimary"
noi mi estimate, saving(my_estimates, replace): stcox `model'
qui mi query
local M=r(M)
scalar cstat=0
qui mi xeq 1/`M': stcox `model'; estat concordance; scalar cstat=cstat+r(C)
scalar cstat=cstat/`M'
noi di "C statistic over imputed data= " cstat

Would you be able to advise how I could go about getting a calibration plot following mi?

[Rich Goldstein]

not sure exactly what you mean by "calibration plot" but this is what I do - after using -mi estimate-, I use -mi predict- which gives me predicted values in only the non-imputed part of the data (data set 0) and I can then use -lowess- in the standard way (see, e.g., Austin, PC and Steyerber, EW (2014), "Graphical assessment of internal and external calibration of logistic regression models by using loess smoothers," Statistics in Medicine, 33: 517-535; note that there may be complications with survival data because of the censoring but see Austin, PC, Harrell, FE, Jr. and van Klaveren, P (2020), "Graphical calibration curves and integrated calibration index (ICI) for survival models," Statistics in Medicine, 39: 2714-2742