I've used the twang package to calculate propensity score weights in my observational data set. The twang vignette uses the survey package to apply the weighted PS in analyses, but does not provide an example for Cox proportional hazard estimation.
The survey package provides a function called svycoxph() which which seems to do the trick. However, when i use coxph() from the survival package an set the weight= argument to the extracted weighted PS from twang's ps() object i get an indentical HR estimate, albeit a much more narrow 95%CI.
Can I use the extracted propensity score weights from a ps object from the twang package as a vector in the weight= argument of the survival package's coxph()?
Should i include an example?
Update: When i added a cluster(patient ID) term in the the coxph formula, both the HR and the 95% CI became identical to the syvcoxph output. So I guess that did the trick. But can i conclude that this is ok to use?
Best Answer
I just came across your post because I used twang to calculate inverse probability weights and wanted to use a coxph model afterwards. I had used svyglm before, but wasn't aware of svycoxph. So thanks for this hint.
I tried what you did and got almost identical results with coxph and cluster(patient ID) as with svycoxph.
In the description of svycoxph, they say:
In the description of coxph, they give the following explanation:
From these two explanations I conclude that coxph using cluster(id) is fine for large numbers of clusters (what you probably had). If you have a small number of clusters and want to account for this, svycoxph will be safer. In my case, with 67 clusters, it made a small difference so I will use svycoxph.