Hi everyone,
I am performing a survival analysis comparing outcomes between two treatment groups on a multiple imputed dataset.
To account for imbalance in baseline covariates between groups, I have calculated a propensity score and standardised IPTW.
Note that the covariates used to calculate the propensity score include some, but not all of the covariates used in the final cox proportional hazards model, such that I have;
// generating propensity score
mi estimate, saving(psmiest): logit treatment_var treatment_covar_list
mi predict xb_mi using psmiest
quietly mi xeq: generate propscore = invlogit(xb_mi)
// generating numerator for subsequent iptw stabilisation
mi estimate, saving(ps_nmiest): logit treatment_var
mi predict xb_mi_num using ps_nmiest
quietly mi xeq: generate iptw_num = invlogit(xb_mi_num)
// generating iptw
generate iptw = 1/propscore if treatment_var == 1
replace iptw = (1-(1/propscore)) if treatment_var == 0
// generating stabilised iptw
generate siptw = iptw_num/propscore if treatment_var == 1
replace siptw = (1-iptw_num)/(1-propscore) if treatment_var == 0
// generating stabilised IPTW CPH model
mi stset surv [pweight=siptw], failure(dead)
mi estimate, hr dots: stcox treatment_var outcome_covar_list
As part of my analysis I would like to show an adjusted KM curve stratified by treatment group and compare the two curves using a weighted log-rank test as described by Xie & Liu(1) and cited in a number of other articles.
So far, I think I can generate the KM curves adjusted for the IPTW by:
mi stset surv [pweight=siptw], failure(dead)
sts graph, by(treatment_var)
Questions:
1) I am not sure where the relevance of the sts graph, adjustfor(varlist) option comes in. Any explanation or help very welcome.
2) Furthermore, as I unfortunately do not have a statistics background, I cannot interpret Xie & Liu's equations for a weighted log-rank test (shown below, where Gw is the weighted LR) into Stata code to test the differences between the two curves generated above. If anyone can describe the code that would produce a weighted log-rank test, that would be much appreciated.
where dj =dj1 +dj0 and Yj =Yj1 +Yj0 be the number of deaths and the number at risk in the whole sample at time tj, j=1,....,D
3) Finally, the above curve is adjusted for covariates that predict the treatment assignment, but not the outcome covariates identified as significant in the CPH model. Is there any way to produce a survival curve stratified by treatment grouping that reflects my final IPTW weighted CPH model, ie. also adjusted for the significant outcome_covar_list?
Apologies for any obvious omissions or misunderstanding on my part.
Kind regards,
Kai
(1) Xie J, Liu C (2005) Adjusted Kaplan-Meier estimator and log-rank test with inverse probability of treatment weighting for survival data. Stat Med 24:3089–3110. doi: 10.1002/sim.2174
http://onlinelibrary.wiley.com/doi/10.1002/sim.2174/pdf
I am performing a survival analysis comparing outcomes between two treatment groups on a multiple imputed dataset.
To account for imbalance in baseline covariates between groups, I have calculated a propensity score and standardised IPTW.
Note that the covariates used to calculate the propensity score include some, but not all of the covariates used in the final cox proportional hazards model, such that I have;
// generating propensity score
mi estimate, saving(psmiest): logit treatment_var treatment_covar_list
mi predict xb_mi using psmiest
quietly mi xeq: generate propscore = invlogit(xb_mi)
// generating numerator for subsequent iptw stabilisation
mi estimate, saving(ps_nmiest): logit treatment_var
mi predict xb_mi_num using ps_nmiest
quietly mi xeq: generate iptw_num = invlogit(xb_mi_num)
// generating iptw
generate iptw = 1/propscore if treatment_var == 1
replace iptw = (1-(1/propscore)) if treatment_var == 0
// generating stabilised iptw
generate siptw = iptw_num/propscore if treatment_var == 1
replace siptw = (1-iptw_num)/(1-propscore) if treatment_var == 0
// generating stabilised IPTW CPH model
mi stset surv [pweight=siptw], failure(dead)
mi estimate, hr dots: stcox treatment_var outcome_covar_list
As part of my analysis I would like to show an adjusted KM curve stratified by treatment group and compare the two curves using a weighted log-rank test as described by Xie & Liu(1) and cited in a number of other articles.
So far, I think I can generate the KM curves adjusted for the IPTW by:
mi stset surv [pweight=siptw], failure(dead)
sts graph, by(treatment_var)
Questions:
1) I am not sure where the relevance of the sts graph, adjustfor(varlist) option comes in. Any explanation or help very welcome.
2) Furthermore, as I unfortunately do not have a statistics background, I cannot interpret Xie & Liu's equations for a weighted log-rank test (shown below, where Gw is the weighted LR) into Stata code to test the differences between the two curves generated above. If anyone can describe the code that would produce a weighted log-rank test, that would be much appreciated.
where dj =dj1 +dj0 and Yj =Yj1 +Yj0 be the number of deaths and the number at risk in the whole sample at time tj, j=1,....,D
3) Finally, the above curve is adjusted for covariates that predict the treatment assignment, but not the outcome covariates identified as significant in the CPH model. Is there any way to produce a survival curve stratified by treatment grouping that reflects my final IPTW weighted CPH model, ie. also adjusted for the significant outcome_covar_list?
Apologies for any obvious omissions or misunderstanding on my part.
Kind regards,
Kai
(1) Xie J, Liu C (2005) Adjusted Kaplan-Meier estimator and log-rank test with inverse probability of treatment weighting for survival data. Stat Med 24:3089–3110. doi: 10.1002/sim.2174
http://onlinelibrary.wiley.com/doi/10.1002/sim.2174/pdf
