Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
51,141 characters · 10 sections · 69 citation commands
$$-GNF: A Copula-based Sensitivity Analysis to Unobserved Confounding Using Normalizing Flows
\looseness=-1 Epidemiologists, sociologists, economists, and other applied scientists, often leverage randomized controlled trials (RCTs), as RCTs provide the safest methodological route to disentangle cause and effect. RCTs are the gold-standard as they require the least assumptions wright1921correlationcausation, fisher1936designofexperiments, cox1958planning, Imbens2015CIinSocialscience. Most importantly, by randomizing which experimental subjects (e.g., people, villages, schools) should take the treatment and which subjects should abstain, an RCT provides unconditional ignorability or exchangability or unconfoundedness. Unconfoundedness implies no unobserved confounders in the causal system of interest: no unobserved common causes of the treatment and the outcome. When unconfoundedness is satisfied, scholars can calculate the causal effect of interest from collected data; that means that the causal quantity is identified ROBINS1986gcom, Rubin1990unconfoundedness, robins2008ipw, pearl2009causality, hernan2009ipw. However, despite the importance of the RCT design, it remains infeasible for a slew of applied settings. It may be costly to implement (e.g., testing a population-wide medicine); it may be unethical (e.g., testing a new drug); or, it may be impractical to implement (e.g., testing a social policy across the world). Therefore, applied researchers often rely on observational data -- which are often secondary data sources with no treatment randomization and where the experimenter had no control over the data generating process. Yet when using observational data, scholars make themselves susceptible for failing to satisfy the unconfoundedness assumption, even when some confounders are observed.
\looseness=-1 Because the unconfoundedness assumption is so critical and at the same time untestable in observational studies Rubin1990unconfoundedness, d2020accounting, methodologists (statisticians, computer scientist, and others) have developed various frameworks for stress testing how causal effect estimates change under varying strength of unconfoundedness failure. These sort of tests are named sensitivity analysis schlesselman1978assessing, manski1990sensitivitybounds, imbens2003sensitivity, brumback2004sensitivity, vanderweele2011bias, also known as bias analysis in epidemiology cornfield1959smoking, cochran1973controlling, rothman2008modern, lash2009applying. Nonetheless, existing sensitivity analysis frameworks are limited in at least three ways, and our proposed $\rho$-GNF improve on these limitations, thereby moving the state-of-the-art forward.
First, with the idea of interventional equivalence of structural causal models (SCMs), we propose a deep-learning method for sensitivity analysis based on the graphical normalizing flow (GNF) wehenkel2020GNF, causal graphical normalizing flow (c-GNF) balgi2022cgnf and Causal Normalizing Flow (CNF) JavaloySV23, because of GNF's attractive properties of non-linearity and invertibility for counterfactual inference and the similarities to the most well studied and popular elliptical copula, i.e., Gaussian copula. Hence, we aptly name the model $\rho$-GNF, where $\rho{ \in }[-1,+1]$ is the bounded sensitivity parameter of the Gaussian copula that controls the degree of unconfoundedness between the observed treatment and outcome. Unlike most sensitivity analysis methods where the sensitivity parameters are unbounded and difficult to specify and interpret, $\rho$ is bounded in the range $[-1,+1]$ and it represents the non-causal dependence between the treatment and outcome. Second, we show that with interventional equivalence, $\rho$-GNF enables us to estimate the interventional causal effects such as average causal effect (ACE) as a function of $\rho$. We call this the $\rho_{curve}$, and this curve enables us to identify the ACE bounds and analyze the confounding strength required to explain away the causal effect. Thus, the $\rho_{curve}$ enables us to provide bounds for the causal effect given a specific interval of $\rho{ \in }[-1,+1]$ that the domain expert considers appropriate. We also define $\rho_{value}$ as the value of $\rho$ that explains away the causal effect and demonstrate its similarities to the widely used E-value vanderweele2017sensitivity. Third, we empirically demonstrate sharp and narrower bounds compared to the widely popular assumption-free bounds through simulated as well as real-world experiments. Unlike existing sensitivity analysis methods that provide bounds either only for discrete or continuous outcomes, $\rho$-GNF accommodates both discrete and continuous outcomes.
\looseness=-1 Our work proceeds as follows. After briefly reviewing related literature in Section (ref), we define our notation and the main objective of ACE estimation, using the idea of interventional equivalence of SCMs with bivariate Gaussian copula in Section (ref). We perform the identification and estimation of the ACE and analyse the sensitivity to the different degree of unconfoundedness in Section (ref)-(ref), via the sensitivity parameter $\rho$. In Section (ref), we present our results with simulated as well as real-world data under different settings of outcome variable (i.e., continuous or binary or categorical), and compare them with the popular assumption-free (AF) bounds. Finally, in Section (ref), we conclude with discussing the key contributions of our $\rho$-GNF method in encouraging the use of sensitivity analysis when working with non-randomized observational data.
\looseness=-1 The sensitivity analysis literature can be roughly categorized into two streams: (i) identify the bounds of the causal effect as functions of some sensitivity parameters that encode the strength of the unobserved confounders robins1989analysis, manski1990sensitivitybounds, vanderweele2011bias, ding2016sensitivity, sjolander2020note, sjolander2021novel, pena2022simple; and (ii) identify how large the influence of the unobserved confounders needs to be to explain-away the causal effect imbens2003sensitivity, vanderweele2017sensitivity, veitch2020sense, sjolander2022values. While robins1989analysis, manski1990sensitivitybounds provide assumption-free (AF) bounds of the causal effect for binary outcome, more recent methods such as ilse2021efficient extend the bounds to categorical outcomes. Other methods provide bounds as functions of sensitivity parameters to be tuned by a domain expert sjolander2020note, sjolander2021novel, pena2022simple. cinelli2019sensitivity, cinelli2020making study sensitivity to unobserved confounding in a linear SCM setting. In contrast to the bounds stream, there are methods that fall under the explain-away stream. For example, imbens2003sensitivity, vanderweele2017sensitivity reason on the lines of the minimum strength of the unmeasured confounder that is needed, conditional on the measured confounders, to explain-away the estimated causal effect. Similar to imbens2003sensitivity and the E-value vanderweele2017sensitivity, the recently developed Austen plots veitch2020sense identify the influence of the confounding needed to explain a specific amount of bias in the causal effect estimate.
\looseness=-1 While there exists a wide spectrum of sensitivity analysis methods with unique advantages, they are not without limiting assumptions. For example, robins1989analysis, manski1990sensitivitybounds, sjolander2020note, sjolander2021novel, pena2022simple require the outcome variable being binary, and other methods assume a specific type of parametric model, e.g., cinelli2019sensitivity, cinelli2020making assume a linear parametric model. While some methods offer sharp ACE bounds robins1989analysis, manski1990sensitivitybounds, sjolander2020note, sjolander2021novel, pena2022simple, method such as vanderweele2017sensitivity may result in wider bounds than the AF bounds, as shown by ioannidis2019limitations and sjolander2020note. Some methods are exclusively suited for specific causal estimands, e.g., ACE or conditional ACE (CACE) or mediation effects tchetgen2012semiparametric, lindmark2018sensitivity. While robins1989analysis, manski1990sensitivitybounds offer no sensitivity parameters that can explain a certain causal effect, vanderweele2017sensitivity, veitch2020sense offer multiple parameters that are unbounded and hard to specify for the domain analyst ioannidis2019limitations.
To summarize, even though there exists several sensitivity analysis frameworks, there is still a lack of unifying method that is flexible enough that can suit many different types of observational data, with easy-to-use sensitivity parameters, and that can be applied to not only binary outcome variables but also categorical and continuous outcomes. Moreover, it's imperative to establish a method enabling researchers to specify distributional assumptions regarding the unobserved causes within the causal system under study. Such assumptions enable tighter ACE bounds, enhancing the certainty of outcomes. Our $\rho$-GNF method targets all these lacks. Our method use deep neural networks, allowing for maximum flexibility and non-linearity. $\rho$-GNF provides a single, bounded sensitivity parameter $\rho{ \in }[-1,+1]$ that is easily interpreted as the measure of non-causal association due to unobserved confounders. Thus, our method enhances an applied researcher's causal toolbox, and shows how deep learning can further causal inference.
\looseness=-1 Let us consider the standard canonical representation (in Chapter 3.4 of peters2017eci) for the structural causal model (SCM) wright1921correlationcausation, Haavelmo1943SEM, Goldberger1972SEMecon, Fienberg1975Introduction2SEM, TabarJYCFL2022ML4PP, where $A$ is the treatment (or cause) and $Y$ is the outcome (or effect), and $\varepsilon_A$ and $\varepsilon_Y$ are their respective unobserved causes such that
where $(\varepsilon_A, \varepsilon_Y)$ follows the joint CDF $\textcolor{orange}{\mathbb{F}_{\varepsilon_A,\varepsilon_Y}(\varepsilon_A,\varepsilon_Y)}$ and is not limited to the usual assumptions of standard normal or uniform random variables as in most common SCM definitions. Using the universality of the uniform (also know as the probability integral transform) angus1994probability, the noise variables $\varepsilon_A$ and $\varepsilon_Y$ of the SCM in Eq. (ref) can equivalently be written in terms of uniform variables $U_A$ and $U_Y$ in the interval $[0,1]$ resulting in the SCM for Figure (ref) (\textcolor{orange}{orange}+\textcolor{blue}{blue}) as
where $\mathbb{F}_{\varepsilon_A}$ and $\mathbb{F}_{\varepsilon_Y}$ respectively denote the marginal CDFs of $\varepsilon_A$ and $\varepsilon_Y$, and $(U_A, U_Y)$ follows the joint CDF $\textcolor{blue}{\mathbb{F}_{U_A,U_Y}(U_A,U_Y)}$ with uniform marginals in [0,1]. From the universality of the uniform, Eq. (ref) can further be simplified in terms of $\mathbb{F}_{A}$ and $\mathbb{F}_{Y|A}$ that denote the marginal CDFs of $A$ and $Y$ conditioned on $A$, respectively as below.
\looseness=-1 Further, we represent $U_A$ and $U_Y$ in Eq. (ref) as transformation of standard normal variables using the CDF of the standard normal $\Phi$, as shown in Figure (ref) (\textcolor{orange}{orange}+\textcolor{blue}{blue}+\textcolor{red}{red}), as
where $(Z_A, Z_Y)$ follows the joint CDF $\textcolor{red}{\mathbb{F}_{Z_A,Z_Y}(Z_A,Z_Y)}$ with standard normal marginals. Eqs. (ref)-(ref) represent observationally and interventionally equivalent SCMs from Figure (ref), but with different unobserved noises and corresponding joint CDFs. Observational equivalence means that the researcher specified model yields the same distribution as observed by the true (nature or God) data-generating process {$\mathbb{F}_{A,Y}(A,Y)$}; interventional equivalence means that the researcher specified model follows the same interventional distribution {$\mathbb{F}_{Y}(Y|do(a))$} and {$\mathbb{F}_{A}(A|do(y))$} mooij2016distinguishing. Since these noises (distributions) are unknown or unobserved, the completeness of $do$-calculus tian2002unconfoundedchildrencriteria, HuangV06docalculus, pearl2012docalculus states that the causal effects of interest are not identifiable from any of the equivalent SCMs in Eqs. (ref)-(ref), without further assumptions. To achieve causal effect/estimand identification and bounds, in the subsequent sections, we propose to use a bivariate Gaussian copula, to model the back-door non-causal noise distributions so that the causal effects can be parametrically estimated using deep-neural-network-inspired normalizing flows trained only on observational data.
A copula is a multivariate distribution function defined on the unit hypercube with uniform marginals sklar1959fonctions, sklar1973random. As the name suggests, a copula `ties' or `links' or `couples' a multidimensional joint distribution to its marginals nelsen2007introduction. From the result of Sklar's Theorem, we have that the bivariate joint CDF $\textcolor{blue}{\mathbb{F}_{U_A,U_Y}(U_A,U_Y)}$ in Eqs. (ref) and (ref) with uniform marginals in [0,1] can be represented using a bivariate copula $\textcolor{blue}{\textcolor{blue}{\mathbf{\mathbb{C}}}(U_A,U_Y)}$.
Eqs. (ref) and (ref) follow from the scale-invariance property of copula $\textcolor{blue}{\textcolor{blue}{\mathbf{\mathbb{C}}}(U_A,U_Y)}$ to strictly increasing transformations/CDFs $\textcolor{orange}{\mathbb{F}_{\varepsilon_A}}$, $\textcolor{orange}{\mathbb{F}_{\varepsilon_Y}}$, and $\textcolor{red}{\Phi}$. \looseness=-1 The unknown copula $\textcolor{blue}{\mathbf{\mathbb{C}}(U_A,U_Y)}$ in Figure (ref) essentially models the non-causal back-door association, and the degree of the non-causal association between $U_A$ and $U_Y$ may be quantified using measures of association such as spearman1987proof, spearman2010proof's $\rho_S$ or kendall1938new's $\tau_K$, where $\rho_S,\tau_K{\in}[-1,+1]$. From the scale-invariance property of $\rho_S$ and $\tau_K$ to the strictly increasing transformations $\textcolor{orange}{\mathbb{F}_{\varepsilon_A}}$, $\textcolor{orange}{\mathbb{F}_{\varepsilon_Y}}$, and $\textcolor{red}{\Phi}$, this measure of association between $U_A$ and $U_Y$ is the same between $\varepsilon_A$ and $\varepsilon_Y$, and between $Z_A$ and $Z_Y$, i.e., $\textcolor{orange}{\rho_S(\varepsilon_A, \varepsilon_Y}){=}\textcolor{blue}{\rho_S(U_A, U_Y)}{=}\textcolor{red}{\rho_S(Z_A,Z_Y)}{=}\textcolor{blue}{\rho_{S_{\textcolor{blue}{\mathbf{\mathbb{C}}}}}}$. In other words, the \textcolor{orange}{orange}, \textcolor{blue}{blue} and \textcolor{red}{red} back-door paths in Figure (ref) induce the same measure of non-causal association $\textcolor{blue}{\rho_{S_{\textcolor{blue}{\mathbf{\mathbb{C}}}}}}$ that can be represented using the copula $\textcolor{blue}{\textcolor{blue}{\mathbf{\mathbb{C}}}(U_A,U_Y)}$. This result intuitively follows as the SCMs are observationally and interventionally equivalent, i.e., same measures of total observed and causal associations are expected, implying the non-causal associations also to be the same. As discussed above, for identifying the causal effect of interest, it is sufficient if we observe/known/hypothesize the noises or the copula to adjust for the non-causal back-door path. However, the copula $\textcolor{blue}{\mathbf{\mathbb{C}}(U_A,U_Y)}$, although uniquely exists because of the strictly increasing continuous transformations $\textcolor{orange}{\mathbb{F}_{\varepsilon_A}}$, $\textcolor{orange}{\mathbb{F}_{\varepsilon_Y}}$, and $\textcolor{red}{\Phi}$, remains unknown and cannot be estimated from observational data as the noises are unobserved.
Since the copula $\textcolor{blue}{\mathbf{\mathbb{C}}(U_A,U_Y)}$ is unknown and unlearnable, it is inevitable to make assumptions about it to achieve causal effect identification. Specifically, $\textcolor{blue}{\mathbf{\mathbb{C}}(U_A,U_Y)}$ may be chosen from any of the vast families of copulas such as Archimedean copulas ling2020deep (Clayton, Frank, Gumbel, etc.), elliptical, or empirical copulas nelsen2007introduction, salvadori2007extremes, durante2016principles, benali2021mtcopula. Recent works such as zheng2021copula, zheng2022sensitivity have proposed sensitivity analysis with one of the most well studied and used elliptical copula, namely the Gaussian copula, but without the use of normalizing flows. In our current work, we present our analysis by assuming and approximating the unknown copula $\textcolor{blue}{\mathbf{\mathbb{C}}(U_A,U_Y)}$ with the Gaussian copula, while extending the normalizing flows with monotonic transformers that are recently shown to be universal non-linear SCM approximators Huang2018NAF, wehenkel2019UMNN, wehenkel2020GNF, balgi2022cgnf.
\looseness=-1 The Gaussian copula has been widely used in the fields of quantitative finance cherubini2004copula, salmon2009recipe, mackenzie2014formula, hydrology research renard2007use, zhang2019copulas, logistics kumar2019copula, astronomy takeuchi2010constructing, and similar fields nelsen2007introduction, salvadori2007extremes, durante2016principles. This is one of the several motivations for the particular selection of the Gaussian copula. The assumption and approximation of the unknown copula $\textcolor{blue}{\mathbf{\mathbb{C}}(U_A,U_Y)}$ in Eqs. (ref)-(ref) with a Gaussian copula $\textcolor{red}{\Phi_\rho(\Phi^{-1}(}\textcolor{blue}{U_A}\textcolor{red}{),\Phi^{-1}(}\textcolor{blue}{U_Y}\textcolor{red}{))}$ achieves causal effect identification as the back-door non-causal association between $Z_A$ and $Z_Y$ may be identified as \\ $\textcolor{red}{\mathbb{F}_{Z_A,Z_Y}(Z_A,Z_Y)}{=}\textcolor{blue}{\textcolor{blue}{\mathbf{\mathbb{C}}}}(\textcolor{red}{\Phi(Z_A),\Phi(Z_Y)}){\approx}\textcolor{red}{\Phi_\rho(Z_A,Z_Y)}$, where $\textcolor{red}{\rho}{\in}[-1,+1]$ is the Pearson's correlation between $Z_A$ and $Z_Y$. As the copula and the noises $\varepsilon_X$ and $U_X$ in Figure (ref) may exhibit non-linear dependence, it is more appropriate to equivalently represent the linear Pearson's correlation parameter $\textcolor{red}{\rho}$ in the Gaussian copula $\textcolor{red}{\Phi_\rho}$ in terms of the non-linear measure of association $\textcolor{blue}{\rho_{S_{\mathbf{\mathbb{C}}}}}$ (or $\textcolor{blue}{\tau_{K_{\mathbf{\mathbb{C}}}}}$) using the following results for bivariate Gaussian copula kruskal1958ordinal, meyer2013bivariate. This measure of back-door non-causal association $\textcolor{blue}{\rho_{S_{\mathbf{\mathbb{C}}}}}$ due to the Gaussian copula assumption equates to the Gaussian copula parameter $\textcolor{red}{\rho}$ as $\textcolor{red}{\rho}{=}2\sin(\pi\textcolor{blue}{\rho_{S_{\mathbf{\mathbb{C}}}}}/6)$ (i.e., $\textcolor{blue}{\rho_{S_{\mathbf{\mathbb{C}}}}}{\approx}\textcolor{red}{\rho}$). With the Gaussian copula $\textcolor{red}{\Phi_\rho}$ assumption in Eq. (ref)-(ref), the back-door non-causal association in Eqs. (ref)-(ref) is known, signifying that the causal effects are now identifiable under the assumed copula. Moreover, they can be estimated from a given observational dataset $\{(A^\ell,Y^\ell)\}^{N_{train}}_{\ell{ = }1}$ to train a parametric model as the proposed $\rho$-GNF in Figure (ref) by rewriting Eq. (ref) as
where $\mathbb{T}_A(\bullet;\theta_A)$ and $\mathbb{T}_{Y|A}(\bullet;\theta_Y)$ represent monotonic transformations parameterized by deep neural networks $\theta{=}(\theta_A,\theta_Y)$ using the integration-based unconstrained monotonic neural network (UMNN) transformer wehenkel2019UMNN. The conditioning of $Y$ on its parent variable $A$ in $\mathbb{T}_{Y|A}(\bullet;\theta_Y)$ is done with the graphical conditioner from graphical normalizing flow (GNF) wehenkel2020GNF using the directed acyclic graph (DAG) $A{\rightarrow}Y$ assumed in Figure (ref). We aptly refer to our model with the Gaussian copula assumption in Eq. (ref) as $\rho$-GNF due to the similarities to GNF and c-GNF that propose normalizing flows for observational density estimation and causal inference, but not for sensitivity analysis (herein lies our novelty over GNF or c-GNF). As any normalizing flow, the UMNN transformers and graphical conditioners are trained by maximizing the $\log$-likelihood of the observational training dataset $\{(A^\ell,Y^\ell)\}^{N}_{\ell{ = }1}$ wehenkel2020GNF, balgi2022cgnf, for a fixed $\rho$.
\looseness=-1 In principle, any copula maybe assumed in place of Gaussian copula in Eqs. (ref)-(ref). The Gaussian copula assumption, similar to normalizing flows tabak2010nf,tabak2013nf,rezende2015variationalNF, Papamakarios2017MAF, papamakarios2021NF_pmi, kobyzev2020NF, facilitates efficient computation of the $\log$-likelihood of the observational training dataset. Thus, enabling computationally efficient training of $\rho$-GNF in Eq. (ref). The Gaussian copula assumption further enables sampling $(Z_A, Z_Y)$ from $\textcolor{red}{\mathbb{F}_{Z_A,Z_Y}(Z_A,Z_Y)}{\approx}\textcolor{red}{\Phi_\rho(Z_A,Z_Y)}$ efficiently for the estimation of Monte-Carlo expectation in Eq. (ref), thus enabling computationally efficient inference. Most importantly, the Gaussian copula assumption provides a single bounded sensitivity parameter $\rho{\in}[-1,+1]$ for sensitivity analysis that can be used to control/model/block/adjust the back-door non-causal association under which the ACE is identifiable. Thus, enabling simple and efficient sensitivity analysis. Our subsequent experiments and results show that the Gaussian copula assumption works well empirically, as ilse2021efficient also observe for multiple unobserved confounders, the joint distribution of the observed variables becomes increasingly Gaussian due to the central limit theorem.
\looseness=-1 The main objective is to estimate the average causal effect (ACE), which can be expressed as
where $Y_a$ denotes the potential outcome under the intervention $A{:=}a$. In practise, the estimation of the ACE is done by Monte-Carlo expectation estimation by drawing the samples from the interventional distributions to approximate $\mathbf{E}[Y_1]$ and $\mathbf{E}[Y_0]$ as indicated below in Eqs. (ref)-(ref), after having trained the $\rho$-GNF for a specific measure of unobserved confounding $\rho$ on the given observational dataset. The First Law of Causal Inference pearl1999probabilitiesofcausation, pearl2009abductionactionprediction, pearl2009causality, pearl2018bookofwhy provides three steps, i.e, abduction, action and prediction, to estimate $\mathbf{E}[Y_a]$ and thus the $ACE_{\rho}$ in Eq. (ref).
\looseness=-1 vanderweele2017sensitivity propose the E-value as the minimum strength of association, on the risk ratio scale, that an unmeasured confounder would need to have with both the treatment and the outcome to fully explain away a specific causal treatment–outcome association, conditional on the measured covariates. A large E-value implies that considerable unmeasured confounding would be needed to explain away an effect estimate. A small E-value implies that little unmeasured confounding would be needed to explain away an effect estimate. Similar to the E-value, we propose the $\rho_{value}$ which represents the Gaussian copula parameter value that explains away the causal association between the observed treatment $A$ and the observed outcome $Y$. In other words, setting the Gaussian copula parameter $\rho{=}\rho_{value}$ results in $ACE_{\rho}{=}0$, i.e., $\mathbf{E}[Y_1]{=}\mathbf{E}[Y_0]$ in Eqs. (ref) and (ref), i.e., the potential outcomes are independent of the treatments/interventions. This implies that the strictly increasing transformation $\mathbb{T}^{-1}_{Y|A}$ modeling $Y$ in Eq. (ref) is independent of $A$, i.e., we have $Y{=}\mathbb{T}^{-1}_{Y}(Z_Y)$.
From the scale-invariance property of Spearman's correlation to strictly increasing transformations $\mathbb{T}^{-1}_{A}$ and $\mathbb{T}^{-1}_{Y}$, we have $\textcolor{red}{\rho_S(Z_A,Z_Y)}{=}\textcolor{blue}{\rho_{S_\mathbf{\mathbb{C}}}}{=}\textcolor{darkgreen}{\rho_S(A,Y)}{=}\textcolor{darkgreen}{\rho_{S_{Obs}}}$, where $\textcolor{darkgreen}{\rho_{S_{Obs}}}$ represents the observed Spearman's correlation between the observational data $\textcolor{darkgreen}{A}$ and $\textcolor{darkgreen}{Y}$. For a given observational dataset $A$ and $Y$, if the total observed association $\textcolor{darkgreen}{\rho_{S_{Obs}}}$ between $A$ and $Y$ is completely due to the non-causal association $\textcolor{blue}{\rho_{S_\mathbf{\mathbb{C}}}}$, there can be no contribution from the causal association, i.e., $ACE_{\rho}{=}0 \Leftrightarrow \rho{=}\textcolor{darkgreen}{\rho_{value}}$ where we have $\textcolor{darkgreen}{\rho_{value}}{=}2sin\left( \pi \textcolor{darkgreen}{\rho_{S_{Obs}}}/6\right)$. In Section (ref), we also experimentally verify that $\rho{=}\textcolor{darkgreen}{\rho_{value}} \Leftrightarrow ACE_{\rho}{=}0$. Hence, the E-value and $\rho_{value}$ represent the same concept, rather interpreted in different scales. Similar to the E-value that can be computed directly from the observational dataset mathur2018website, the $\rho_{value}$ is also identified as a single, bounded intuitive parameter that can be computed directly from Spearman's rho $\textcolor{darkgreen}{\rho_{S_{Obs}}}$ between $A$ and $Y$, i.e., $\rho_{value}{=}2\sin(\pi \textcolor{darkgreen}{\rho_{S_{Obs}}}/6)$. This further shows the similarities of the E-value and $\rho_{value}$. Since the total observed association between $A$ and $Y$ remains constant for a given observational dataset, the $ACE_{\rho}$ varies with $\rho \in [-1,+1]$. We aptly refer this sensitivity plot of $ACE_{\rho}$ to $\rho$ for a given observational dataset as the $\rho_{curve}$.
\looseness=-1 The $\rho_{value}$ equips an analyst/domain-expert to determine the sign of the ACE (i.e., whether the treatment is harmful or beneficial), which is arguably the most important part of the ACE and one of the ultimate goals of causal inference. Specifically, suppose the domain expert hypothesizes a measure of confounding in the interval $[\rho_{min},\rho_{max}]$. Thus, the $\rho_{curve}$ enables us to bound the true ACE to the narrower interval $[ACE_{\rho_{max}},ACE_{\rho_{min}}]$, which may in turn help us identify the most important insight of the causal inference, i.e., the sign of the true ACE. In particular, if $\rho_{min}{>}\rho_{value}$ (or $\rho_{max}{<}\rho_{value}$) then we may conclude that the true ACE is negative (or positive), as observed from the $\rho_{curve}$ in Figures (ref), (ref), and (ref). Furthermore, unlike E-value that only indicates the strength of the association due to unobserved confounding, $\rho_{value}$ indicates both strength and the sign of the association due to unobserved confounding, i.e., positive or negative association between $A$ and $Y$.
\looseness=-1 Our $\rho$-GNF is implemented in PyTorch paszke2017pytorch\footnote{The $\rho$-GNF code is available at \url{https://github.com/sobalgi/rhoGNF}.} by adapting the baseline code of UMNN wehenkel2019UMNN and GNF wehenkel2020GNF\footnote{GNF: \url{https://github.com/AWehenkel/Graphical-Normalizing-Flows}, UMNN: \url{https://github.com/AWehenkel/UMNN}.}. As normalizing flows are developed for continuous variables, we use the Gaussian dequantization trick from c-GNF balgi2022cgnf to model discrete variables into $\rho$-GNF. For the experiments, we present two different settings: (i) simulated dataset with continuous outcomes in Section (ref), and (ii) simulated dataset with binary outcomes in Section (ref). We additionally present the experiments using a real-world dataset with categorical outcomes analysing the impact of the IMF (International Monetary Fund) program on the degree of child poverty in the Global-South region by denoting the outcome as the total degree of child poverty ranging in 0-7 classes.
In practise, the $\rho_{curve}$ for a given observational dataset is obtained by varying the values of $\rho {\in}[-1,+1]$, i.e., we select $\rho{=}\{-0.99, -0.8, -0.6, -0.4, -0.2, 0.0, +0.2, +0.4, +0.6,\\ +0.8,$ $+0.99\}$ to train $\rho$-GNFs. We estimate the respective $ACE_{\rho}$ from Eqs. (ref)-(ref) to plot the $\rho_{curve}$ as interpolation of the points ($\rho, ACE_{\rho}$) shown in Figures (ref), (ref) and (ref). Our empirical ACE bounds are obtained as the infimum and supremum of the $\rho_{curve}$, i.e., $\inf[ACE_{\rho}] \leq ACE_{true} \leq \sup[ACE_{\rho}]$. We identify the $\rho_{value}$ that explains away the causal association for the benefit of domain expert/analyst.
In our first set of simulated experiments, we consider the SCM with continuous treatment $A$ and outcome $Y$ variables by hoover2006causality. This is a well-studied SCM from economics and econometrics, and is defined as below. The true causal, non-causal, and total observed (causal+non-causal) associations in the SCM are known and are used for verification purposes.
For a given set of SCM parameters $\alpha,\beta,\delta$ in Eq. (ref), we have the corresponding Pearson's correlation $\textcolor{darkgreen}{\rho_{P_{Obs}}}{=}\frac{\sigma_{A,Y}}{\sigma_{A}\sigma_{Y}}$, $\sigma^2_A{=}1$, $\sigma^2_Y{=}\alpha^2{+}\delta{+}2\alpha\beta$, $\sigma_{A,Y}{=}\alpha{+}\beta$, $\rho_{true}{=}\beta/\delta$ and $ACE_{true}{=}\alpha$ as tabulated in Table (ref). Figure (ref) indicates the scatter plot of each of the six observational datasets from the SCM in Eq. (ref) as parameterized in Table (ref). The SCMs that are observationally equivalent with the same $\textcolor{darkgreen}{\rho_{P_{Obs}}}$ present similar observational data distributions even though these data distributions have been generated from interventionally dissimilar SCMs. Figure (ref) shows observationally equivalent datasets result in equivalent $\rho_{curve}$, e.g., $\rho_{curve}$ of $\{SCM_1,SCM_2,SCM_3\}$ are similar since they are observationally equivalent as verifiable in the scatter plots of their respective observational dataset from Figure (ref). From Figure (ref) and Table (ref), we empirically observe the results $\rho{=}\textcolor{darkgreen}{\rho_{value}}$ corresponds to $ACE_{\rho}{=}0$ and $\rho{=}\rho_{true}$ corresponds to $ACE_{\rho}{=}ACE_{true}$ as presented in Section (ref). In summary, the main contribution of our method is the estimation of ACE as a function of the sensitivity parameter $\rho$ that indicates the degree of the confounding, thus distinguishing observationally equivalent SCMs based on the assumed degree of unconfoundedness.
The SCM in Eq. (ref) used to illustrate $\rho$-GNF with continuous outcomes may be considered a simple case, i.e. a linear SCM similar to the ones studied in cinelli2019sensitivity, cinelli2020making. However, as the $\rho$-GNF is parameterized to learn arbitrary non-linear monotonic functions using UMNN transformers wehenkel2019UMNN, our work is a generalization of the sensitivity analysis for simple linear SCMs to complex non-linear SCMs, as we demonstrate in the subsequent experiments with binary and categorical outcomes in Section (ref).
In the second set of experiments, we consider multiple (twenty) randomly generated observational datasets from randomly sampled data generating processes as proposed in sjolander2020note, sjolander2021novel, pena2022simple, where all variables are binary. As proposed, we randomly sample the parameters $\{\mathbb{P}(U), \mathbb{P}(A|U), \mathbb{P}(Y|A,U)\}$ of the binary data generating process (DGP) from the uniform distribution in [0,1], where $U$ denotes the confounder of $A$ and $Y$ which is intentionally hidden during training to simulate unobserved confounding. For this simple binary outcome/treatment/confounder case, the gold-standard AF bounds of a constant width of 1 $(=p_1{+}p_0)$ with 100% certainty of including the true ACE are presented in robins1989analysis, manski1990sensitivitybounds as
where $p_a{=}\mathbb{P}(A{=}a)$ and $q_a{=}\mathbb{P}(Y{=}1|A{=}a)$ are estimated using observational data. The $\rho_{curve}$ in Figure (ref) does include the true ACE in all of the randomly generated DGPs. Further, we obtain a narrower bound of 0.83$\pm$0.03 (i.e., $\approx$14% narrower than the AF bounds) such that our empirical bounds lie within the AF bounds, unlike the bounds by vanderweele2017sensitivity which may be wider than the AF bounds, as shown in sjolander2020note. The former comes as no surprise: Since our bounds are empirically obtained particularly assuming the Gaussian copula, they must lie within the AF bounds. Figure (ref) verifies it, i.e.,
Subsequently, we also extend the experiment from binary outcome to categorical outcome using a real-world dataset with categorical outcomes analysing the IMF (International Monetary Fund) program impact on the degree of child poverty in the Global-South region with the outcome as the total degree of child poverty ranging in 0-7 classes balgi2022counterfactual. Unlike the binary outcomes, the AF bounds are not available for the non-binary outcomes such as the categorical total degree of child poverty (degree 0: no poverty to degree 7: severe poverty). Since the degree of child poverty is formulated as the sum of seven binary individual dimension of child poverty, we may identify the AF bounds for each binary individual dimension of child poverty and extend the AF bounds to the total degree of child poverty. From Eq. (ref) and the observational dataset of the seven binary individual dimensions of child poverty, the assumption-free lower and upper bounds are identified as,
Since, the total degree of child poverty is the sum of these seven binary individual dimensions of child poverty, the lower and upper assumption-free bounds total degree of child poverty can be identified as sum of the respective lower and upper bounds. Thus, obtain the AF bounds for the total degree of child poverty as $[-3.3362,+3.5638]$ as indicated in the $\rho_{curve}$ in Figure (ref). Note that since the binary AF bounds are of width 1, and the total degree of child poverty is a sum of seven binary variables, the AF bounds for the total degree of the child poverty is of width 7. From Figure (ref), our ACE bounds $[-1.73, +2.11]$ are verified empirically with a narrower width of 3.84, i.e., a reduction of the bounds by $45.1\%$ under the Gaussian copula assumption.
\looseness=-1 We proposed a novel copula-based approach for sensitivity analysis termed $\rho$-GNF to bound the causal effect in the risk difference scale where $\rho{\in}[-1,+1]$ is a bounded sensitivity parameter representing the unobserved back-door non-causal association between the observed treatment and outcome. Under the Gaussian copula assumption, we showed that $\rho$-GNF enabled us to estimate the causal effect as a function of $\rho$ in the form of $\rho_{curve}$. The $\rho_{curve}$ enabled us to identify narrower empirical bounds in contrast to the wider AF bounds, irrespective of discrete or continuous outcome variables. We identified $\rho_{value}$ as the measure of the unobserved confounding that explains away the causal effect and presented the similarities to the E-value with experimental validation. Further, the $\rho_{curve}$ enabled us to provide finer bounds for the causal effect given an interval of $\rho$ values that the domain expert considers appropriate to identify the sign of the ACE, thus deducing if the treatment is beneficial or harmful. The adaptability of $\rho$-GNF for both discrete and continuous outcomes should encourage the use of sensitivity analysis when working with non-randomized observational data to draw causal conclusions.