Bias and Sensitivity Analysis When Estimating Treatment Effects from the Cox Model with Omitted Covariates

Summary Omission of relevant covariates can lead to bias when estimating treatment or exposure effects from survival data in both randomized controlled trials and observational studies. This paper presents a general approach to assessing bias when covariates are omitted from the Cox model. The proposed method is applicable to both randomized and non-randomized studies. We distinguish between the effects of three possible sources of bias: omission of a balanced covariate, data censoring and unmeasured confounding. Asymptotic formulae for determining the bias are derived from the large sample properties of the maximum likelihood estimator. A simulation study is used to demonstrate the validity of the bias formulae and to characterize the influence of the different sources of bias. It is shown that the bias converges to fixed limits as the effect of the omitted covariate increases, irrespective of the degree of confounding. The bias formulae are used as the basis for developing a new method of sensitivity analysis to assess the impact of omitted covariates on estimates of treatment or exposure effects. In simulation studies, the proposed method gave unbiased treatment estimates and confidence intervals with good coverage when the true sensitivity parameters were known. We describe application of the method to a randomized controlled trial and a non-randomized study.

Using theorems (2.2), (2.4) and (2.6) in Jiang (2010, page 22-23), the ratio Añ ik /Bñ i converges in probability to Under the assumption T ⊥T + |X, is the survival function conditional on x, H 0 (·) is the cumulative baseline hazard, and S + (T i |x) is the survival function of censoring time given x. It follows that E xc e Xθ * e −H 0 (T i )e Xθ+Cβ S + (T i |X) .
In the special case when S + (T i |X) = 1 −T i /τ X for 0 ≤T i ≤ τ X , D ik will be an indeterminate form 0/0 asT i approaches τ X . In this case, l'Hospital's 2 rule can be used: Other similar cases can be handled by algebraic elimination, L'Hopital's rule, or other methods to manipulate the expression. More generally, 1 .
It can be shown using the definition of convergence in probability that 1 Writing the conditional expectation E (·|T ≤ T + ) as E OBS (·) and assuming P (T ≤ T + ) = 0 (i.e. the probability of observing an event is not zero), Uñ k (θ * ) converges in probability to U k (θ * , θ; β) = 0 = E OBS (X k − D k ), which leads to (4).

Web Appendix B: Proof for (7)
Using the simulation setting in Lin et al (1998) that S + (t|x, c) = 1 − t/τ and K = 1, the equation (4) can be written as To be consistent with Lin et al (1998), we assume the event is sufficiently rare and H 0 (t) sufficiently small that e H 0 (t)e xθ+cβ ≈ 1 (i.e. O(H 0 (t)e xθ+cβ ) can be ignored). The left side of (A-1) can then be approximated by The uncensoring probability can be written as When H 0 (t) is small, (A-2) is approximately Using (6), the right side of (A-1) is therefore and it follows The results (A-4) and (A-5) are the same as the formulae in Lin et al (1998).

Web Appendix C: Proof for (11)
The equation (4) shows that For X ∼ B(1, p), by the result (6), the right side of (A-6) is where ϕ 0 and ϕ 1 are the uncensoring rates in control and treatment groups respectively. The left side of (A-6) is These lead to (9).

5
Web figure 1: Kaplan-Meier plots for complete data S(t|x, c) (left) and the observed data S(t|x) (right). The sample size is 10,000. The data were generated from h(t|x, c) = e x+βc , X ∼ B(1, 0.5) and C ∼ B(1, 0.5). Under this setting, S(t|x) = 0.5S(t|x, c = 1) + 0.5S(t|x, c = 0). As β increases from 0.5 to 5, on the left side, S(t|x, c = 1) tend to zero, reducing the observed difference between S(t|x = 1) and S(t|x = 0) (right). As a consequence, the marginal hazard ratio exp(θ * ) = h(t|x = 1)/h(t|x = 0), attenuates to 1 and the bias increases from -0.14 to -0.45. But as β increases from 5 to 8, there is no apparent change in the difference S(t|x = 1) and S(t|x = 0). As a result, the bias does not increase with β and -0.45 is the limit of bias. 1−e −e γ 0 +(γ 1 +θ)x+βc e γ 0 +(γ 1 +θ)x+βc for x = 0 and 1. It can be seen that γ 1 should be −θ to ensure that the censoring probabilities are the same in the two treatment groups. The overall censoring probability is controlled by the value of γ 0 . 9 Web figure 5: Kaplan-Meier plots for complete data S(t|x, c) (left) and the observed data S(t|x) (right). The sample size is 200,000. The data were generated from h(t|x, c) = e x+5c , X ∼ B(1, 0.5) and C ∼ B(1, 0.5). S(t|x) = 0.5S(t|x, c = 1) + 0.5S(t|x, c = 0). For ease of simulation, we generated t + from uniform (0, τ ), where τ was solved from (13) such that the overall proportions censored could be 0%, 50% or 90%. When 0% data are censored, S(t|x, c = 1) gets close to zero for t < 0.05 and so most of the information about θ (i.e. the difference between S(t|x = 1) and S(t|x = 0)) is supplied by the difference between S(t|x = 1, c = 0) and S(t|x = 0, c = 0). When 50% data are censored, the times with c = 0 are more likely to be censored than those with c = 1. As a consequence, the difference between S(t|x = 1, c = 0) and S(t|x = 0, c = 0) cannot supply much information about θ and so the bias increases. As the censoring percentage increases from 50% to 90%, the difference between S(t|x = 1, c = 0) and S(t|x = 0, c = 0) becomes less informative. At the same time, the difference between S(t|x = 1, c = 1) and S(t|x = 0, c = 1) starts to supply more of the information about θ and the bias decreases. When 90% of the data are censored, almost all of the few deaths observed will occur in subjects with the same value of c (c = 1) and the bias is very low. Web figure 6: Illustrative example showing density plots of 100,000 bootstrap samples ofθ * and the correspondingθ = R(θ * ). A sample of size 200 was generated from h(t|x, c) = e x+3c , X ∼ B(1, 0.5), C|x ∼ B(1, 0.3 + 0.4x) and T + ∼ unif orm(0, 1). In this simulation, the distribution ofθ = R(θ * ) is slightly right skewed. It will be not accurate to use the standard error of the bootstrap sample ofθ = R(θ * ) to construct the C.I. Web figure 7: Sensitivity analysis for the p-value for the treatment effect of folinic acid on age of sitting for children with Down's syndrome. Web