Dear Statalisters,
I would like to automate the calculation of the integrated
discrimination improvement (or IDI) from several nested models, in
order to evaluate the incremental predictive ability of two
biomarkers.
Model 1 = clinical factors
Model 2 = clinical factors + BNP
Model 3 = clinical factors + Ca125
Model 4 = clinical factors + BNP + Ca125
Model 5 = clinical factors + BNP##Ca125
The formula to calculate the IDI is:
Absolute IDI= (p_event2-p_event1)+(p_nonevent1-p_nonevent2),
Relative IDI= [(p_event2-p_nonevent2)/(p_event1-p_nonevent1)]-1
The probabilities for events & nonevents are calculated as follow:
Model 1= logistic death X1 X2 X3
predict double m6m_0, pr
Model 2= logistic death X1 X2 X3 BNP
predict double m6m_bnp, pr
Model 3= logistic death X1 X2 X3 Ca125
predict double m6m_ca, pr
Model 4= logistic death X1 X2 X3 BNP Ca125
predict double m6m_bnpca, pr
Model 5= logistic death X1 X2 X3 BNP##Ca125
predict double m6mbnpca, pr
For Model 2 vs Model 1:
sum m6m_0 if death==1
(0.3197818)
sum m6m_bnp if death==1
(0.336075)
sum m6m_0 if death==0
(0.1323866)
sum m6m_bnp if death==0
(0.1292155)
Absolute IDI .....display (.336075-.3197818)+(.1323866-.1292155)
Relative IDI.......display ((.336075-.1292155)/(.3197818-.1323866))-1
Then, I tried the following code (for relative IDI) with error:
---------------------------------------------------------------------------------------------------------------
tempname hdle
capt erase info.dta
postfile `hdle' str10 idi2_1 idi3_1 idi4_1 idi5_1 /*
*/ using info
foreach var of varlist m6m_0 m6m_bnp m6m_ca m6m_bnpca m6mbnpca {
sum `var' if !death
loc m6m_0_mean0=r(mean)
loc m6m_bnp_mean0=r(mean)
loc m6m_ca_mean0=r(mean)
loc m6m_bnpca_mean0=r(mean)
loc m6mbnpca_mean0=r(mean)
sum `var' if death
loc m6m_0_mean1=r(mean)
loc m6m_bnp_mean1=r(mean)
loc m6m_ca_mean1=r(mean)
loc m6m_bnpca_mean1=r(mean)
loc m6mbnpca_mean1=r(mean)
gen
idi2_1=((m6m_bnp_mean1-m6m_bnp_mean0)/(m6m_0_mean1-m6m_0_mean0))-1
gen idi3_1=((m6m_ca_mean1-m6m_ca_mean0)/(m6m_0_mean1-m6m_0_mean0))-1
gen idi4_1=((m6m_bnpca_mean1-m6m_bnpca_mean0)/(m6m_0_mean1-m6m_0_mean0))-1
gen idi5_1=((m6mbnpca_mean1-m6mbnpca_mean0)/(m6m_0_mean1-m6m_0_mean0))-1
post `hdle' (`idi2_1') (`idi3_1') (`idi4_1') (`idi5_1')
}
postclose `hdle'
use info, clear
l, abbrev(10) noobs sep(0)
---------------------------------------------------------------------------------------------------------------
Can someone help me with this?
Warm regards,
Eduardo
*
* For searches and help try:
* http://www.stata.com/help.cgi?search
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
