A flexible approach for causal inference with multiple treatments and clustered survival outcomes

Abstract

When drawing causal inferences about the effects of multiple treatments on clustered survival outcomes using observational data, we need to address implications of the multilevel data structure, multiple treatments, censoring, and unmeasured confounding for causal analyses. Few off-the-shelf causal inference tools are available to simultaneously tackle these issues. We develop a flexible random-intercept accelerated failure time model, in which we use Bayesian additive regression trees to capture arbitrarily complex relationships between censored survival times and pre-treatment covariates and use the random intercepts to capture cluster-specific main effects. We develop an efficient Markov chain Monte Carlo algorithm to draw posterior inferences about the population survival effects of multiple treatments and examine the variability in cluster-level effects. We further propose an interpretable sensitivity analysis approach to evaluate the sensitivity of drawn causal inferences about treatment effect to the potential magnitude of departure from the causal assumption of no unmeasured confounding. Expansive simulations empirically validate and demonstrate good practical operating characteristics of our proposed methods. Applying the proposed methods to a dataset on older high-risk localized prostate cancer patients drawn from the National Cancer Database, we evaluate the comparative effects of three treatment approaches on patient survival, and assess the ramifications of potential unmeasured confounding. The methods developed in this work are readily available in the $\text{R}$ package $\text{riAFTBART}$ .

Keywords: Bayesian machine learning; multilevel survival data; observational studies; sensitivity analysis.

FIGURE 1

Relative biases among 250 replications for each of four methods, IPW‐riCox, DR‐riAH, PEAMM, and riAFT‐BART, and three treatment effects $CAT{E}_{1,2}$, $CAT{E}_{1,3}$, and $CAT{E}_{2,3}$ based on 5‐year RMST under four data configurations: (proportional hazards vs nonproportional hazards) $\times$ (10% censoring proportion vs 40% censoring proportion). The true treatment effects under proportional hazards are $CAT{E}_{1,2}^{0,PH}=7.7$ months, $CAT{E}_{1,3}^{0,PH}=3.6$ months and $CAT{E}_{2,3}^{0,PH}=-4.1$ months. The true treatment effects under nonproportional hazards are $CAT{E}_{1,2}^{0,nPH}=8.1$ months, $CAT{E}_{1,3}^{0,nPH}=3.9$ months, and $CAT{E}_{2,3}^{0,nPH}=-4.2$ months

FIGURE 2

Relative biases in the estimates of three pairwise treatment effects ${\text{CATE}}_{1,2}$, ${\text{CATE}}_{1,3}$, and ${\text{CATE}}_{2,3}$ among 1000 replications using data simulated for illustrative sensitivity analysis. Three causal analyses were performed: (i) including unmeasured confounder (including UMC), (ii) sensitivity analysis (SA) using priors of increasing width for the confounding function centered around the true confounding function, $\text{Unif}(c^0-\omega \widehat{\sigma },c^0+\omega \widehat{\sigma })$, where $\omega =0,0.5,1$, and (iii) ignoring unmeasured confounder (ignoring UMC)

FIGURE 3

The posterior mean of the counterfactual survival curves for each of three treatment groups in NCDB data. The solid curves are the average by treatment group of the individual‐specific survival curves estimated following equation (7). The dashed survival curves are the Kaplan‐Meier estimates for each treatment group. Solid gray curves are estimates of individual‐specific survival curves for 15 randomly selected patients from three treatment groups

FIGURE 4

The location effect in terms of the log survival time in months represented by the posterior mean and credible intervals of the random intercept $b_k$, $k=1,\dots ,9$

FIGURE 5

The Kaplan‐Meier estimators for each of three treatment groups superimposed on the individual‐specific survival curves predicted from the riAFT‐BART model fitted on the NCDB data. Solid black survival curves are the Kaplan‐Meier estimators. Solid grey curves are predicted individual‐specific survival curves for 50 patients randomly selected from each treatment group. Panel A: radical prostatectomy (RP) group; Panel B: external beam radiotherapy combined with androgen deprivation (EBRT+AD) group; Panel C: external beam radiotherapy plus brachytherapy with or without androgen deprivation (EBRT+brachy$\pm$AD) group