INTRODUCING PROXIMAL CAUSAL INFERENCE FOR EPIDEMIOLOGISTS
Paul N. Zivich, Stephen R. Cole, Jessie K. Edwards, Grace E. Mulholland, Bonnie E. Shook‐Sa, Eric J. Tchetgen Tchetgen
Abstract
Causal inference with observational data has generally proceeded under the assumption of conditional exchangeability. That is, the action (e.g., treatment, exposure, intervention) is independent of the potential outcomes, conditional on a set of covariates (1). However, exchangeability is often questionable. Miao et al. (3) have proposed an alternative approach to identification, a generalization of previous work by Kuroki and Pearl (2), which allows for unmeasured confounding with particular causal structures. Specifically, if there exists a measured variable that is a potential cause of the action and unrelated to the outcome except through measured confounders and a known but unmeasured confounder (i.e., treatment proxy), and another measured variable that is a potential cause of the same outcome and unrelated to the action except through measured confounders and the same unmeasured confounder (i.e., outcome proxy), then the average causal effect (ACE) can be identified nonparametrically under a set of sufficient conditions. Here, we briefly introduce proximal causal inference, where “proximal” denotes that the pair of measured variables is a consequence of the unmeasured confounder, to the epidemiology community and demonstrate its application using a simulation study. To motivate our simulation, we draw from the work of Lu et al. (4). We wish to estimate the ACE of a new treatment (|$A=1$|), versus standard treatment (|$A=0$|), on 12-month human immunodeficiency virus viral load (|$Y$|). The population ACE can be expressed as |$E\!\left[{Y}^1-{Y}^0\right]$|, where |${Y}^a$| is the potential outcome under treatment |$a$|. We consider the confounders age (|$X$|) and CD4 cell count (|$U$|), where age was observed and the true CD4 count was unobserved. The other observed variables, self-rated health (|$W$|) and measured CD4 cell count (|$Z$|), are potential proxies, where measured CD4 count is a treatment proxy (i.e., a measure of the true CD4 count that features in treatment guidelines) and self-rated health is an outcome proxy (i.e., it is affected by true CD4 count and subsequently affects viral load). Several possible variations on the causal structure are depicted in Figure 1. Directed acyclic graphs for the 3 data-generating scenarios considered. |$A$|, action of interest; |$Y$|, outcome of interest; |$Z$|, action or treatment proxy; |$W$|, outcome proxy; |$X$|, traditional observed confounder; |$U$|, unobserved confounder. A) In scenario 1, both standard and proximal g-computation are expected to be unbiased, since |$\left\{X,W\right\}$| and |$\left\{X,W,Z\right\}$| block all backdoor paths. B) In scenario 2, standard g-computation is expected to be biased but proximal g-computation is expected to be unbiased. C) In scenario 3, both standard and proximal g-computation are expected to be biased. Notice that this model does not include the treatment proxy. Since the outcome model is a linear regression without interaction terms between |$A$| and any covariate, the ACE is estimated by |${\hat{\mathrm{\alpha}}}_1$|, the least-squares estimate of |${\mathrm{\alpha}}_1$|. In this simple case where linear regression is used to model the conditional mean values of both |${W}_i$| and |${Y}_i$|, existing software for 2-stage least-squares estimation can be used for point and variance estimation while allowing for the outcome proxy model to be misspecified (5). The described proximal g-computation procedure can be generalized. If the outcome model above had included interaction terms, the following modification could have been used to estimate the ACE (5). The predicted values of the outcome, |${\hat{Y}}_i^a$|, are generated from |$\hat{\boldsymbol{\mathrm{\alpha}}}$| by setting |${A}_i=a$| for all records, where |$a$| is in |$\left\{0,1\right\}$|; then the ACE is estimated as the mean difference between the predictions, |$\frac{1}{n}{\sum}_{i=1}^n\left({\hat{Y}}_i^1-{\hat{Y}}_i^0\right)$|. Proximal g-computation can also be applied with categorical or count outcome proxies. Discrete outcomes are also possible, with a nonlinear regression model for the outcome (e.g., logistic regression), but numerical integration or Monte Carlo simulation is needed for estimation (5). When multiple outcome proxies are present, all outcome proxies must be modeled separately, and all predicted outcome proxy values should be included in the outcome model. The first model corresponds to a minimally sufficient adjustment set for the ACE in Figure 1A and is expected to have the smallest asymptotic variance (8). The latter also includes |${Z}_i$|. As with proximal g-computation, the ACE is the coefficient for treatment in the case of linear regression with no interaction terms for |$A$| (i.e., |${\mathrm{\beta}}_1$| in the first model or |${\mathrm{\gamma}}_1$| in the second). We consider 3 scenarios (Figure 1, Web Appendix 2). In all scenarios, age affects treatment, viral load, CD4 cell count, measured CD4 cell count, and self-rated health. In scenario 1 (Figure 1A), CD4 count only affects treatment through measured CD4 count and only affects viral load through self-rated health (i.e., there is no direct effect of |$U$| on |$A$| or |$Y$|). In this case, both standard approaches and proximal causal inference are identified and expected to allow for unbiased estimation of the ACE. In scenario 2, CD4 count directly affects treatment and viral load (Figure 1B), as well as the indirect effects of scenario 1. Here, the standard approaches are expected to be biased (because treatment and potential outcomes are not independent given the measured covariates), but the ACE is identified with proximal causal inference and expected to allow for consistent estimation. Finally, in scenario 3, measured CD4 count affects both treatment and viral load directly, meaning |$Z$| would no longer be a valid treatment proxy variable (Figure 1C). Therefore, the ACE is not identified under any of the approaches considered, and all estimators ought to yield biased estimates. To assess the performance of the approaches, we compared the bias, empirical standard error, ratio of the estimated standard error to the empirical standard error, root mean squared error, and 95% confidence interval coverage for |$n=500$| with 4,000 iterations (9). For both proximal and standard g-computation, M-estimation was used, with the standard error of the estimated ACE estimated by the empirical sandwich variance estimator (10). Estimated standard errors were used to construct Wald-type 95% confidence intervals. All simulations were conducted in Python (with NumPy (11), SciPy (12), pandas (13), and delicatessen (14)) and independently replicated in both SAS (with IML) and R (with tidyverse (15) and rootSolve (16)). Simulation results are provided in Table 1 and Web Tables 1 and 2. In scenario 1, standard and proximal g-computation provided approximately unbiased estimates of the ACE of |$A$| on |$Y$| and approximate 95% confidence interval coverage. However, the root mean squared error indicated that the standard g-computation approach with the minimally sufficient adjustment set was preferred (due to the smaller empirical standard error). In short, a penalty in terms of the estimated variance was observed for proximal g-computation when accounting for a nonexistent bias. In scenario 2, standard g-computation was biased but proximal g-computation was approximately unbiased with appropriate coverage, as expected. Proximal g-computation was also superior in terms of root mean squared error. In scenario 3, all estimators were biased. Simulation Results From a Comparison of Standard g-Computation With Proximal g-Computation for Estimation of the Average Causal Effecta Abbreviations: ESE, empirical standard error; RMSE, root mean squared error; SER, standard error ratio. a Simulation results were derived using Python software (11–14). b Bias was defined as the estimated average causal effect minus the true average causal effect. c The ESE was calculated as the standard deviation of the estimated average causal effects. d The SER was calculated as the average estimated standard error divided by the ESE. The standard error was estimated by the sandwich variance estimator. SER > 1 indicates that the variance estimator is overestimating the variance of the estimator, and SER < 1 indicates that the variance estimator is underestimating the variance of the estimator. e The RMSE was calculated as the square root of the sum of the squared bias and squared ESE. f Coverage was defined as the proportion of Wald-type 95% confidence intervals that contained the true average causal effect. Wald-type 95% confidence intervals were calculated using the sandwich variance estimate. g Minimal standard g-computation is based on the minimal adjustment set expected to result in the smallest asymptotic variance (i.e., |$\left\{X,W\right\}$|). Standard g-computation included all measured confounding variables (i.e., |$\left\{X,Z,W\right\}$|). Proximal causal inference broadens the set of structures under which the ACE can be identified and estimated. As with other approaches to causal inference with observational data, proximal causal inference relies on unverifiable assumptions to make progress (17). When those assumptions are violated, estimators are biased, as shown in the final scenario with an invalid treatment proxy. For selecting valid treatment and outcome proxy variables, work on proxy variable selection provided in the context of negative controls provides some considerations (18, 19). When considering proximal causal inference, it is also important to note that |$U$| represents a known unmeasured confounder and not an unknown unmeasured confounder. The identification assumptions underlying proximal causal inference stipulate certain conditions to hold true for |$U$|, which are questionable if |$U$| were to represent an unknown unmeasured confounder. Specifically, claims about the completeness condition (i.e., that both proxy variables vary sufficiently relative to |$U$|’s variability) are difficult without some idea about the source of unmeasured confounding. Any identification and estimation strategy should be motivated by the best available substantive knowledge. Proximal causal inference has been extended in several directions. Here, estimation relied on parametric models, which require assumptions regarding correct specification of those models. Other work by Miao et al. (3), Cui et al. (20), and Ghassami et al. (21) further extends proximal causal inference by proposing semiparametric and nonparametric estimators that obviate the need for parametric models or doubly robust methods that provide robustness against partial misspecification (5, 22). Furthermore, there are extensions of proximal causal inference for other parameters, sources of systematic error, and data structures, including mediation (23), time-varying confounding (24), interference (25), front-door criterion with hidden mediators (21), synthetic controls (26), test-negative vaccine effectiveness (27), and data fusion (28). Given pervasive measurement error, perhaps most (if not all) measured confounders are proxies for underlying common causes; yet it remains to be seen how often a pair of proxy variables can be leveraged to improve causal effect estimation in epidemiology. Applied examples of proximal causal inference may help epidemiologists to better assess the benefits of proximal causal inference. This work was supported in part by National Institutes of Health grants T32-AI007001 (P.N.Z.), R01-AI157758 (S.R.C., J.K.E., B.E.S.-S., G.E.M.), and P30-AI050410 (S.R.C., B.E.S.-S.). All data were simulated from the data-generating mechanisms shown in Web Appendix 2. Corresponding software code is available at https://github.com/pzivich/publications-code. The views expressed in this article are those of the authors and do not reflect those of the National Institutes of Health. Conflict of interest: none declared.