Single Proxy Control

Unmeasured confounding is a well-known threat to valid causal inference from observational data. An approach that is sometimes used in practice to assess residual confounding bias, is to check whether known null effects can be recovered free of bias, by evaluating whether the exposure or treatment of interest is found to be associated with a so-called negative control outcome (NCO), upon adjusting for measured confounders (Rosenbaum, 1989; Lipsitch et al., 2010; Shi et al., 2020). An observed variable is said to be a valid negative control outcome or more broadly, an outcome confounding proxy, to the extent that it is associated with hidden factors confounding the exposure-outcome relationship in view, although not directly impacted by the exposure. Therefore, an NCO which is empirically associated with the exposure, might suggest the presence of residual confounding. In the event such an association is present, a natural question is whether the negative control outcome can be used for bias correction.

The most well-established NCO approach for debiasing observational causal effect estimates is the difference-in-differences approach (DiD) (Card and Krueger, 1994; Lechner, 2011; Caniglia and Murray, 2020). In fact, DiD may be viewed as directly leveraging the pre-treatment outcome as an NCO since it cannot logically be causally impacted by the treatment. Identification then follows from an additive equi-confounding assumption that the unmeasured confounder association with the post-treatment outcome of interest matches that with the pre-treatment outcome on the additive scale (Sofer et al., 2016). The baseline outcome in DiD is thus implicitly assumed to be a valid NCO, and equi-confounding is equivalent to the so-called parallel trends assumption, that the average trends in treatment-free potential outcomes for treatment and untreated units are parallel. In practice, equi-confounding or equivalently parallel trends may not be reasonable for a number of reasons, including if the outcome trend is also impacted by an unmeasured common cause with the treatment. Furthermore, additive equi-confounding may not be realistic as a broader debiasing method in non-DiD settings where the NCO is not necessarily a pre-treatment measurement of the outcome of interest, but is instead a post-treatment measurement of a different type of outcome (and might therefore have support on a different scale than the outcome of interest has).

To address these potential limitations of additive equi-confounding, Tchetgen Tchetgen (2013) introduced the Control Outcome Calibration Approach (COCA) as a simple yet formal counterfactual NCO approach to debias causal effect estimates in observational analyses. At its core, COCA essentially treats the NCO variable as a proxy measurement for the treatment-free potential outcome, which therefore is associated with the latter, and which becomes independent of the treatment assignment mechanism, upon conditioning on the treatment-free counterfactual outcome. As the treatment-free potential outcome can be viewed as an ultimate source of unmeasured confounding, this assumption formalizes the idea that, as a relevant proxy for the source of residual confounding, the NCO would be made irrelevant for the treatment assignment mechanism if one were to hypothetically condition on the underlying potential outcome.

For identification and inference for a continuous outcome, the original COCA approach of Tchetgen Tchetgen (2013) involves the correct specification of a regression model for the NCO, conditional on the treatment-free potential outcome and measured confounders, together with a rank-preserving structural model which effectively assumes a constant individual-level treatment effect. In this paper, we develop a nonparametric COCA identification framework for the average causal effect for the treated, which equally applies irrespective of the nature of the primary outcome, whether binary, continuous, or polytomous. Importantly, as we show, the proposed COCA identification framework completely obviates the need for rank preservation, therefore accommodating an arbitrary degree of effect heterogeneity across units. Relatedly, an alternative counterfactual approach named proximal causal inference has recently developed in causal inference literature (Miao et al., 2016, 2018; Tchetgen Tchetgen et al., 2024), which leverages a pair of negative treatment and outcome control variables, or more broadly treatment and outcome confounding proxies, to nonparametrically identify treatment causal effects subject to residual confounding without invoking a rank-preservation assumption. Importantly, while proximal causal inference relies on two proxies for causal identification, in contrast, COCA is a single-proxy control approach which therefore may present practical advantages. For estimation and inference, we introduce three strategies to implement COCA which improve on prior methods: (i) an extended propensity score approach, (ii) a so-called outcome calibration bridge function approach, and (iii) a doubly robust approach which carefully combines approaches (i) and (ii) and remains unbiased, provided that either approach is also unbiased, without necessarily knowing which method might not be unbiased. Finally, we illustrate the methods with an application evaluating the causal effect of a Zika outbreak on birth rate in Brazil, and we conclude with possible extensions to our methods and a brief discussion.

Consider an observational study where, as represented in Figure 1, one has observed an outcome variable $Y$, a binary treatment $A$ whose causal effect on $Y$ is of interest, and measured pre-treatment covariates $X$. We are concerned that as displayed succinctly in the figure with the bow arc, the association between $A$ and $Y$ is confounded by hidden factors.

Throughout, $Y^{a}$ denotes the potential outcome or counterfactual, had possibly contrary to fact, the exposure been set to $a$ by an external hypothetical intervention. Furthermore, throughout, we also make the consistency assumption:

Hereafter, we aim to make inferences about the causal effect of treatment on the treated (ETT), denoted by $\psi^{*}=E\big{(}Y^{a=1}-Y^{a=0}\,\big{|}\,A=1\big{)}$. Under consistency, $\psi_{1}^{*}=E\big{(}Y^{a=1}\,\big{|}\,A=1\big{)}$ is identified by $E\big{(}Y\,\big{|}\,A=1\big{)}$; to identify the counterfactual mean $\psi_{0}^{*}=E\big{(}Y^{a=0}\,\big{|}\,A=1\big{)}$ requires additional assumptions. Standard methods often resort to the no unmeasured confounding assumption, i.e., $Y^{a=0}\,\rotatebox[origin={c}]{90.0}{$\models$}\,A\,\big{|}\,X$, a strong assumption we do not make. Instead, we suppose that one has measured a valid NCO $W$, possibly multi-dimensional, which is known a priori to satisfy the following conditions:

Assumption 2-(i) encodes the key assumption of a known null causal effect of the treatment on the NCO in potential outcome notation. Assumption 2-(ii) encodes that $W$ is relevant for predicting the treatment-free potential outcome of interest. Assumption 2-(iii) states that $W$ is independent of $A$ conditional on treatment-free potential outcome and covariates. These conditions formally encode the assumption that $W$ is a valid proxy for the treatment-free potential outcome, a source of residual confounding bias; $W$ is only associated with the treatment mechanism to the extent that it is associated with the confounding mechanism captured by the potential outcome. We illustrate these NCO assumptions with the causal graph displayed in Figure 2. The thick arrows in the graph indicate the deterministic relationship defining the observed outcome in terms of potential outcomes and treatment variables by the consistency assumption. The missing arrows on the graph formally encode the core conditional independence conditions implied by Assumption 2.

It is enlightening to consider a data generating mechanism that is compatible with Assumption 2. As an example, we consider the following latent variable model for a continuous outcome:

where $U$ is a continuously distributed unobserved variable and $h_{y}(u,x)$ is a function that is strictly monotone in $u$ for all $x$, but otherwise completely unrestricted. In the Supplementary Material A.1, we establish that expressions (1a)-(1b) imply Assumption 2. Expression (1a) means that the treatment-free outcome is a monotonic transformation of an unobserved variable $U_{0}$, a specific instance of a so-called changes-in-changes model (Athey and Imbens, 2006). Although the latter would also assume under the causal graph in Figure 2 that $W=h_{w}(U_{0},X)$ where $h_{w}(u,x)$ is strictly monotone in $u$ for all $x$, an assumption we do not make. Expression (1b) corresponds to Assumption 2-(ii) and (iii). In fact, our formulation accommodates an additional measurement error in $W$, say $\epsilon_{w}$, so that $W=h_{w}(U_{0},X,\epsilon_{w})$ where $h_{w}(u,x,\epsilon_{w})$ varies in $u$ (potentially non-monotonic) and $\epsilon_{w}\,\rotatebox[origin={c}]{90.0}{$\models$}\,A\,\big{|}\,(U_{0},X)$, (1b) is satisfied. Figure 3 provides a graphical representation compatible with, although not necessarily with, expressions (1a)-(1b).

For identification and estimation in the case of a continuous outcome, Tchetgen Tchetgen (2013) further assumed the rank-preserving structural model:

which, by consistency, implies a constant individual-level causal effect $\psi^{*}=Y^{a=1}-Y^{a=0}$. Under this model, he noted that upon defining $Y(\psi)=Y-\psi A$, then one can deduce from Assumption 2 that $W\,\rotatebox[origin={c}]{90.0}{$\models$}\,A\,\big{|}\,(Y(\psi),X)$ if and only if $\psi=\psi^{*}$, in which case, given (2), $Y^{a=0}=Y(\psi^{*})=Y-\psi^{*}A$ which motivates a regression-based implementation of COCA, that entails searching for the parameter value of $\psi=\psi^{*}$ such that $E\big{\{}W\,\big{|}\,A,Y(\psi),X\big{\}}=E\big{\{}W\,\big{|}\,Y(\psi),X\big{\}}$. A straightforward implementation of the approach uses linear models whereby for each value of $\psi$ on a sufficiently fine grid, one obtains an estimate of the regression model $E\big{\{}W\,\big{|}\,Y(\psi),X\big{\}}=\beta_{1}+\beta_{2}A+\beta_{3}Y(\psi)+\beta_{4}X$ using ordinary least squares (OLS), with estimated coefficients $\big{(}\widehat{\beta}_{1}(\psi),\widehat{\beta}_{2}(\psi),\widehat{\beta}_{3}(\psi),\widehat{\beta}_{4}(\psi)\big{)}$. Then a 95% confidence interval for $\psi^{*}$ consists of all values of $\psi$ for which a valid test of the null hypothesis $\beta_{2}(\psi)=0$ fails to reject at the 0.05 type 1 error level. Such hypothesis test might be performed by verifying whether the interval $\widehat{\beta}_{2}(\psi)\pm 1.96\widehat{\mathrm{SE}}\big{(}\widehat{\beta}_{2}(\psi)\big{)}$ covers $0,$ with $\widehat{\mathrm{SE}}\big{(}\widehat{\beta}_{2}(\psi)\big{)}$ the OLS estimate of the standard error of $\widehat{\beta}_{2}(\psi)$. Tchetgen Tchetgen (2013) also describes a potentially simpler one-shot approach which fits a single regression $E\big{(}W\,\big{|}\,A,Y,X\big{)}=\beta_{1}+\beta_{2}^{*}A+\beta_{3}Y+\beta_{4}X$ via OLS where $\beta_{2}^{*}=-\beta_{3}/\psi_{0},$ in which case $\widehat{\psi}=-\widehat{\beta}_{2}^{*}/\widehat{\beta}_{3}$ where $\big{(}\widehat{\beta}_{2}^{*},\widehat{\beta}_{3}\big{)}$ are OLS estimates; a corresponding standard error estimator of $\widehat{\psi}$ is given in the Supplementary Material A.2 for convenience. Though practically convenient, validity of either approach relies on both correct specification of the linear model for $W$ given $(Y^{a=0},X)$, and on the rank-preserving structural model, which may be biologically implausible. In the following, we describe alternative methods aimed at addressing these limitations.

Our estimator is derived from the influence function IF in Result 5 and adopts the cross-fitting approach (Schick, 1986; Chernozhukov et al., 2018), which is implemented as follows. We randomly split $N$ study units, denoted by $\mathcal{I}=\{1,\ldots,N\}$, into $K$ non-overlapping folds, denoted by $\{\mathcal{I}_{1},\ldots,\mathcal{I}_{K}\}$. For each $k=1,\ldots,K$, we estimate the EPS and COCA confounding bridge functions using observations in $\mathcal{I}_{k}^{c}=\mathcal{I}\setminus\mathcal{I}_{k}$, and then evaluate the estimated nuisance functions using observations in $\mathcal{I}_{k}$ to obtain an estimator of $\psi^{*}$. We refer to $\mathcal{I}_{k}^{c}$ and $\mathcal{I}_{k}$ as the estimation and evaluation folds, respectively. To use the entire sample, we take the simple average of the $K$ estimators.

We introduce the following additional notation in order to facilitate the discussion. Let $\mathcal{H}(V)$ be the Reproducing Kernel Hilbert Space (RKHS) of $V$ endowed with a universal kernel function $\mathcal{K}$, such as the Gaussian kernel $\mathcal{K}$, i.e., $\mathcal{K}(v,v^{\prime})=\exp\big{\{}-\big{\|}v-v^{\prime}\|_{2}^{2}/\kappa\big{\}}$ where $\kappa\in(0,\infty)$ is a bandwidth parameter; see Chapter 4 of Steinwart and Christmann (2008) for the definition and examples of the universal kernal function. For each $k=1,\ldots,K$, let $\mathbbm{P}^{(-k)}(V)=|\mathcal{I}_{k}^{c}|^{-1}\sum_{i\in\mathcal{I}_{k}^{c}}V_{i}$ and $\mathbbm{P}^{(k)}(V)=|\mathcal{I}_{k}|^{-1}\sum_{i\in\mathcal{I}_{k}}V_{i}$. For a function $g(O)$, let $\big{\|}g\big{\|}_{P,2}=\big{[}E\big{\{}g^{2}(O)\big{\}}\big{]}^{1/2}$ be the $L_{2}(P)$-norm of $g$.

We estimate the EPS and COCA bridge functions by adopting a recently developed minimax estimation approach (Ghassami et al., 2022). We remark that other approaches (e.g., Mastouri et al. (2021)) can also be adopted with minor modification. Note that $\omega^{*}(Y,X)$ and $b^{*}(W,X)$ satisfy

Therefore, following Ghassami et al. (2022), minimax estimators of $\omega^{*}$ and $b^{*}$ are given by

where $\big{\|}\cdot\big{\|}_{\mathcal{H}}$ is an RKHS norm and $\lambda_{p}$, $\lambda_{\omega}$, $\lambda_{q}$, and $\lambda_{b}$ are positive regularization parameters.

We make a few remarks about the minimax estimation approach, of which details are relegated to the Supplementary Material A.6. First, despite the complicated formulas, closed-form representations of $\widehat{\omega}^{(-k)}$ and $\widehat{b}^{(-k)}$ are available from the representer theorem (Kimeldorf and Wahba, 1970; Schölkopf et al., 2001). Second, the bandwidth and regularization parameters can be selected via cross-validation. Lastly, $\omega^{*}$ may vary widely because it is a ratio of two probabilities. In such cases, the proposed minimax estimator may result in significantly small or negative estimates. To mitigate this issue, one may consider a practical approach to regularize the minimax estimator when it appears to be ill-behaved.

Using the minimax estimators of the nuisance functions, a semiparametric estimator $\widehat{\psi}$ of $\psi^{*}$ is then obtained as follows:

Under regularity conditions, the semiparametric estimator $\widehat{\psi}$ is consistent and asymptotically normal for $\psi^{*}$.

Assumption 5-(i) states that nuisance functions and the corresponding estimators are uniformly bounded. Assumption 5-(ii) states that the estimated nuisance functions are consistent for the true nuisance functions in the $L_{2}(P)$ norm sense. Assumption 5-(iii) states that the cross-product rate of nuisance function estimators are $o_{P}(N^{-1/2})$. Assumption 5-(iii) is satisfied if $b^{*}$ and $\omega^{*}$ are sufficiently smooth, the conditional expectation operators $f(W,X)\mapsto E\{f(W,X)\,\big{|}\,Y,A=0,X\}$ and $f(Y,X)\mapsto E\{f(Y,X)\,\big{|}\,W,A=0,X\}$ are sufficiently smooth, and $\widehat{\omega}^{(-k)}$ and $\widehat{b}^{(-k)}$ are estimated over an RKHS with fast enough eigendecay; see Section 5 of Ghassami et al. (2022) for details. Importantly, if one nuisance function is estimated at sufficiently fast rates, the other nuisance function is allowed to converge at a substantially slower rate provided that the cross-products remain $o_{P}(N^{-1/2})$. This is an instance of the mixed-bias property described by Rotnitzky et al. (2020) and Ghassami et al. (2022). It is also worth highlighting the structure of the mixed bias in the current context which is the minimum of two product biases, each containing a bias term for a nuisance function and a projected bias term for the other nuisance function; this property was first reported in Ghassami et al. (2022) for a large class of functionals including ours.

Result 6 establishes that $\widehat{\psi}$ is consistent and asymptotically normal (CAN) for $\psi^{*}$.

Using the variance estimator $\widehat{\sigma}^{2}$, valid $100(1-\alpha)$% confidence intervals for the ETT are given by $\big{(}\widehat{\psi}+z_{\alpha/2}\widehat{\sigma}/\sqrt{N},\widehat{\psi}+z_{1-\alpha/2}\widehat{\sigma}/\sqrt{N}\big{)}$ where $z_{\alpha}$ is the $100\alpha$th percentile of the standard normal distribution. Alternatively, one may construct confidence intervals using the multiplier bootstrap (van der Vaart and Wellner, 1996, Chapter 2.9); see Section A.6 for details.

Lastly, the cross-fitting estimator depends on a specific sample split, and thus, may produce outlying estimates if some split samples do not represent the entire data. To mitigate this issue, Chernozhukov et al. (2018) proposes to use median adjustment from multiple cross-fitting estimates; the detail can be found in the Supplementary Material A.6.

The Zika virus, which can be transmitted from a pregnant woman to her fetus, can cause serious brain abnormalities, including microcephaly (i.e., an abnormally small head) (Rasmussen et al., 2016). Brazil is one of the countries hardest hit by the Zika virus. In particular, the outbreak in 2015 resulted in over 200,000 cases in Brazil by 2016 (Lowe et al., 2018). As a result, many prior works (Castro et al., 2018; Diaz-Quijano et al., 2018; Taddeo et al., 2022; Tchetgen Tchetgen et al., 2024) asked whether the Zika virus outbreak caused a drop in birth rates.

We re-analyzed the dataset analyzed in Taddeo et al. (2022) and Tchetgen Tchetgen et al. (2024). In the dataset, we focused on 673 municipalities in two states of Brazil, Pernambuco and Rio Grande do Sul, which are northeastern and southernmost states. Out of the 1248 cases of microcephaly that occurred in Brazil by November 28, 2015, 51.8% (646 cases) were reported in Pernambuco (PE), less than 10 cases of Zika-related microcephaly were reported in Rio Grande do Sul (RS) (Gregianini et al., 2017), which shows that PE was severely impacted by the Zika virus outbreak, while RS was minimally affected. Based on their epidemiologic histories, we defined 185 and 488 municipalities in PE and RS as treated and control groups, respectively.

For each municipality, we included the following variables in the analysis. As pre-treatment covariates, we included municipality-level population size, population density, and proportion of females measured in 2014. We used the post-epidemic municipality-level birth rate in 2016 as the outcome $Y$, where the birth rate is defined as the total number of live human births per 1,000 persons. We used the pre-epidemic municipality-level birth rates in 2013 and 2014 as the outcome proxies (i.e., NCO), denoted by $W_{1}$ and $W_{2}$, respectively. To be valid proxies, the birth rates in 2013 and 2014 must satisfy Assumption 2: (i) birth rates in 2013 and 2014 cannot be causally impacted by the Zika virus epidemic which occurred in 2015, (ii) birth rates in 2013 and 2014 are correlated with what the birth rate in 2016 would have been had there not been a Zika virus epidemic, and (iii) birth rates in 2013 and 2014 are independent of a municipality’s Zika epidemic status, upon conditioning on its Zika virus epidemic-free potential birth rate in 2016. The first two conditions are uncontroversial, while the third condition largely relies on the extent to which pre-epidemic birth rates can accurately be viewed as a proxy for the counterfactual birth rate had the pandemic not occurred, and as such would not further be predictive of whether the municipality experienced a high rate of Zika virus incidence, conditional on the region’s epidemic-free counterfactual birth rate in 2016. Although one might consider this last assumption reasonable, ultimately, it is empirically untestable without making an alternative assumption. Nevertheless, in the Supplementary Material A.8, we describe a straightforward sensitivity analysis to evaluate the extent to which violation of the assumption might impact inference.

Using the dataset, we estimated $\psi^{*}=E\big{(}Y^{a=1}-Y^{a=0}\,\big{|}\,A=1\big{)}$, i.e., the difference between the observed average birth rate of Pernambuco and a forecast of what it would have been had the Zika outbreak been prevented. Therefore, the ETT quantifies the average treatment effect of the Zika outbreak on the birth rate within the Pernambuco region. Of note, the crude estimand $E(Y\,\big{|}\,A=1)-E(Y\,\big{|}\,A=0)$ was estimated to be equal to $3.384$, suggesting that municipalities in the PE region (with higher incidence of Zika virus) experienced a higher birth rate than RS regions in 2016 during the Zika virus outbreak. An immediate concern is that this crude association between $A$ and $Y$ might be subject to significant confounding bias, leading us to conduct two separate analyses geared at addressing residual confounding bias; the proposed COCA methods, which we compared with a standard difference-in-differences analysis. Thus, we estimate the ETT using the approach outlined in Section 4 where the NCOs are specified as either (i) $W_{1}$, birth rate in 2013, or (ii) $W_{2}$, birth rate in 2014, or (iii) $(W_{1},W_{2})$. For comparison, we also obtained doubly-robust parametric estimators of the ETT using the three NCO specifications; see the Supplementary Material A.7 for details on how these estimators were constructed.

Table 1 summarizes corresponding results. We find that the six COCA estimates vary between $-1.833$ and $-2.410$, meaning between $1.833$ and $2.410$ birth per 1,000 persons were reduced in PE due to the Zika virus outbreak, an empirical finding better aligned with the scientific hypothesis that Zika may likely adversely impact the birth rates of exposed populations. Compared to the crude estimate of 3.384, the negative effect estimates indeed provide compelling evidence of potential confounding. We also obtain an estimate using the difference-in-difference estimator under a standard parallel trends assumption (e.g., Card and Krueger (1994); Angrist and Pischke (2009)), which yields a considerably smaller effect estimate varying between $-1.156$ and $-1.041$; noting that the DiD estimator requires the assumption of equi-confounding of the $A-Y$ association in the pre- and post-periods, while our proposed estimator does not (but instead requires conditions (i)-(iii) outlined above). Regardless of the estimator, all estimates appear to be consistent with the anticipated adverse causal impact of the Zika Virus epidemic. Consequently, we conclude that based on inferences aimed at accounting for confounding (DiD and COCA), the Zika virus outbreak likely led to a decline in the birthrate of affected regions in Brazil, which agrees with similar findings in the literature (Castro et al., 2018; Diaz-Quijano et al., 2018; Taddeo et al., 2022; Tchetgen Tchetgen et al., 2024).

We have described a COCA nonparametric identification framework, therefore extending previous results of Tchetgen Tchetgen (2013) to a more general setting accommodating outcomes of arbitrary nature and obviating the need for an assumption of constant treatment effects, i.e. rank preservation. We have proposed three estimation strategies including a doubly robust method which has appealing robustness properties. Interestingly, the COCA central identifying assumption, that conditioning on the treatment-free counterfactual would in principle shield the treatment assignment from any association with the NCO is isomorphic to an analogous assumption in the missing data literature where an outcome might be missing not at random; however, a fully observed so-called shadow variable (the missing data analog of an NCO) reasonably assumed to be conditionally independent of the missing data process given the value of the potentially missing outcome. For example, Zahner et al. (1992) considered a study of the children’s mental health evaluated through their teachers’ assessments in Connecticut. However, the data for the teachers’ assessments are subject to nonignorable missingness. As a proxy of the teacher’s assessment, a separate parent report is available for all children in this study. The parent report is likely to be correlated with the teacher’s assessment, but is unlikely to be related to the teacher’s response rate given the teacher’s assessment and fully observed covariates. Hence, the parental assessment is regarded as a shadow variable for the teacher’s assessment in this study. The literature on shadow variables is fast-growing (d’Haultfoeuille, 2010; Kott, 2014; Wang et al., 2014; Miao and Tchetgen Tchetgen, 2016; Li et al., 2023; Miao et al., 2023), the methods developed in this paper have close parallels to shadow variable counterparts in this literature. This connection to shadow variables is particularly salient when the COCA confounding bridge function is not uniquely defined, which can easily occur for instance when the shadow variable (or analogously the NCO) is multivariate, therefore significantly complicating inference. Fortunately, the methods developed by Li et al. (2023) for the analogous shadow variable setting directly apply to the corresponding COCA setting and thus provide a complete solution for identification and inference for the average treatment effect for the treated without relying on completeness conditions nor on unique identification of either the EPS or the COCA confounding bridge function. We refer the interested reader to this latter work for further details. It is worth noting that the doubly robust estimator proposed in this paper appears to be completely new, and different from those of Miao and Tchetgen Tchetgen (2016), Li et al. (2023), and Miao et al. (2023) and therefore may also be of use in shadow variable applications. Likewise, the doubly robust estimators proposed in the latter works can equally be applied to the current COCA setting as an alternative inferential approach.

Additionally, as mentioned in the previous Section, the key assumption that conditioning on the treatment-free potential outcome, would in principle make the NCO or outcome proxy irrelevant to treatment mechanism is ultimately untestable, and may in certain settings not hold exactly. In fact, this would be the case if the NCO were in fact explicitly used in assigning the treatment in which case the assumption might be violated. In order to address such eventuality, the analyst might consider several candidate proxies/NCOs when available, and may even perform an over-identification test, by inspecting the extent to which the estimated causal effect depends on the choice of proxy. Alternatively, a sensitivity analysis might also be performed to evaluate the potential impact of a hypothesized departure from the assumption. In the context of the Zika virus application, an over-identification test and a sensitivity analysis were carried out as illustrative examples, with the corresponding results and discussion provided in the Supplementary Material A.8 and A.9.

Finally, as previously mentioned in the Introduction, COCA offers an alternative approach to proximal causal inference for debiasing observational estimates of the ETT by leveraging negative control or valid confounding proxies. A key difference highlighted earlier between these two frameworks is that COCA relies on a single valid NCO which directly proxies the treatment-free potential outcome, while proximal causal inference requires both valid NCO and negative control treatment variables that proxy an underlying unmeasured confounder. Importantly, COCA takes advantage of the fact that the treatment-free potential outcome is observed in the untreated, while in proximal causal inference, the unmeasured confounders for which proxies are available, are themselves never observed, arguably a more challenging identification task. Despite the practical advantage of needing one rather than two proxies, it is important to note that though COCA identifies the ETT, it fails to nonparametrically identify the population average treatment effect (ATE), without an additional assumption. In contrast, proximal causal inference provides nonparametric identification of both causal parameters and thus can be interpreted as providing richer identification opportunities. A key reason for this difference in the scope of identification is the fact that in the current paper, we have emphasized an interpretation of the NCO as a proxy for the treatment-free potential outcome but not for the potential outcome under treatment, in the sense that under our conditions, it must be that the treatment-free potential outcome does not only shield the treatment from the NCO, but also shields the potential outcome under treatment from the latter. Two potential strategies to recover COCA identification of the population ATE, might be either (i) evoke a rank-preservation assumption which if appropriate would imply that the ETT and the ATE are equal (this is the assumption made in Tchetgen Tchetgen (2013)); or (ii) identify a second proxy NCO which is a valid proxy for the counterfactual outcome under treatment. The second condition would be needed if $Y^{a=1}$ can also be viewed as a hidden confounder. By a symmetry argument, one can show that (ii) in fact would provide identification of the average counterfactual outcome under treatment for the untreated. A weighted average of both counterfactual means would then provide identification of the ATE. Details are not provided, but can easily be deduced from the presentation.

The authors would like to thank James Robins, Thomas Richardson, and Ilya Shpitser for helpful discussions.

The data and the analysis R code are accessible on the GitHub repository located at http://github.com/qkrcks0218/SingleProxyControl.