Doubly robust omnibus sensitivity analysis of externally controlled trials with intercurrent events

ABSTRACT

1 INTRODUCTION

1.1 Why consider external controls

In medical research, the gold standard for evaluating new treatments has been randomized controlled trials (RCTs). Regulatory bodies often require these rigorously controlled clinical studies to validate the effectiveness and safety of these treatments for specific patient groups. Despite the high regard for randomized, double-blind trials, they may not always be viable or ethical, particularly in rare or severe diseases with limited treatment options. Recruiting enough participants for such trials is often difficult, and using placebos in situations when alternative effective treatments are available may be unethical or impractical.

In these cases, single-arm trials (SATs) can be a practical alternative, though they come with limitations. Without direct comparison data for untreated subjects, researchers are compelled to infer these outcomes from external sources such as previous studies or real-world databases. These are referred to as externally controlled trials (ECTs). However, data from past trials might not always be relevant due to changes in the treatment landscape and patient demographics, reducing its evidence capacity compared with RCTs. Recently, real-world data have gained popularity for external control arms due to its accessibility and contemporaneity with treatment groups.

1.2 Estimand considerations with external controls

The Food and Drug Administration’s (FDA) latest draft guidelines for natural history studies in rare disease drug development emphasized 5 critical concerns with using external controls (Food and Drug Administration, 2023). These concerns, outlined in Table 1, may introduce biases in research findings when real-world data are utilized as a control arm in ECTs. For a more structured analysis of these biases, we classified them into 2 principal categories: baseline incomparability and outcome mean non-exchangeability. Each category is characterized by unique mechanisms for introducing bias, as detailed in Table 1.

Besides the challenges of non-exchangeability of the external controls in clinical trials, managing intercurrent events such as participant dropout, non-compliance, or premature termination of therapy adds more complexities. Typically, the missingness at random (MAR) framework is adopted after the absence of outcomes following these events. While the MAR assumption is often deemed plausible, it remains unverifiable and may lack practical applicability since it assumes all participants persist with the study medication, addressing only a theoretical aspect of the treatment effect, as noted in the guidelines (ICH, 2021). A more plausible assumption would be that the treatment effect may quickly fade away, leading to the missing not at random (MNAR) pattern for the intercurrent events. The MNAR assumption can also be put forward for the sensitivity analysis, suggesting that any treatment effect observed, while a participant was active in the study is nullified upon their discontinuation. To evaluate this effect, the implementation of control-based imputation (CBI) models has been proposed, offering a nuanced and realistic assessment of treatment outcomes in clinical trials.

1.3 Primary analysis and doubly robust omnibus sensitivity analysis

Numerous statistical methodologies have been developed to address the potential biases from using external controls in ECTs. Most of these methods utilize techniques such as propensity score stratification, matching, or weighting to mitigate selection bias by balancing baseline covariates between external controls and SATs. Nonetheless, the outcome exchangeability of external controls cannot be verified with observed data due to the absence of concurrent controls. Therefore, sensitivity analyses become essential to evaluate the impact of potential violations of outcome non-exchangeability and to understand the effects of intercurrent events under different realistic scenarios.

In this context, we introduce a sensitivity analysis framework tailored for ECTs with intercurrent events. This framework encompasses models such as the tilting models, which control the degree of outcome mean non-exchangeability in external controls and alterations in outcome distributions after intercurrent events. Central to our framework is the establishment of identification results, the development of efficient influence functions (EIFs), and EIF-motivated tilting estimators under sensitivity models, which jointly capture the effect of changes due to outcome mean non-exchangeability and intercurrent events. These EIF-motivated estimators have several advantageous statistical properties, such as local efficiency, double robustness, and asymptotic normality. By analytically establishing the conditions for desirable asymptotic properties, our estimator allows flexible models for nuisance parameters while maintaining root-n consistency (Theorem 4). Therefore, our major contribution is to derive the locally efficient estimator for evaluating the treatment effect under the sensitivity models and to jointly assess the robustness and reliability of multiple assumptions in ECTs with intercurrent events in a more efficient and flexible manner.

Our paper is organized as follows: Section 2 presents a brief overview of sensitivity analysis. Section 3 introduces the notation and develops the semi-parametric efficient estimator for the primary analysis. Section 4 details the tilting sensitivity models and the efficient inference for the sensitivity analysis. Another efficient estimator for the CBI sensitivity analysis under the jump-to-reference (J2R) model is presented in the Supplementary Materials. Section 5 discusses one practical method for choosing the sensitivity parameters by bounding their impacts. Extensive simulation studies for both primary and sensitivity analyses are presented in Section 6. Section 7 illustrates our approach with an antidepressant trial, and Section 8 concludes with a discussion.

2 RELATED WORKS

Before we delve into the proposed sensitivity analysis framework, we provide a review of sensitivity analysis methods. Causal inference involving observational studies usually requires no unmeasured confounding assumption, that is, the treatment assignment is ignorable conditional on a set of covariates (Rosenbaum and Rubin, 1983b; Imbens and Rubin, 2015). Yet, claiming the absence of confounders in the treatment-outcome relationship is untestable and often implausible in practice. Thus, it is advised to conduct a series of sensitivity analyses assessing how robust the causal findings are against the plausible violations of the unconfoundedness (Faries et al., 2024). The problem of sensitivity analyses has been studied in a variety of fields with the earliest work in Cornfield et al. (1959), which is later extended in Imbens (2003), Rosenbaum (1987), and Rosenbaum and Rubin (1983a). However, one concerning issue regarding this framework is that it demands a specific parametric assumption on the unmeasured confounder U. In many cases, sensitivity analyses can be quite fragile against the model misspecification of U as shown in Zhang and Tchetgen (2019).

To circumvent modeling the conditional (or marginal) distribution for the unmeasured confounder, a plethora of sensitivity approaches have been proposed, which preserve the critical elements in Imbens (2003) and Rosenbaum and Rubin (1983a) and allow for flexible strategies to model the distribution of U. Zhang and Tchetgen (2019) leverage the modern semi-parametric theory to obtain a consistent estimation in a model while placing no distributional assumptions on U. The idea of partial misspecification of the nuisance parameters endows the framework with an unrestricted law of latent variables and has been similarly considered in the context of measurement error (Tsiatis and Ma, 2004), mixed models (Garcia and Ma, 2016), and statistical genetics (Allen et al., 2005). Another line of work to tackle this problem is based on the idea of “omitted-variable” bias, which can be computed easily without needs to specify the parametric form of the potential unobserved confounding. The idea of general bias formula is introduced in VanderWeele and Arah (2011), and is further extended in Cinelli and Hazlett (2020) with more flexibility and robustness.

As illustrated so far, most of the work blurs the line between sensitivity analysis and model checking by introducing sensitivity parameters under stringent parametric assumptions. For example, Little and Rubin (2019) suggest that the ignorability assumption in Heckman (1979) can be tested as the result of their Gaussian parametric assumption, thereby inducing testable implications of the untestable ignorability assumption. Moreover, the lack of such “clean” separation requires that the observed model to be refit for each setting of the sensitivity parameters, which will be an onerous task as modern non-parametric models may be adopted to fit the potential outcomes. To address this concern, Robins et al. (2000) propose and extended by Franks et al. (2020) and Nabi et al. (2024) to use the “tilting,” or “selection” function to decouple the sensitivity analyses from the observed data model. Typically, such sensitivity analysis specification does not impose any parametric assumptions on the distribution of the observed data or the unmeasured confounder, but only on a relaxed version of the unconfoundedness assumption:

where |$Y(a)$| is the potential outcome under treatment a, A is the treatment assignment, and X is the baseline covariates. Here, the first term constitutes the observed data density, while the second term is the selection function governed by the sensitivity parameter |$\psi$|⁠. Other approaches in Blackwell (2014) and Yang and Lok (2018) take a different track by representing the confounding as a function of the observed covariates, and describe the conditional potential outcomes difference varied by the treatment assignment as |$q(a,X;\alpha )= {\mathbb {E}} \lbrace Y(a)\mid A=a,X\rbrace - {\mathbb {E}} \lbrace Y(a)\mid A=1-a,X\rbrace$|⁠, where the confounding function q is characterized by the single sensitivity parameter |$\alpha$|⁠. As the observed data distribution is free of the sensitivity parameters (⁠|$\psi$| or |$\alpha$|⁠), it achieves the “clean” factorization of the identified and unidentified parts of the sensitivity analysis framework. Veitch and Zaveri (2020), extending from Imbens (2003), posit a probabilistic model to bypass the need to specify any distributional assumption on U, which also decouples the sensitivity analyses from the observed data and leads to tractable bias calculations.

3 SINGLE-ARM TRIAL DATA WITH EXTERNAL CONTROLS: PRIMARY ANALYSIS

To ground ideas, we first focus on cross-sectional studies. Let |$S=1$| denote trial participation, and the trial data consist of |$\lbrace V_{i}=(X_{i},A_{i},R_{i},Y_{i},S_{i}=1):i\in \mathcal {R}\rbrace$|⁠, where |$R=1$| indicates the absence of intercurrent events (eg, treatment discontinuation) during follow-up and 0 otherwise. Let |$S=0$| denote the external participation, and the external data consist of |$\lbrace V_{i}=(X_{i},A_{i},R_{i},Y_{i}, S_{i}=0):i\in \mathcal {E}\rbrace$|⁠. Assume |$(X_{i},A_{i},R_{i},Y_{i},S_{i})$| are independent and identically distributed, and we omit the subscript i for the simplicity of notation. Let |$V=(X,A,R,Y,S)$| be the random vector of all observed variables and follow the observed data distribution |$P_{0}$|⁠. The treatment policy strategy designates a treatment effect estimand that measures the total effect of the treatment assignment and the intercurrent event on the outcome. Therefore, to define the estimand unambiguously, we extend the causal framework in Lipkovich et al. (2020) and introduce the potential outcomes framework for R and Y. Denote |$R(a)$| as the potential indicator for the absence of intercurrent events under treatment a, |$Y(a, r)$| as the potential outcome under treatment a with intercurrent event status r, and |$Y(a) = Y\lbrace a, R(a)\rbrace$|⁠.

The estimand of interest is the average treatment effect (ATE) for the trial, defined as

which is a treatment policy estimand, as the occurrence of intercurrent events is considered irrelevant. However, due to the missing outcomes following these events, the treatment policy strategy cannot be implemented for intercurrent events that are terminal, such as treatment discontinuation. Table 2(A) outlines several key causal assumptions to identify the treatment effect for the primary analysis.

Assumptions 1 and 2 are standard causal assumptions for identification (Rosenbaum and Rubin, 1983b). Assumption 1 maps the potential outcomes to the observed data, and Assumption 2 ensures that each participant has a positive probability of being recruited into the SAT or external controls, and not having intercurrent events. Assumption 3 states that the conditional expectation of |$Y(0)$| is the same for the trial and the external controls. Assumptions 4 and 5 imply that the intercurrent events occur at random for SAT and the external controls, respectively. These assumptions are satisfied if the covariates X capture all the confounding variables. Take a weight-loss trial for an example, where the intercurrent events are the non-compliance with the prescribed diet. If we assume that the covariates, such as age, baseline weight, and other lifestyle factors, capture all the confounding variables, it follows that given the same covariates for 2 participants, they have the same likelihood of being non-compliant with the diet, regardless of their weight loss. Therefore, we can conclude that the weight loss will not be affected by the occurrence of the non-compliance conditional on these covariates, and the post-intercurrent event outcomes are exchangeable to the observed outcomes.

Theorem 1 provides 3 identification formulas for |$\tau$| under the assumptions in Table 2(A), and Table 2(B) summarizes the necessary models for the identification.

We give some intuitions behind these identification formulas. Theorem 1(a) describes that the individual treatment effect given the covariates X is |$\mu _{1}(X)-\mu _{0}(X)$|⁠. Taking the expectation over the trial population yields the identification for ATE. Theorem 1(b) can be understood as a transportability problem. The first term, corresponding to |$\pi _{S}(X)\mu _{1}(X)$|⁠, adjusts the outcomes from SAT |$SY$| by |$R/\pi _{R_{1}}(X)$|⁠, which weights the observed subjects by their response propensity. The second term, corresponding to |$\pi _{S}(X)\mu _{0}(X)$|⁠, adjusts the outcomes |$(1-S)Y$| from the external controls by |$Rq_{S}(X)/\pi _{R_{0}}(X)$|⁠, which transports the external controls to the trial via the density ratio |$q_{S}(X)$| after response propensity weighting. In Theorem 1(c), it predicts the outcomes for the treatment group by |$\pi _{R_{1}}(X)Y+\lbrace 1-\pi _{R_{1}}(X)\rbrace \mu _{1}(X)$|⁠, which imputes the post-intercurrent event outcomes by |$\mu _{1}(X)$|⁠. Similarly, it predicts the outcomes of the external controls by |$\mu _{0}(X)$|⁠. The difference between these 2 predictions marginalized over the trial population quantifies the ATE.

Based on these identification formulas, infinitely many estimators can be constructed. To develop a more principled estimator, we derive the EIF for |$\tau$|⁠, and the resulting EIF-motivated tilting estimator achieves the rate double robustness, local efficiency, and asymptotic normality. The details of the estimator and these properties are relegated to Theorems S1 and S2 in the Supplementary Materials.

4 SENSITIVITY ANALYSIS UNDER TILTING MODELS

4.1 Assumptions and a graphical representation

Assumptions 3-5 in Table 2 are critical for the identification of |$\tau$|⁠. However, these assumptions may be subject to violations in practice and are unverifiable based on the observed data. Here, we develop the tilting sensitivity models to assess whether the primary analysis result is sensitive to the violation of these assumptions.

Model 1 posits that each unobserved outcome distribution is a “tilted version” of the observed outcomes, where |$\gamma _{S}$|⁠, |$\gamma _{R_{0}}$|⁠, and |$\gamma _{R_{1}}$| are treated as the sensitivity parameters, entailing the level of EC outcome non-exchangeability and the effect of intercurrent events within each arm; see Figure 1 for an illustration when Assumptions 3-5 in Table 2 are violated due to unmeasured confounders |$U_{S}$|⁠, |$U_{R_{0}}$|⁠, and |$U_{R_{1}}$| under the tilting sensitivity models.

FIGURE 1

Schematic plot of the tilting sensitivity models subject to unmeasured confounders |$U_{S}$|⁠, |$U_{R_{0}}$|⁠, and |$U_{R_{1}}$|⁠.

Open in new tabDownload slide

Given negative (or positive) sensitivity parameters, the unobserved outcome distribution is tilted to the left (or right) relative to the distribution of observed outcomes, with smaller (or larger) values receiving greater weight. For example, if |$\gamma _S$| is smaller (or larger) than 0, Model 1 implies that trial participants, if untreated, tend to have smaller (or larger) outcomes compared to the external controls, given the same covariates. Similarly, if |$\gamma _{R_0}$| and |$\gamma _{R_1}$| are both smaller (or larger) than 0, Model 1 implies that participants with intercurrent events tend to have smaller (or larger) outcomes compared to those without such events, given the same covariates. In particular, |$\gamma _{S}=0$| leads to Assumption 3, |$\gamma _{R_{0}}=0$| leads to Assumption 4, and |$\gamma _{R_{1}}=1$| leads to Assumption 5.

Remark 1 implies that our tilting sensitivity model, also known as the exponential tilting model, is connected to the logistic selection specification (2) with flexible formulation, including many non-parametric models, such as sieve approximation, Dirichlet process mixtures, and Bayesian additive regression trees. The logistic selection model, despite its drawbacks noted in Copas and Li (1997), is widely used to assess selection bias in missing data (Robins et al., 2000; Dahabreh et al., 2023) and to conduct sensitivity analyses for unmeasured confounding in causal inference (Franks et al., 2020; Nabi et al., 2024).

4.2 Identification and EIF

The following theorem establishes the non-parametric identification of |$\tau$| when the sensitivity parameters |$\gamma _S$|⁠, |$\gamma _{R_0}$|⁠, and |$\gamma _{R_1}$| are fixed.

The identification formula in Theorem 2 is derived under the same logic as Theorem 1. Under the tilting sensitivity model, we have

where (3) reduces to |$\mu _{1}(X)$| when |$\gamma _{R_{1}}=0$| as the intercurrent events occur at random for SAT, and (4) reduces to |$\mu _{0}(X)$| when |$\gamma _{R_{0}}=\gamma _{S}=0$| as external controls are exchangeable to SAT and the intercurrent events occur at random for external controls. Similarly, we can derive the EIF for |$\tau$| under Model 1 to motivate the semi-parametric efficient estimator.

The EIF in Theorem 3 is constituted by 3 parts. The first part (5) is contributed by the SAT with no intercurrent events (ie, |$S=1,R=1$|⁠); the second part (6) is contributed by the SAT with intercurrent events. When |$\gamma _{R_{1}}=0$|⁠, indicating the intercurrent event occurs at random for SAT, |$g(V;\gamma _{R_{1}})$| equals to |$Y-\mu _{1}(X)$|⁠, and part (6) reduces to (S10) for primary analysis; the third part (7) is contributed by the external controls. When |$\gamma _{R_{0}}=\gamma _{S}=0$|⁠, indicating external controls are exchangeable to SAT and the intercurrent event occurs at random for external controls, |$h(V;\gamma _{R_{0}},\gamma _{S})$| equals to |$R\lbrace Y-\mu _{0}(X)\rbrace /\pi _{R_{0}}(X)$|⁠, and the part (7) reduces to (S11). Next, we construct an estimator for |$\tau$| by solving the empirical mean of |$\phi _{\mathrm{eff}}^{t}(V;P_{0},\gamma _S, \gamma _{R_{0}}, \gamma _{R_{1}})$| with |$P_{0}$| replaced by its estimated counterpart, and present the asymptotic properties for |$\widehat{\tau }^{t}$| in Theorem 4.

Theorem 4 shows that |$\widehat{\tau }^{t}$| is root-n consistent and asymptotically normal for fixed sensitivity parameters |$\gamma _S$|⁠, |$\gamma _{R_{0}}$|⁠, and |$\gamma _{R_{1}}$| when the remainder term |$N^{1/2}\Vert {\rm Rem}^{t}(\widehat{P},P_{0})\Vert _{L_{2}}$| is |$o_{\mathbb {P} }(1)$|⁠. Intuitively, when |$\gamma _{S}=\gamma _{R_{1}}=\gamma _{R_{0}}=0$|⁠, we have |$\Vert \widehat{c}(X;\gamma )-c(X;\gamma )\Vert _{L_{2}}=0$| for any |$\gamma$|⁠, and

Thus, the remainder term reduces to

which is at the same convergence rate as |$\Vert {\rm Rem}(\widehat{P},P_{0})\Vert _{L_{2}}$| for the primary analysis in Theorem S2. The remainder term |$\Vert {\rm Rem}^t(\widehat{P},P_0)\Vert _{L_{2}}$| suggests that the error of |$\widehat{\tau }^{[t]}$| is only affected by the estimation errors of the nuisance models in second-order terms. Therefore, |$\widehat{\tau }^{[t]}$| is more robust and remains consistent when flexible machine learning methods used for nuisance estimation converge at rates faster than |$N^{-1/4}$|⁠. This condition is satisfied by some machine learning methods (Kennedy, 2016; Bradic et al., 2019), which is the so-called rate double robustness (Chernozhukov et al., 2018). Note that the conditional expectations |$b(X;\gamma )$| and |$c(X;\gamma )$| are critical to obtain accurate tilting estimates for the sensitivity analysis. In general, these conditional expectations require heavy computation or strong restrictions on the outcome model. Fortunately, these terms are analytically tractable when the observed outcome model belongs to the class of exponential family mixtures, for example, the Dirichlet processes mixture models (Dorie et al., 2016).

5 CALIBRATING SENSITIVITY PARAMETERS

The magnitude of the sensitivity parameters indicates the strength of the non-ignorability of the indicators |$(R,S)$| given the covariates, which is commonly caused by the existence of unmeasured confounders. However, it is practically infeasible to identify the sensitivity parameters with the observed data. Furthermore, assessing whether a confounder with such strength plausibly exists, given the prior knowledge and domain expertise, is arguably challenging as well. While sensitivity parameters are not directly identifiable from the data, it is reasonable to bound their relative strength using the observed data.

Following the calibration approach proposed by Franks et al. (2020), we introduce a method to determine the plausible quantities for sensitivity analysis based on the observed data. Assume the logistic selection specification outlined in Remark 1 holds. Next, we assume that the relative strength of the unmeasured confounder cannot exceed that of the observed covariates, meaning it should not account for more variation of the indicators as the most important covariate. To measure the relative strength, we adopt the “implicit |$R^{2}$|” concept from Imbens (2003), which generalizes variance-explained measures to the case of binary outcomes. For example, the partial variances |$\rho _{Y(0)\mid X}^{2}$| explained by |$Y(0)$| given X is

where |$\sigma _{Y}^{2}= {\mathbb {E}} [\text{var}\lbrace Y(0)\mid X,S=0\rbrace ]$|⁠. Then, we propose a target value |$(\rho ^{*})^2$| for the unidentified |$\rho _{Y(0)\mid X}^{2}$| using the observed data. In specific, we compute the partial variance explained by each covariate |$X_{j}$| given all other covariates |$X_{-j}$|⁠, and set |$(\rho ^{*})^2=\max _{j}\rho _{X_{j}\mid X_{-j}}^{2}/(1-\max _{j}\rho _{X_{j}\mid X_{-j}}^{2})$|⁠. Here, |$(\rho ^{*})^2$| represents the maximum partial variance explained by adding one covariate |$X_{j}$| to the others, relative to the baseline variance that needs to be explained, referred to as the partial Cohen’s f in Cinelli and Hazlett (2020). Setting |$\rho _{Y(0)\mid X}^{2}=(\rho ^{*})^2$| allows us to calibrate |$\gamma _{S}$| by |$\gamma _{S}^{*}$|⁠, which implies that the information gained by adding |$Y(0)$| to X as a predictor of S is comparable to the maximum information gain by the most important covariate. To calibrate the sensitivity parameter |$\gamma _{S}$|⁠, the following one-to-one mapping is adopted:

Similar bounding procedures apply to the calibration of |$\gamma _{R_{0}}$| and |$\gamma _{R_{1}}$|⁠.

6 SIMULATION STUDY

We first conduct a set of simulations to evaluate the operating characteristics of the proposed estimators under possible model misspecification when Assumptions 3-5 in Table 2 are satisfied. Set the sample sizes of the SAT and EC to be around |$N_{\mathcal {R}}=200$| and |$N_{\mathcal {E}}=500$| with total size |$N=700$|⁠. The covariates |$X\in \mathbb {R}^{5}$| are generated by |$X_{j}\sim N(0.25,1)$| for |$j=1,\cdots ,4$| and |$X_{5}\sim \text{Bernoulli}(0.5)$|⁠. Consider a nonlinear transformation of the covariates and denote |$Z_{j}=\lbrace X_{j}^{2}+2\sin (X_{j})-1.5\rbrace /\sqrt{2}$| for |$j=1,\cdots ,4$| and |$Z_{5}=X_{5}$|⁠. We generate the indicator of being selected to SAT or EC by |$S\mid X\sim \text{Bernoulli}\lbrace \pi _{S}(X)\rbrace ,$| where |$\pi _{S}(X)=\text{logit}^{-1}(\alpha _{S}+0.1\sum _{j=1}^{5}Z_{j})$| and |$\alpha _{S}$| is chosen adaptively to ensure the average of S is about |$N_{\mathcal {R}}/N$|⁠. Next, we generate the indicators of the intercurrent events and the outcomes for SAT and EC by

|$Y\mid X,S=1\sim \mathcal {N}(\sum _{j=1}^{5}Z_{j}/2,1)$|⁠, and |$Y\mid X,S=0\sim \mathcal {N}(\sum _{j=1}^{5}Z_{j}/3,1)$|⁠, where |$(\alpha _{R_{S1}},\alpha _{R_{S0}})$| are adaptively chosen to ensure the average propensity of R is around 0.5. With a large sample size of Monte Carlo simulation, we compute true ATE |$\tau =0.13$|⁠. First, we assess the robustness of the proposed estimator |$\widehat{\tau }^{t}$| when the sensitivity parameters |$\gamma _S = \gamma _{R_0}=\gamma _{R_1}=0$|⁠. Denote the tilting estimator |$\widehat{\tau }^t$| with fixed zero-valued sensitivity parameters as |$\widehat{\tau }^{[0]}$|⁠, we consider 2 model specifications of the propensity (PS) of the participation |$\pi _{S}(X)$| and the intercurrent events |$\pi _{R_{s}}(X)$|⁠, and the outcome models (OM). In particular, we fit the corresponding parametric models with the covariates Z as the correctly specified models or with the covariates X as the misspecified models. We compare our proposed EIF-motivated tilting estimator with other 2 estimators, which are constructed solely based on PS or OM, denoted by |$\widehat{\tau }^{PS}$| and |$\widehat{\tau }^{OM}$|⁠, respectively. Figure 2A shows the point estimation results based on 500 Monte Carlo experiments. When PS and OM are correctly specified, the considered estimators are all unbiased. However, |$\widehat{\tau }^{PS}$| and |$\widehat{\tau }^{OM}$| are biased when their required models are misspecified. Our proposed tilting estimator is shown to be doubly robust with fixed zero-valued sensitivity parameters as it is consistent if either PS or OM is correctly specified.

FIGURE 2

Performance of the PS-based, OR-based, and EIF-motivated tilting estimator with fixed sensitivity parameters: (A) zero-valued under 4 different model specifications, and (B) non-zero-valued when all the nuisance models are correctly specified, based on 500 Monte Carlo simulations.

Open in new tabDownload slide

Table 3(A) shows the absolute bias, standard errors (SE), mean squared error (MSE), coverage rates (CRs), and the average CI lengths of each estimator. We construct the corresponding 95% Wald-type CIs for inference, where the variances are estimated by non-parametric bootstrap with size |$B=50$|⁠. We observe that both |$\widehat{\tau }^{[0]}$| and |$\widehat{\tau }^{OM}$| exhibit the smallest average confidence interval lengths. However, when the OM is misspecified, the CR for |$\widehat{\tau }^{[0]}$| is closer to the nominal level compared to |$\widehat{\tau }^{OM}$|⁠. This observation underscores the double robustness of |$\widehat{\tau }^{[0]}$|⁠, aligning with the results shown in Figure 2 and supporting our claims in Theorem S2 for the primary analysis.

Next, we assess the performance of our proposed tilting estimators |$\widehat{\tau }^{t}$| when Assumptions 3-5 in Table 2 are violated. We keep the same data-generating process for the X, Z, and Y but generate the indicators using Bernoulli sampling with different propensities:

where |$(\gamma _{S},\gamma _{R_{1}},\gamma _{R_{0}})$| control the strength of the unmeasured confounder, describing how the propensity depends on the potential outcome after accounting for the covariates. We consider a range of sensitivity parameters where the indicators are equally confounded by the potential outcomes, that is, |$\gamma =\gamma _{S}=\gamma _{R_{1}}=\gamma _{R_{0}}$|⁠. Under the hypothetically fixed sensitivity parameters |$(\gamma _{S},\gamma _{R_{1}},\gamma _{R_{0}})$| and assuming that all the nuisance models are correctly specified, Figure 2(B) shows the point estimation for all 3 estimators, and Table 3(B) presents the finite-sample performances of |$\widehat{\tau }^{t}$| in details. One important finding is that our proposed estimator is always consistent across a range of fixed sensitivity parameters, as indicated by its small bias and approximately correct CRs. In contrast, the other 2 estimators are unable to handle the joint sensitivity analysis for multiple assumptions, suggested by their non-negligible biases for fixed non-zero sensitivity parameters.

7 REAL-DATA APPLICATION

We examine a study designed to estimate an antidepressant drug effect on the scores of the Hamilton Depression Rating Scale for 17 items (HAMD-17). This study was conducted under the Auspices of the Drug Information Association, which collects the data at baseline and weeks 1, 2, 4, 6, and 8 for |$N=196$| patients with 99 in the control group and 97 in the treatment group. However, some patients may drop out during the study for various reasons. Our primary interest is the ATE on the change of HAMD-17, irrespective of the intercurrent events such as dropout-related missing data. According to the guidelines in ICH (2021), the ATE is defined as the mean difference of the change in the HAMD-17 scores from the baseline to the final time point in week 8. We adhere to the same analysis plan as Liu et al. (2024) with covariates X, including the investigation sites and baseline HAMD-17 scores. Let |$Y(a)$| and R be the change of HAMD-17 scores under treatment a and the indicator of whether a patient stayed in the study at week 8.

In our application, we consider the original trial as a single-arm trial, with the concurrent control group being considered as external controls to illustrate the proposed sensitivity analysis. We assume the potential outcomes |$Y(a)$| follow a Gaussian mixture model fitted using the R package flexmix. Next, we bound the magnitude of the sensitivity parameters using the approach described in Section 5. We illustrate this approach with the baseline HAMD-17 score, which is the most important predictors in terms of partial variance explained, with |$(\rho _{S}^{*})^2\approx 0.02$|⁠, |$(\rho _{R_{1}}^{*})^2\approx 0.11$|⁠, and |$(\rho _{R_{0}}^{*})^2\approx 0.04$|⁠. To map these values to the sensitivity parameters, we apply the one-to-one mapping formula (8), and obtain the calibrated sensitivity parameters |$|\gamma _{R_{0}}^{*}|\approx 0.02$|⁠, |$|\gamma _{R_{1}}^{*}|\approx 0.02$|⁠, and |$|\gamma _{S}^{*}|\approx 0.01$|⁠. Figure 3 illustrates the ATE estimates across a range of hypothetical sensitivity parameters, adjusting for the potential outcomes as the unmeasured confounder in the logistic selection specification. The shaded area indicates the unmeasured confounder with impacts up to the values of the calibrated sensitivity parameters. Here, “NS” denotes “not significant,” meaning the 95% confidence interval of the ATE contains 0.

FIGURE 3

Average treatment effects of an antidepressant drug effect on the HAMD-17 scores over a grid of hypothesized sensitivity parameters under the tilting sensitivity models.

Open in new tabDownload slide

When the assumptions are satisfied, that is, all the sensitivity parameters equal 0, the ATE estimates are |$\widehat{\tau }^{t} = -1.42$| with the 95% bootstrap CI as |$(-2.80, -0.05)$|⁠, which is statistically significant. Next, we assume the impact of the confounders acts toward hurting our preferred hypothesis, that is, |$\gamma _{S}^{*}< 0$|⁠. Here, a negative value of |$\gamma _{S}$| suggests that the unobserved change in HAMD-17 scores in the concurrent control group tends to be lower (ie, better) than the observed external controls, reducing the absolute value of the effect size. This could occur if patients are more likely to participate in the single-arm trial when less depressed. When |$\gamma _{S} = -0.01$| and |$\gamma _{R_{0}} = \gamma _{R_{1}} = -0.02$|⁠, the estimated treatment effect of the antidepressant drug becomes |$\widehat{\tau }^{t} = -1.28$|⁠, where the unmeasured confounder is as strong as the baseline HAMD-17 scores. Although the tilted estimate is below 0, suggesting the effectiveness of the drug on the HAMD-17 scores, it is no longer statistically significant, as its 95% CI is |$(-2.65, 0.09)$|⁠. Thus, following our sensitivity analyses, we show that the magnitude of the possible drug effects on the HAMD-17 scores is robust, but the significance of such effects is not robust against unmeasured confounding at the strength of baseline HAMD-17 scores. However, domain knowledge is still required to consider the plausibility of unmeasured confounders of such strength level under this situation.

8 DISCUSSION

In this paper, we develop a semi-parametric efficient framework for sensitivity analysis under the tilting models. Motivated by Tukey’s factorization, this framework effectively separates the model checking from the sensitivity analysis, which does not rely on any modeling assumption and fits perfectly well with the EIF-motivated tilting estimators. Decoupling the sensitivity analysis from the model fit assessment is crucial and ubiquitous within the model-based sensitivity analysis. However, joint sensitivity analysis for multiple assumptions remains largely unexplored to the best of our knowledge. By simultaneously assessing the EC outcome mean non-exchangeability and the effects of intercurrent events, our framework hopes to shed more light on the advancements of joint modeling for sensitivity analysis.

Future work could extend our framework to longitudinal trials with intercurrent events, particularly those with irregular and informative observation patterns, as discussed by Yang (2021) and Smith et al. (2024). Such an extension may increase the number of sensitivity parameters across multiple time points, and introduce additional challenges in deriving the EIFs conditioned on the historical information. Another potential extension could focus on the choices of sensitivity parameters, which is profoundly useful in practice. Our approach relies on bounding the magnitude of sensitivity parameters using observed data, and substantive domain expertise should be consulted to examine whether an unmeasured confounder with such strength is plausible. For some hybrid control designs, the sensitivity parameters can be partially identified with the help of concurrent controls; similar ideas have been explored in Gao et al. (2024) to adjust for EC outcome mean non-exchangeability. Thus, the internal validity from the hybrid controls can be leveraged to inform the choices of sensitivity parameters. In summary, our proposed semi-parametric sensitivity analysis is both efficient and flexible as it is rate-doubly robust, locally optimal, and can be incorporated with a range of models with a modern causal inference workflow.

ACKNOWLEDGMENTS

We thank the editor, associate editor, and two anonymous referees for their constructive suggestions and valuable feedback on our manuscript.

FUNDING

None declared.

CONFLICT OF INTEREST

None declared.

DATA AVAILABILITY

The data that support the findings in this paper are available in the Drug Information Association Missing Data at https://www.lshtm.ac.uk/research/centres-projects-groups/missing-data#dia-missing-data, collected by Mallinckrodt et al. (2014), and are also provided at https://github.com/Gaochenyin/SensDR.

REFERENCES