I am building a Cox PH model in r to analyze my data and I am unsure what the best way to do it is. It is observational data and the idea is that we want to know the effect of a treatment A used at two different hospitals, site 1 and 2. We examine time-to-death and we have a few covariates (binary and continous) that we would like to include as well. However, due to differences in standard of care it is not given that often at site 2 whereas it is more common at site 1. This means that when it is NOT given at site 1 there is something very different going on with the patient. What we want is therefore to compare all patients treated with A (site 1 and 2) with patients not treated with A at site 2 only. So, we want to exclude patients from site 1 not receiving the treatment, because we believe there is another reason for not recieving treatment A. How would you build this model? Should I construct a dummy variable with a binary for TreatmentA vs. Site2_NoTreatmentA and exclude the Site1_NoTreatmentA or should I build the model with a strata or cluster for site?

Thank you in advance.


1 Answer 1


Presumably the overall patient population of interest is the same at both institutions, or there would be no point in combining data between them. If the reasons for treatment choice within that patient population differ between the two institutions, then those reasons need to be included in the model. Simply ignoring individuals not given treatment A at one institution would risk substantial bias in any attempt to apply results to the overall population of interest.

Any model would need to include as many outcome-associated covariates as is reasonable without overfitting. "A few covariates" probably won't be adequate to capture the reasons for death beyond those associated specifically with treatment, so an estimate of treatment effect based on just a few covariates would be unreliable. You would also seem to need a treatment-by-institution interaction, whether you include the institution as a modeled predictor or stratify by institution. Covariates might be included directly in the model or via probability-of-treatment weighting. Counterfactual modeling would probably be needed to identify the true treatment effect from such observational data. See this page and its links for an introduction to the issues and to the associated difficulties in causal inference.

This situation is so tricky that you need to consult with a local statistician with whom you can discuss the intricacies of the data and work through which particular modeling strategy will best test your hypothesis of interest.

  • $\begingroup$ Thank you for taking the time to answer. If we just consider the problem with having two sites. If I include a cluster(site) It will be the same as including the dummy variable as interaction term, right? $\endgroup$
    – User LG
    Oct 10, 2022 at 5:48
  • 1
    $\begingroup$ @UserLG that's not correct. A cluster term doesn't affect coefficient estimates, it only adjusts standard error estimates, and it doesn't allow for different treatment effects between sites. A strata(site) term with a treatment interaction or a site dummy coefficient with a treatment interaction is needed to allow for different treatment effects between sites. $\endgroup$
    – EdM
    Oct 10, 2022 at 13:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.