Flexible sensitivity analysis for causal inference in observational studies subject to unmeasured confounding

Consider the causal impact of smoking on individuals’ health as a motivating example. Randomly assigning individuals to smoke is neither practical nor ethical, thus we need to rely on observational data to learn the causal effect of interest. Specifically, the observational study in Bazzano et al. (2003) explored whether cigarette smoking positively affects homocysteine levels, a known risk factor for cardiovascular disease. Drawing upon observations from the U.S. National Health and Nutrition Examination Survey (NHANES) 2005–2006, the study compares homocysteine levels between smokers and non-smokers while adjusting for observed covariates like age, gender, race, and body mass index (BMI) under the unconfoundedness assumption. However, in the presence of suspected underlying factors associated with both increased probability of smoking and elevated homocysteine levels, such as certain gene expressions, estimates adjusted solely for observed covariates are prone to bias. In such instances, conducting sensitivity analyses becomes essential to discern how the estimated causal effect changes in response to the presence of unmeasured confounders. Such analyses help determine the degree to which we need to deviate from the unconfoundedness assumption to attribute the estimated positive effect to unobserved confounders. We will revisit this example in Section 5.1 to illustrate our proposed methods.

The remainder of the paper is organized as follows. In Section 2, we review causal inference in observational studies under the unconfoundedness assumption. In Section 3, we introduce our sensitivity analysis framework and provide identification and estimation procedures. In Section 4, we provide practical suggestions for calibration and specification of the sensitivity parameters. In Section 5, we revisit the motivating example, apply the proposed method to this real-world application, and conduct a simulation study to evaluate the finite sample performance of our proposed methods. In Section 6, we discuss alternative sensitivity parameters as well as the extensions to other causal inference settings. In Section 7, we conclude. In the supplementary material, we present additional technical details.

Assumption 1 is strong and untestable. Observational studies are often plagued with unmeasured confounding, that is, the existence of some hidden variables $U$ that affect both the treatment and outcome simultaneously. Under these circumstances, we need to conduct sensitivity analyses to assess the impact of $U$ on the causal estimates. Existing methods only serve purposes in specific statistical estimation or testing strategies, as reviewed in Section 1. To provide a unified framework that can deal with the standard estimators $\hat{\tau}^{\textup{reg}}$, $\hat{\tau}^{\textup{ht}}$ and $\hat{\tau}^{\textup{dr}}$ reviewed in Section 2 simultaneously, we adopt the following parametrization of confounding.

In Definition 1, $\varepsilon_{1}(X)$ and $\varepsilon_{0}(X)$ are two sensitivity parameters. In sensitivity analysis, we can first fix them to obtain the corresponding estimators and then vary them within a range to obtain a sequence of estimators. Definition 1, as well as the theory below, allows $\varepsilon_{1}(X)$ and $\varepsilon_{0}(X)$ to depend on the covariates. For simplicity of implementation, we can further assume that they are not functions of the covariates.

The parametrization in Definition 1 compares the observable conditional means of the potential outcomes, $E\{Y(1)\mid Z=1,X\}$ and $E\{Y(0)\mid Z=0,X\}$, with the corresponding counterfactual conditional means of the potential outcomes, $E\{Y(1)\mid Z=0,X\}$ and $E\{Y(0)\mid Z=1,X\}$, respectively. It quantifies the violation of unconfoundedness in Assumption 1. Technically, the form of unconfoundedness in Assumption 1 is stronger than we need if the parameter of interest is $\tau$. When $\varepsilon_{1}(X)=\varepsilon_{0}(X)=1$ in Definition 1, we can recover all identification and estimation results for $\tau$ in Section 2. We revisit the running example in Section 1.1. If suspected factors associated with both increased probability of smoking and elevated homocysteine levels exist, under Definition 1, $\varepsilon_{1}(X)>1$ and $\varepsilon_{0}(X)>1$, leading to an upward bias of estimators of $\tau$ constructed under the unconfoundedness assumption. The sensitivity parameters in Definition 1 are also related to the parametrization based on the treatment selection model $\textup{pr}\{Z=z\mid X,Y(z)\}$ for $z=0,1$. We provide more discussion in Section S1.8 in the supplementary material.

Robins (1999) and Scharfstein et al. (1999) discussed the potential use of the parametrization in Definition 1. Similar parametrizations also appeared in other settings for causal inference (Tchetgen and Shpitser 2012; VanderWeele et al. 2014; Ding and Lu 2017; Yang and Lok 2018; Jiang et al. 2022). The confounding function in Kasza et al. (2017) also adopts ratios between the marginal means, $E\{Y(1)\mid Z=1\}/E\{Y(1)\mid Z=0\}$, in the special case of binary outcomes, without specifying it as a flexible function of the observed covariates. Franks et al. (2020) used a similar parameterization in Bayesian inference and pointed out that the parametrization in Definition 1 does not impose any testable restrictions on the observed data distributions. Nevertheless, despite the attractiveness of this parametrization, the statistical theory for sensitivity analysis is still missing for the standard estimators in the canonical setting of observational studies.

This subsection will present the key identification and estimation results for $\tau$, in parallel with the results under the unconfoundedness assumption. In particular, by setting $\varepsilon_{1}(X)=\varepsilon_{0}(X)=1$ in the formulas below, we can recover the corresponding results under the unconfoundedness assumption. We will provide theory and methods for the outcome regression, the inverse propensity score weighting, and the doubly robust estimators for $\tau$ under Definition 1.

In Theorem 1, we essentially get the identification for $E\{Y(1)\}$ and $E\{Y(0)\}$ from the identification formulas for the two counterfactual conditional means in (4) and (5), which also lead to the two identification formulas for $\tau$ in (6) and (7). Under the unconfoundedness assumption, $\varepsilon_{1}(X)=\varepsilon_{0}(X)=1$, (7) reduces to the first identification formula in (2). For general $\varepsilon_{z}(X)$, we need to reweight the conditional means $\mu_{z}(X)$ in the identification of the counterfactual means $E\{Y(z)\mid Z=1-z\}$.

With the fitted outcome model, we can construct two plug-in estimators corresponding to (6) and (7) called the predictive and projective estimators for $\tau$, respectively:

$\displaystyle\hat{\tau}^{\textup{pred}}$

$\displaystyle=$

$\displaystyle n^{-1}\sum_{i=1}^{n}\left\{Z_{i}Y_{i}+(1-Z_{i})\hat{\mu}_{1}(X_{i})/\varepsilon_{1}(X_{i})\right\}-n^{-1}\sum_{i=1}^{n}\left\{Z_{i}\hat{\mu}_{0}(X_{i})\varepsilon_{0}(X_{i})+(1-Z_{i})Y_{i}\right\},$

and

$\displaystyle\hat{\tau}^{\textup{proj}}$

$\displaystyle=$

$\displaystyle n^{-1}\sum_{i=1}^{n}\left\{Z_{i}\hat{\mu}_{1}(X_{i})+\frac{(1-Z_{i})\hat{\mu}_{1}(X_{i})}{\varepsilon_{1}(X_{i})}\right\}-n^{-1}\sum_{i=1}^{n}\left\{Z_{i}\hat{\mu}_{0}(X_{i})\varepsilon_{0}(X_{i})+(1-Z_{i})\hat{\mu}_{0}(X_{i})\right\}.$

The terminologies “predictive” and “projective” are from the survey sampling literature (Firth and Bennett 1998; Ding and Li 2018). The estimators $\hat{\tau}^{\textup{pred}}$ and $\hat{\tau}^{\textup{proj}}$ differ slightly: the former uses the observed outcomes when available; in contrast, the latter replaces the observed outcomes with the fitted values.

More interestingly, we can also identify $\tau$ by an inverse propensity score weighting formula although Definition 1 only involves the conditional means of the potential outcomes.

Theorem 2 modifies the classic inverse probability weighting formulas with two extra factors $w_{1}(X)$ and $w_{0}(X)$ depending on both the propensity score and the sensitivity parameters. With the fitted propensity score model, (8) motivates the following estimator for $\tau$:

$\displaystyle\hat{\tau}^{\textup{ht}}=n^{-1}\sum_{i=1}^{n}\hat{w}_{1}(X_{i})\frac{Z_{i}Y_{i}}{\hat{e}(X_{i})}-n^{-1}\sum_{i=1}^{n}\hat{w}_{0}(X_{i})\frac{(1-Z_{i})Y_{i}}{1-\hat{e}(X_{i})},$

where $\hat{w}_{1}(X_{i})=\hat{e}(X_{i})+\{1-\hat{e}(X_{i})\}/\varepsilon_{1}(X_{i})$ and $\hat{w}_{0}(X_{i})=\hat{e}(X_{i})\varepsilon_{0}(X_{i})+1-\hat{e}(X_{i})$ are the estimated weights.

The identification formulas based on outcome regression and inverse propensity score weighting motivate us to develop a combined strategy that leads to doubly robust estimation. Following Bang and Robins (2005), we calculate the efficient influence function for $\tau$. The details of the efficient influence function are in Bickel et al. (1993), and the proof of Theorem 3 is in Section S3.3 in the supplementary material. Readers focusing only on applying the method do not need to understand the efficient influence function to understand our proposed estimation procedure.

The efficient influence function for $\tau$ has mean $0$ by definition, so $\tau$ has the following representation:

$$\tau=E\left[w_{1}(X)\frac{Z}{e(X)}Y-\frac{\{Z-e(X)\}\mu_{1}(X)}{e(X)\varepsilon_{1}(X)}\right]-E\left[w_{0}(X)\frac{1-Z}{1-e(X)}Y-\frac{\{e(X)-Z\}\mu_{0}(X)\varepsilon_{0}(X)}{1-e(X)}\right]$$

which motivates the following estimator based on the fitted propensity score and outcome models:

$$\hat{\tau}^{\textup{dr}}=n^{-1}\sum_{i=1}^{n}\left\{\hat{w}_{1}(X_{i})\frac{Z_{i}Y_{i}}{\hat{e}(X_{i})}-\frac{\check{Z}_{i}\hat{\mu}_{1}(X_{i})}{\hat{e}(X_{i})\varepsilon_{1}(X_{i})}-\hat{w}_{0}(X_{i})\frac{(1-Z_{i})Y_{i}}{1-\hat{e}(X_{i})}-\frac{\check{Z}_{i}\hat{\mu}_{0}(X_{i})\varepsilon_{0}(X_{i})}{1-\hat{e}(X_{i})}\right\}.$$

It is numerically identical to the following forms

	$\displaystyle\hat{\tau}^{\textup{dr}}$	$\displaystyle=$	$\displaystyle\hat{\tau}^{\textup{ht}}-n^{-1}\sum_{i=1}^{n}\check{Z}_{i}\left\{\frac{\hat{\mu}_{1}(X_{i})}{\hat{e}(X_{i})\varepsilon_{1}(X_{i})}+\frac{\hat{\mu}_{0}(X_{i})\varepsilon_{0}(X_{i})}{1-\hat{e}(X_{i})}\right\}$
		$\displaystyle=$	$\displaystyle\hat{\tau}^{\textup{pred}}+n^{-1}\sum_{i=1}^{n}\left\{\frac{1-\hat{e}(X_{i})}{\varepsilon_{1}(X_{i})}\frac{Z_{i}\check{Y}_{i}}{\hat{e}(X_{i})}-\hat{e}(X_{i})\varepsilon_{0}(X_{i})\frac{(1-Z_{i})\check{Y}_{i}}{1-\hat{e}(X_{i})}\right\}$
		$\displaystyle=$	$\displaystyle\hat{\tau}^{\textup{proj}}+n^{-1}\sum_{i=1}^{n}\left\{\hat{w}_{1}(X_{i})\frac{Z_{i}\check{Y}_{i}}{\hat{e}(X_{i})}-\hat{w}_{0}(X_{i})\frac{(1-Z_{i})\check{Y}_{i}}{1-\hat{e}(X_{i})}\right\},$

which, similar to (3), augments $\hat{\tau}^{\textup{ht}}$ by imputed outcomes and modifies $\hat{\tau}^{\textup{pred}}$ and $\hat{\tau}^{\textup{proj}}$ by inverse propensity score weighted residuals. Importantly, Theorem 4 below shows that $\hat{\tau}^{\textup{dr}}$ is doubly robust.

Under standard regularity conditions, the doubly robust estimator $\hat{\tau}^{\textup{dr}}$ is asymptotically normally distributed with mean $\tau$ and variance achieving the variance of the efficient influence function. Therefore, we can conduct statistical inference based on the closed-form plug-in variance estimator or the nonparametric bootstrap variance estimator of $\hat{\tau}^{\textup{dr}}$. We can also construct double machine learning estimators and apply the cross-fitting technique to achieve the same asymptotic distribution under a relaxed condition on the complexity of the spaces of the nuisance parameters (Chernozhukov et al. 2018).

Specifying the ranges of the sensitivity parameters is fundamentally challenging in causal inference because the observed data do not provide direct information on the sensitivity parameters. We provide a method to calibrate the sensitivity parameters using the observed covariates. The identification and estimation procedures rely on the unidentifiable sensitivity parameters $(\varepsilon_{1}(X),\varepsilon_{0}(X))$. Calibration based on the observed covariates provides us guidance on their sizes, helping us to understand their relative magnitudes in the analytical framework. Zhang and Small (2020) proposed a calibration method based on the observed covariates as a benchmark to understand how sensitive the estimates are with respect to an unobserved confounder in matched observational studies. In our sensitivity analysis framework, we do not directly model the unobserved confounders or their association with the treatment or potential outcome. Also, our proposed sensitivity parameter depends on the set of observed covariates $X$ we have controlled for, thus it is challenging to directly characterize the corresponding calibration of the observed covariates. Below we provide a leave-one-covariate-out approach to calibrate the sensitivity parameters. Similar approaches appeared in other sensitivity analysis methods (Chapter 17 in Ding 2024).

For each observed covariate $X_{j}$ in the $p$-dimensional $X=(X_{1},\ldots,X_{p})$, $j=1,\ldots,p$, we first drop it as if it is an unobserved confounder, then estimate the value of sensitivity parameters with $X_{j}$ unobserved. Under the unconfoundedness Assumption 1, we have

$\displaystyle\varepsilon_{z}(X_{-j})$

$\displaystyle=$

$\displaystyle\frac{E\{Y(z)\mid Z=1,X_{-j}\}}{E\{Y(z)\mid Z=0,X_{-j}\}}\ =\ \frac{E\{\mu_{z}(X)\mid Z=1,X_{-j}\}}{E\{\mu_{z}(X)\mid Z=0,X_{-j}\}}$

for $z=0,1$, where $X_{-j}$ denotes the $p-1$ dimensional observed covariates after dropping covariate $X_{j}$. The $\varepsilon_{z}(X_{-j})$ measures how much the observed data deviates from the unconfoundedness assumption when we delete an observed confounder $X_{j}$, and can be estimated using observed data. Thus, we can summarize the distribution of $\varepsilon_{z}(X_{-j})$ to characterize how strong a contribution each covariate $X_{j}$ has as a confounder. To summarize $\varepsilon_{z}(X_{-j})$ from a function of $X_{-j}$ to a scalar, we can further marginalize over the distribution of $X_{-j}$ to compute the mean of sensitivity parameter ignoring $X_{j}$, or use other summary statistics such as maximum, minimum, upper or lower quantiles of the estimated distribution to guide us on the interpretation of the magnitude of the sensitivity parameter $\varepsilon_{z}(X)$. We will illustrate this calibration method in Section 5.1. In addition, we can also calibrate using a subset of covariates and provide a leave-multiple-covariates-out approach. The definition, identification, and estimation of the sensitivity parameters $\varepsilon_{z}(X_{-\mathcal{S}})$ similarly follow from the unconfoundedness assumption, where $\mathcal{S}\subset\{1,\ldots,p\}$ denotes the indices of the subset of covariates we use for calibration.

In this subsection, we provide practical suggestions on specifying the sensitivity parameters and reporting the results. When the sensitivity parameters are constants not depending on $X$, i.e., $\varepsilon_{1}(X)=\varepsilon_{1}$ and $\varepsilon_{0}(X)=\varepsilon_{0}$, the proposed estimators of $\tau$ monotonically decrease in both $\varepsilon_{1}$ and $\varepsilon_{0}$. Therefore, it often suffices to present one side of the parameters for sensitivity analysis if we are interested in how strong the unobserved confounding should be to explain away the estimated causal effect under the unconfoundedness assumption. For instance, if the estimated causal effect under unconfoundedness is positive, researchers can report a table of estimated average causal effects for different $(\varepsilon_{1},\varepsilon_{0})$ with $\varepsilon_{1}\geq 1$ and $\varepsilon_{0}\geq 1$. Such table shows how robust the positive result is when we increase the sensitivity parameters, and presents how large the sensitivity parameters should be to overturn the estimated significant positive effect, while the other direction with small sensitivity parameters is usually of less interest since it will only increase the estimated positive effect. For a similar reason, when the estimated average causal effect is negative under the unconfoundedness assumption, we suggest reporting results with $\varepsilon_{1}\leq 1$ and $\varepsilon_{0}\leq 1$. However, in the case when we would like to learn how strong the unobserved confounding is needed to match some previously estimated causal effects with larger magnitude, the other direction may be of interest as well.

In our proposed framework, we allow for the dependence of the sensitivity parameters on observed covariates. When the sensitivity parameters depend on $X$, two alternative approaches can be pursued: modeling the sensitivity parameters explicitly or adopting a worst-case interpretation. In the first approach, we can assume $\varepsilon_{z}(X)=\exp(\alpha_{z}+\beta_{z}^{{\mathrm{\scriptscriptstyle T}}}X)$ for $z=0,1$, and specify the value of $\varepsilon_{z}(X)$ based on domain knowledge about the parameters $(\alpha_{z},\beta_{z})$. In the second approach, under the special scenarios when the outcome is non-negative (e.g., binary), treating the estimated results with constant values $(\varepsilon_{1},\varepsilon_{0})$ serves as a worst-case bound for the true average causal effect. Here, $\varepsilon_{z}$ is determined as either $\min(\varepsilon_{z}(X))$ or $\max(\varepsilon_{z}(X))$, depending on the direction of the average causal effect. We give a formal result below.

Proposition 1 gives bounds based on the three identification formulas in Theorems 1–3. Define $\hat{\tau}^{*}_{\textsc{l}}$ and $\hat{\tau}^{*}_{\textsc{u}}$ as estimators using the same formula as $\hat{\tau}^{*}$ by replacing $(\varepsilon_{1}(X_{i}),\varepsilon_{0}(X_{i}))$ by $(\varepsilon_{1,\textsc{u}},\varepsilon_{0,\textsc{u}})$ and $(\varepsilon_{1,\textsc{l}},\varepsilon_{0,\textsc{l}})$, respectively, for $*=\textup{pred, proj, ht, dr}$. Following from Proposition 1, we can treat $\hat{\tau}^{*}_{\textsc{l}}$ and $\hat{\tau}^{*}_{\textsc{u}}$ as estimators of lower and upper bounds of $\tau$, respectively. When the potential outcomes are non-positive, we have analogous results to Proposition 1.

Suppose our estimated average causal effect under the unconfoundedness assumption is positive, we focus on $\varepsilon_{1}(X)\geq 1$, $\varepsilon_{0}(X)\geq 1$ to answer the question of how large the sensitivity parameters have to be to explain away the estimated positive causal effect. Looping over the pre-specified lists of constant $\varepsilon_{1}$ and $\varepsilon_{0}$, we are finding how large the lower bounds of $\varepsilon_{1}(X)$ and $\varepsilon_{0}(X)$ are to explain away the positive effect.

We apply the proposed sensitivity analysis method to a real-world application. We re-analyze the observational study in Bazzano et al. (2003) to study whether cigarette smoking has a causal effect on homocysteine levels, the elevation of which is considered a risk factor for cardiovascular disease. We compare the homocysteine levels in daily smokers and non-smokers based on the NHANES 2005–2006 data. The NHANES is a series of surveys whose target population is the civilian, noninstitutionalized U.S. population. Since 1999, the NHANES has been conducted annually by interviewing individuals in their homes and asking the responders to complete the health examination competent of the survey. We use the cross-sectional data collected in the NHANES 2005–2006, which includes observed covariates such as gender, age, education level, BMI, and a measure of poverty. Consider the daily smokers to be in the treated group and the never smokers to be in the control group. The baseline homocysteine level in the control group is 7.86 umol/L. Under the unconfoundedness assumption with $\varepsilon_{1}(X)=\varepsilon_{0}(X)=1$, the estimated $\tau$ using the doubly robust estimator is 1.48 with a 95% confidence interval (0.78, 2.18). We use the logit model for the propensity score and the linear model for the outcome, respectively, and use the nonparametric bootstrap method to estimate the variance. Despite a large collection of observed covariates collected in the survey, unmeasured factors such as the lifestyle of their spouses or parents, work intensity, stress level, and certain gene expressions possibly contribute to both a higher probability of smoking and higher homocysteine level. With the concern about the presence of unmeasured confounders, we conduct sensitivity analysis using the doubly robust estimator $\hat{\tau}^{\textup{dr}}$. In this example, it is challenging to specify a correct model of sensitivity parameters as a function of $X$, and the outcome of interest is non-negative. Thus, we adopt the conservative worst-case interpretation in Proposition 1, assume the sensitivity parameters $\varepsilon_{1}(X)$ and $\varepsilon_{0}(X)$ do not depend on $X$, and vary both of them from 0.75 to 1.25. Table 1 reports the results of the doubly robust estimator and its 95% confidence interval for different $(\varepsilon_{1},\varepsilon_{0})$ where we again use nonparametric bootstrap to estimate the variance of $\hat{\tau}^{\textup{dr}}$. We omit the detailed results for small $\varepsilon_{1}$ and $\varepsilon_{0}$ since they all give significantly positive results at the 0.05 level. As $\varepsilon_{1}$ and $\varepsilon_{0}$ increase, the significance level decreases, and the direction of the point estimator changes when both are very large.

However, the length of the confidence interval is not a monotone function of $\varepsilon_{1}(X)$ and $\varepsilon_{0}(X)$ since the asymptotic variance of $\hat{\tau}^{\textup{dr}}$ under our sensitivity analysis framework is a complicated function of the sensitivity parameters. In Section S1.7 in the supplementary material, we give the specific form of the semiparametric efficiency bound for $\tau$ based on its efficient influence function in Theorem 3.

Figure 1 shows the contour plot of the estimated average causal effect using the doubly robust estimator $\hat{\tau}^{\textup{dr}}$ with various values of $\varepsilon_{1}$ and $\varepsilon_{0}$, where the contour lines are generated through grid search. We also implement the leave-one-covariate-out calibration method proposed in Section 4.1 in the example. For each observed covariate $X_{j}$, we label the maximum estimated value of $\varepsilon_{1}(X_{-j})$ and $\varepsilon_{0}(X_{-j})$ in the contour plot. This suggests that to explain away the estimated positive causal effect, the unobserved confounder has to be stronger than all observed covariates. Figure 1 also shows that under Assumption 1, covariates such as BMI and poverty are weak confounders with $\varepsilon_{1}(X)$ and $\varepsilon_{0}(X)$ in the range of $\varepsilon_{1}(X_{-\textup{BMI}})\in[1.002,1.020]$, $\varepsilon_{0}(X_{-\textup{BMI}})\in[0.998,1.001]$ and $\varepsilon_{1}(X_{-\textup{poverty}})\in[0.972,1.026]$, $\varepsilon_{0}(X_{-\textup{poverty}})\in[0.998,1.002]$, respectively.

In this subsection, we conduct a simulation study to evaluate the finite sample performance of our proposed sensitivity analysis framework with the doubly robust estimator.

We generate the data using the following data-generating process with sample size $n=500$. We first generate observed covariates $X=(X_{1},X_{2},X_{3})^{\mathrm{\scriptscriptstyle T}}$ with independent $X_{j}\sim\mathcal{N}(0,0.5^{2})$ for $j=1,2$ and $X_{3}\sim\textup{Bernoulli}(0.5)$ and a binary unobserved confounder $U\sim\textup{Bernoulli}(0.5)$. Next generate the treatment $Z\mid(X,U)\sim\textup{Bernoulli}(0.25+0.5U)$. Then generate the potential outcomes with log-normal conditional distribution, $\log\{Y(1)\}=X_{1}+X_{2}+bU+e_{1}$ and $\log\{Y(0)\}=X_{1}+X_{2}+bX_{3}+e_{0}$, where $e_{1}$ and $e_{0}$ are independent $\mathcal{N}(0,0.5^{2})$. Under this data-generating process, the population average causal effect $\tau$ equals 0. Since $U$ is an unobserved confounder of $Y(1)$ and $Z$ conditional on $X$, the unconfoundedness assumption is violated, thus estimators constructed using merely observed variables will be asymptotically biased, leading to invalid inference. Under the log-normal model of the potential outcomes and the Bernoulli treatment assignment, the true $\varepsilon_{1}(X)=\{0.75\exp(b)+0.25\}/\{0.25\exp(b)+0.75\}$, and $\varepsilon_{0}(X)=1$ because the distribution of $\log\{Y(0)\}$ does not depend on $U$. We evaluate the performance of $\hat{\tau}^{\textup{dr}}$ with various values of $b\in(0,0.2,0.3,0.5,1,1.5)$. For $b>0$, the true $\varepsilon_{1}(X)>1$ and assuming the unconfoundedness assumption will over-estimate $\tau$. Thus, we also take a pre-specified list of values larger than 1, $\varepsilon_{1}(X)\in\{1,1.10,1.16,1.28,1.60,1.93\}$ and $\varepsilon_{0}(X)=1$. We use a nonparametric bootstrap with 200 bootstrap samples to estimate the asymptotic variance. For each value of $b$ and each pair $(\varepsilon_{1}(X),\varepsilon_{0}(X))$, we compute the coverage probability of the 95% confidence interval.

When we observe a positive $\hat{\tau}^{\textup{dr}}$ under the unconfoundedness assumption, at each $\varepsilon_{1}$ level, we reject the null hypothesis that $\tau=0$ if the constructed lower bound of the 95% confidence interval is positive. With a positive $b$, $\tau$ will be overestimated if $\varepsilon_{1}$ takes a value smaller than the truth, thus it is likely to make a false rejection. While for $\varepsilon_{1}$ larger than the truth, we will make a conservative inference thus the false rejection is less likely. To evaluate this, we also compute the frequency of making a false rejection of a significantly positive average causal effect. The number of Monte Carlo samples is 500.

Table 2 reports the simulation results including the coverage and false rejection rates of our proposed estimation procedure under various data-generating processes. First, when the sensitivity parameters are set to be equal to the true values, the coverage rate is correct. Second, the false rejection rates are valid if the value of $\varepsilon_{1}$ is specified to be larger than the truth. The proposed sensitivity analysis procedure has valid finite sample performance.

We expressed the sensitivity parameters as ratios between the observed outcome means and the corresponding counterfactual outcome means. Alternatively, we can also express the sensitivity parameters on the difference scale as Definition 2 below (Robins 1999).

In Definition 2, $\delta_{1}(X)$ and $\delta_{0}(X)$ are the two sensitivity parameters.

Under Definition 2, we can identify the two counterfactual means by $E\{Y(1)\mid Z=0\}=E\left\{\mu_{1}(X)-\delta_{1}(X)\mid Z=0\right\}$ and $E\{Y(0)\mid Z=1\}=E\left\{\mu_{0}(X)+\delta_{0}(X)\mid Z=1\right\}$. Therefore, we have the following identification formulas for $E\{Y(1)\}$ and $E\{Y(0)\}$ based on outcome regression, inverse probability weighting, and doubly robust estimation:

	$\displaystyle E\{Y(1)\}$	$\displaystyle=$	$\displaystyle E\{\mu_{1}(X)\}-E\{(1-Z)\delta_{1}(X)\}\ =\ E\left\{\frac{ZY}{e(X)}-(1-Z)\delta_{1}(X)\right\}$
		$\displaystyle=$	$\displaystyle E\left[\mu_{1}(X)+\frac{Z\{Y-\mu_{1}(X)\}}{e(X)}-(1-Z)\delta_{1}(X)\right],$

and

	$\displaystyle E\{Y(0)\}$	$\displaystyle=$	$\displaystyle E\{\mu_{0}(X)\}+E\{Z\delta_{0}(X)\}\ =\ E\left\{\frac{(1-Z)Y}{1-e(X)}+Z\delta_{0}(X)\right\}$
		$\displaystyle=$	$\displaystyle E\left[\mu_{0}(X)+\frac{(1-Z)\{Y-\mu_{1}(X)\}}{1-e(X)}+Z\delta_{0}(X)\right],$

where the two doubly robust formulas come from the efficient influence functions for $E\{Y(1)\}$ and $E\{Y(0)\}$, respectively. Taking the difference between these two expectations, we have the identification formulas for the average causal effect $\tau$, which can be written as the standard identification formulas under unconfoundedness assumption minus $\delta=E\{Z\delta_{0}(X)\}+E\{(1-Z)\delta_{1}(X)\}.$ When the sensitivity parameters $\delta_{0}(X)$ and $\delta_{1}(X)$ do not depend on $X$, the correction simplifies to $\delta=e\delta_{0}+(1-e)\delta_{1}$, a constant that depends on the probability of the treatment and $(\delta_{1},\delta_{0})$.

Similarly, for the average causal effect on the treated, the correction is $\delta_{\textsc{t}}=E\{Z\delta_{0}(X)\}/e.$ When the sensitivity parameter $\delta_{0}(X)$ does not depend on $X$, the correction simplifies to $\delta_{0}$, a constant that equals the sensitivity parameter itself.

From the above discussion, sensitivity analysis based on the difference scale is much simpler than that based on the ratio scale. Its simplicity is an advantage. However, it also has disadvantages. First, the sensitivity analysis estimators reduce to the corresponding classic estimators minus prespecified constants. It makes sensitivity analysis a tautology of specifying how different the estimators and the estimands are. Second, the sensitivity parameters in Definition 2 depend on the baseline expectations of the outcome. They may be harder to specify than those in Definition 1. This is a classic reason to focus on the ratio scale in sensitivity analysis (Cornfield et al. 1959; Poole 2010; Ding and VanderWeele 2014).

Sjölander et al. (2022) generalized the sensitivity parameters on the difference scale to allow for the transformation of the conditional means by a smooth link function and proposed an outcome regression approach to conduct sensitivity analyses with generalized linear models.

We can extend the current sensitivity analysis framework to analyze a more general class of nonlinear causal parameters, $g(\mu_{1},\mu_{0})$, a general function of the marginal means of the potential outcomes. Previous results focus on the special case when $g(\mu_{1},\mu_{0})=\mu_{1}-\mu_{0}$. More generally, many other function forms are of interest. Specifically, when the potential outcomes are binary, the causal risk ratio and the causal odds ratio, $\textsc{rr}=\mu_{1}/\mu_{0}$ and $\textsc{or}=\{\mu_{1}/(1-\mu_{1})\}/\{\mu_{0}/(1-\mu_{0})\}$ are two canonical causal parameters. We can utilize the nonparametric identification and proposed estimators for $\mu_{1}$ and $\mu_{0}$ and construct estimators for a general causal parameter by simply plugging in estimators $\hat{\mu}_{z}^{*}$ into function $g(\mu_{1},\mu_{0})$ for $*\in\{\textup{pred, prod, ht, dr}\}$ and $z=0,1$. The plug-in outcome and inverse probability weighting estimators are straightforward. The plug-in estimator $g(\hat{\mu}_{1}^{\textup{dr}},\hat{\mu}_{0}^{\textup{dr}})$ remains to be doubly robust and achieves the semiparametric efficiency bound. We also implement estimators for rr, or, and their logarithms in the R package saci.