(Conditional) Independence

Simulation notebook

Goals:

• Illustrate (conditional) independence in general and for causal inference with binary treatments in particular

Acknowledgements: I thank Henri Pfleiderer for his assistance in preparing this notebook.

Formal definition of (conditional) independence

Discrete Variables

Independence of two discrete variables \(Y\) and \(W\):

\[ Y \perp\!\!\!\!\perp W \iff P(Y=y, W=w) = P(Y = y) \cdot P(W = w)\]

This implies by Bayes’ law:

\[ P(Y=y \mid W=w) = \frac{P(Y=y, W=w)}{P(W=w)} = \frac{P(Y = y) \cdot P(W = w)}{P(W=w)} = P(Y = y) \]

Continuous Variables

Equivalently, two continuous variables \(Y\) and \(W\) being independent means:

\[ Y \perp\!\!\!\!\perp W \iff f_{Y,W}(y,w) = f_Y(y) \cdot f_W(w) \] where \(f\) denotes probability density functions.

Mixed Case

In addition to these two cases there is a mixed case. Consider for example a continuous \(Y\) and a discrete \(W\). This is a very common setup in causal inference. We may think of \(Y\) as a continuous outcome and \(W\) a binary treatment indicator. In this case:

\[ Y \perp\!\!\!\!\perp W \iff f_{Y,W}(y,w) = f_Y(y) \cdot P(W=w) \]

Using the definition of the conditional density or probability respectively, this implies:

\[ f_{Y \mid W}(y \mid w) = \frac{f_{Y,W}(y,w)}{P(W=w)} = \frac{f_Y(y) \cdot P(W=w)}{P(W=w)} =f_Y(y) \] and

\[ P(W = w \mid Y = y) = \frac{f_{Y,W}(y,w)}{f_Y(y)} = \frac{f_Y(y) \cdot P(W=w)}{f_Y(y)} = P(W =w) \]

Conditional Independence

For compactness, we focus on the mixed case here. We have three random variables, a continuous \(Y\) and binary \(X\) and \(W\). We say that \(Y\) and \(W\) are independent conditional on \(X\) iff the following holds:

\[ Y \perp\!\!\!\!\perp W \mid X \iff f_{Y,W \mid X}(y, w \mid x) = f_{Y \mid X}(y \mid x) \cdot P(W = w \mid X = x) \]

This also implies for all \(w\) with \(P(W = w \mid X = x)\) that

\[f_{Y\mid W, X}(y \mid w,x) = \frac{f_{Y,W \mid X}(y,w \mid x)}{P(W = w \mid X = x)} = \frac{ f_{Y \mid X}(y \mid x) \cdot P(W = w \mid X = x)}{P(W = w \mid X = x)} = f_{Y \mid X} (y \mid x)\] i.e., the conditioning on \(W\) doesn’t matter once we condition on \(X\). As the (unconditional) independence above, this works similarly for the cases of only continuous and only discrete variables using only densities or probabilities, respectively.

Conditional independence implies conditional mean independence

Conditional independence of \(Y\) and \(W\) given \(X\) implies mean independence, meaning: \(E[Y \mid X,W] = E[Y \mid X]\):

\[ E[Y \mid W,X] = \int_{-\infty}^{+\infty} y \underbrace{f_{Y \mid W,X}(y\mid w,x)}_{=f_{Y \mid X} (y \mid x), \ cond. \ indep.} dy = \int_{-\infty}^{+\infty} y \ \ f_{Y \mid X} (y \mid x) \ dy = E[Y \mid X] \]

All of this is illustrated in the following examples.

General Example

Let’s consider the following DGP:

\[\begin{align} Y= X + U \ \ \ with \ \ \ U \sim \mathcal{N}(0,0.5) \\ X = \mathbf{1}(V\geq 0) \\ W = \mathbf{1}(Z\geq 0) \\ [V, Z]'\sim \mathcal{N}\left([0~~0]',\left[ \begin{array} {rrr} 1 & \rho \\ \rho & 1 \end{array}\right]\right) \end{align}\]

First, consider the case where \(Y\) and \(W\) are independent (i.e. \(\rho =0\)). This is illustrated by drawing a large sample from the described DGP with \(\rho = 0\):

if (!require("tidyverse")) install.packages("tidyverse", dependencies = TRUE); library(tidyverse)
if (!require("ggridges")) install.packages("ggridges", dependencies = TRUE); library(ggridges)
if (!require("MASS")) install.packages("MASS", dependencies = TRUE); library(MASS)

set.seed(1234)
n = 1000000
mu = c(0,0)
rho = 0
sigma = matrix(c(1,rho,rho,1), nrow = 2)
draw = mvrnorm(n, mu, sigma)
X = 1*(draw[,1]>0)
W = 1*(draw[,2] > 0)
Y = X + rnorm(n, sd = 0.5)

Let’s plot the distributions of \(Y\) for \(W=0\) and \(W=1\):

tibble(Y,W,X) %>% mutate(W = as.factor(W)) %>% 
  ggplot(aes(x = Y, y = fct_rev(W), fill = W))+
  geom_density_ridges(alpha = 0.6)+
  xlab("Y") + ylab("Density")

The distribution of \(Y\) is the same when conditioning on different values of \(W\). Consequently, \(Y\) and \(W\) are independent. This is a feature of the (conditional) distribution of \(Y\) and thus also holds for the (conditional) expected value:

\[ E[Y \mid W=1] = E[Y \mid W=0] = E[Y] = 0.5 \]

mean(Y)

[1] 0.5002274

mean(Y[W==1])

[1] 0.49979

mean(Y[W==0])

[1] 0.5006651

Now, let’s introduce some correlation between the variables \(V\) and \(Z\) by setting \(\rho = 0.6\). This creates dependence between \(Y\) and \(W\):

set.seed(1234)
n = 1000000
mu = c(0,0)
rho = 0.6
sigma = matrix(c(1,rho,rho,1), nrow = 2)
draw = mvrnorm(n, mu, sigma)
X =  1*(draw[,1]>0) 
W = 1*(draw[,2] > 0)
Y = X + rnorm(n, sd = 0.5)

tibble(Y,W,X) %>% mutate(W = as.factor(W)) %>% 
  ggplot(aes(x = Y, y = fct_rev(W), fill = W))+
  geom_density_ridges(alpha = 0.6)+
  xlab("Y")

The distributions are not the same anymore, \(Y\) is not independent of \(W\).

It is however, once we condition on \(X\). In the present context, this can be seen by looking at the distribution of \(Y\) in the 4 different subgroups formed by the 2 binary variables \(W\) and \(X\).

tibble(Y,W,X) %>%  mutate(W = as.factor(W), X = as.factor(X), Subgroup = interaction(W,X)) %>%
  ggplot(aes(x = Y, y = fct_rev(Subgroup), fill = Subgroup))+
  geom_density_ridges(alpha = 0.6)+
  scale_fill_discrete(labels = c("W=0, X=0","W=1, X=0","W=0, X=1","W=1, X=1"))+
  xlab("Y") + ylab("Density")

Once we condition on \(X\), i.e. look at the distribution of \(Y\) in the 2 subgroups with \(X=0\) and \(X=1\), the distributions don’t change when conditioning on different values of \(W\). \(Y\) and \(W\) are independent conditional on \(X\): \(Y \perp\!\!\!\!\perp W \mid X\)

In this case, \(E[Y \mid W=1] \neq E[Y \mid W=0] \neq E[Y]\):

mean(Y)

[1] 0.4995484

mean(Y[W==1])

[1] 0.7044881

mean(Y[W==0])

[1] 0.2947923

Conditional on X, we can see that: \(E[Y \mid W= 0, X=x] = E[Y \mid W= 1, X=x] = E[Y \mid X=x]\) for \(x \in \{0,1\}\):

tibble(Y,W,X) %>% group_by(X,W) %>% summarise(mean_Y = mean(Y))

tibble(Y,W,X) %>% group_by(X) %>% summarise(mean_Y = mean(Y))

Example: Potential Outcomes and causal inference

Let’s look at the importance of the above concept in a causal inference context. \(Y\) denotes the outcome, \(W\) the treatment and \(X\) some confounding variable. \(Y(1)\) and \(Y(0)\) are the potential outcomes in the treated and untreated case respectively, so \(Y = W \cdot Y(1) + (1-W) \cdot Y(0)\).

set.seed(1234)
n = 1000000
mu = c(0,0)
rho = 0
sigma = matrix(c(1,rho,rho,1), nrow = 2)
draw = mvrnorm(n, mu, sigma)
X =  1*(draw[,1] > 0)
W = 1*(draw[,2] > 0)
Y0 = X + rnorm(n, sd = 0.5)
Y1 = 0.5 + X + rnorm(n, sd = 0.5)
Y = W*Y1 + (1-W)*Y0

Consider the following DGP:

• \(Y(0) = X + U_0\) with \(U_0 \sim \mathcal{N}(0,0.5)\).

• \(Y(1) = 0.5 + X + U_1\) with \(U_1 \sim \mathcal{N}(0,0.5)\).

and \[\begin{align} X = \mathbf{1}(V\geq 0) \\ W = \mathbf{1}(Z\geq 0) \\ [V, Z]'\sim \mathcal{N}\left([0~~0]',\left[ \begin{array} {rrr} 1 & \rho \\ \rho & 1 \end{array}\right]\right) \end{align}\]

First, consider again the case with \(\rho = 0\) to illustrate independence in the sense that \(Y(w) \perp\!\!\!\!\perp W\).

To see that this holds in the first example, plot the distributions of \(Y(1)\) for the treated and the untreated separately. As we have independence, they should look identical:

tibble(Y1,Y0,Y,W,X) %>% mutate(W = as.factor(W)) %>% 
  ggplot(aes(x = Y1, y = fct_rev(W), fill = W))+
  geom_density_ridges(alpha = 0.6)+
  xlab("Y(1)")+ ylab("Density")

Same for \(Y(0)\):

tibble(Y1,Y0,Y,W,X) %>% mutate(W = as.factor(W)) %>% 
  ggplot(aes(x = Y0, y = fct_rev(W), fill = W))+
  geom_density_ridges(alpha=.6)+
  xlab("Y(0)")+ ylab("Density")

These graphs illustrate that the potential outcomes are independent of the treatment indicator.

The individual treatment effects are given as: \(Y(1) - Y(0) = 0.5 + X + U_1 - X - U_0 = 0.5 + U_1 - U_0\). The ATE is \(E[Y(1)-Y(0)] = 0.5\).

As the potential outcomes are independent of \(W\), we can estimate the ATE using a simple mean comparison between the treated and the untreated group:

TE_mean_comparison = mean(Y[W==1]) - mean(Y[W==0])
print(TE_mean_comparison)

[1] 0.4984061

As above, let’s introduce some correlation between \(V\) and \(Z\) by setting $= 0.6. Now, \(X\) and \(W\) are not independent:

rho = 0.6
sigma = matrix(c(1,rho,rho,1), nrow = 2)
draw = mvrnorm(n, mu, sigma)
X =  1*(draw[,1] > 0)
W = 1*(draw[,2] > 0)
Y0 = X + rnorm(n, sd = 0.5)
Y1 = 0.5 + X + rnorm(n, sd = 0.5)
Y = W*Y1 + (1-W)*Y0

Intuitively, observations with a high \(X\) (and thus a high \(Y\)) are now more likely to be treated, as \(V\) and \(Z\), the two variables determining \(X\) and \(W\), are positively correlated. Comparing means between treated and untreated will give an incorrect, too large, estimate of the ATE. In a first step, let’s again illustrate the lack of independence between the potential outcomes and \(W\):

tibble(Y1,Y0,Y,W,X) %>% mutate(W = as.factor(W)) %>% 
  ggplot(aes(x = Y1, y = fct_rev(W), fill = W))+
  geom_density_ridges(alpha=.6)+
  xlab("Y(1)")+ylab("Density")

Same for \(Y(0)\):

tibble(Y1,Y0,Y,W,X) %>% mutate(W = as.factor(W)) %>% 
  ggplot(aes(x = Y0, y = fct_rev(W), fill = W))+
  geom_density_ridges(alpha=.6)+
  xlab("Y(0)")+ylab("Density")

The potential outcomes are not independent of \(W\). Estimate the ATE by mean comparison:

TE_mean_comparison = mean(Y[W==1]) - mean(Y[W==0])
print(TE_mean_comparison)

[1] 0.9111342

Let’s condition on \(X\) and see whether they are independent. For this purpose, look at the distributions of \(Y(w)\) for the different subgroups formed by \(W\) and \(X\):

tibble(Y,Y1,Y0,W,X) %>%  mutate(W = as.factor(W), X = as.factor(X), Subgroup = interaction(W,X)) %>% 
  ggplot(aes(x = Y1, y = fct_rev(Subgroup), fill = Subgroup))+
  geom_density_ridges(alpha = 0.6)+
  scale_fill_discrete(labels = c("W=0, X=0","W=1, X=0","W=0, X=1","W=1, X=1"))+
  xlab("Y(1)")+ylab("Density")

tibble(Y,Y1,Y0,W,X) %>%  mutate(W = as.factor(W), X = as.factor(X), Subgroup = interaction(W,X)) %>% 
  ggplot(aes(x = Y0, y = fct_rev(Subgroup), fill = Subgroup))+
  geom_density_ridges(alpha = 0.7)+
  scale_fill_discrete(labels = c("W=0, X=0","W=1, X=0","W=0, X=1","W=1, X=1"))+
  xlab("Y(0)")+ylab("Density")

When we look at conditional distributions of the potential outcomes, they are the same for treated and untreated. This shows that conditional on \(X\), the potential outcomes are independent of \(W\). What does this imply for the estimation of the ATE?

We can estimate the ATE by comparing conditional means:

TE_cond_mean_comparison_1 = mean(Y[W==1 & X==1]) - mean(Y[W==0 & X==1])
print(TE_cond_mean_comparison_1)

[1] 0.5019498

TE_cond_mean_comparison_0 = mean(Y[W==1 & X==0]) - mean(Y[W==0 & X==0])
print(TE_cond_mean_comparison_0)

[1] 0.4993736

TE_cond_mean_comparison = mean(X)*TE_cond_mean_comparison_1 + mean(1-X)*TE_cond_mean_comparison_0
print(TE_cond_mean_comparison)

[1] 0.5006614

Alternatively, use a regression model with \(X\) as a control variable:

summary(lm(Y~ W + X))

Call:
lm(formula = Y ~ W + X)

Residuals:
     Min       1Q   Median       3Q      Max 
-2.40876 -0.33723 -0.00068  0.33697  2.46799 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.001105   0.000778  -1.421    0.155    
W            0.500661   0.001097 456.274   <2e-16 ***
X            1.000412   0.001097 911.718   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.5003 on 999997 degrees of freedom
Multiple R-squared:  0.6241,    Adjusted R-squared:  0.6241 
F-statistic: 8.301e+05 on 2 and 999997 DF,  p-value: < 2.2e-16

This recovers the true ATE. Note, that here the true structural model is actually linear and OLS can be used to recover the true coefficients on \(X\) and \(W\):

\[ Y = W\cdot (0.5 + X + U_{1}) + (1-W) \cdot (X + U_{0}) = 0.5W + X +\underbrace{WU_{1} + (1-W)U_{0}}_{=U} = 0.5W + X + U\]