Causal Inference VIII — DAGs, Mediation & the SCM (R companion)

dagitty reads identification off the graph; ggdag draws it; mediation decomposes the effect

R has the reference tooling for the graphical approach. Where the Python notebook coded each adjustment by hand, dagitty (Textor et al.) derives the answer directly from the graph — given a DAG, it returns the valid adjustment sets, finds instruments, and lists the testable conditional independencies. This is the practical payoff of the SCM framework: you draw your assumptions, and the software tells you what is identified and what to condition on.

• dagitty — adjustmentSets() for the back-door criterion, latent-variable handling for the front-door case;
• ggdag — publication-quality DAG drawing;
• mediation (Imai, Keele, Tingley, Yamamoto) — mediate() for the ACME/ADE decomposition with proper (quasi-Bayesian / bootstrap) confidence intervals.

Results match the from-scratch Python notebook: adjust for the confounder, never for the collider, identify through the front door, and split the effect into direct and indirect parts.

1. The back-door criterion — dagitty::adjustmentSets

Encode the confounding DAG ($Z\to X$, $Z\to Y$, $X\to Y$) and ask dagitty what to adjust for. It returns { Z } — the back-door criterion, computed from the graph. ggdag draws it, and the regression confirms that adjusting for $Z$ recovers the true effect of 2 while the naive regression is biased.

suppressMessages({library(dagitty); library(ggdag); library(ggplot2)})
g1<-dagitty('dag{ Z->X Z->Y X->Y }')
cat("Back-door adjustment set for X -> Y:\n"); print(adjustmentSets(g1, exposure="X", outcome="Y"))
set.seed(0); n<-6000; Z<-rnorm(n); X<-Z+rnorm(n); Y<-2*X+3*Z+rnorm(n)
cat(sprintf("\nnaive   Y~X   = %.2f  (biased)\n", coef(lm(Y~X))["X"]))
cat(sprintf("back-door Y~X+Z = %.2f  (true 2)\n", coef(lm(Y~X+Z))["X"]))
options(repr.plot.width=6, repr.plot.height=4.2)
coordinates(g1)<-list(x=c(X=0,Y=2,Z=1), y=c(X=0,Y=0,Z=1))
print(ggdag(g1) + theme_dag() + ggtitle("Confounder Z: adjust to close the back door"))

Warning message:
"package 'dagitty' was built under R version 4.6.1"

Warning message:
"package 'ggdag' was built under R version 4.6.1"

Back-door adjustment set for X -> Y:

{ Z }

naive   Y~X   = 3.50  (biased)

back-door Y~X+Z = 2.01  (true 2)

2. Collider bias — dagitty says adjust for nothing

For the collider DAG ($X\to Y$, $X\to C$, $Y\to C$), dagitty returns the empty set { } — you should condition on nothing, because $C$ is a collider that must be left alone. Conditioning on it opens a non-causal path and biases the estimate, which the regression confirms: adding $C$ drags the estimate of a true effect of 2 far off. dagitty encodes the rule that catches this automatically.

g2<-dagitty('dag{ X->Y X->C Y->C }')
cat("Adjustment set for X -> Y (collider present):\n"); print(adjustmentSets(g2, exposure="X", outcome="Y"))
set.seed(1); X<-rnorm(n); Y<-2*X+rnorm(n); C<-X+Y+rnorm(n)
cat(sprintf("\ncorrect  Y~X    = %.2f  (true 2)\n", coef(lm(Y~X))["X"]))
cat(sprintf("WRONG    Y~X+C  = %.2f  (conditioning on the collider CREATES bias)\n", coef(lm(Y~X+C))["X"]))
options(repr.plot.width=6, repr.plot.height=4.2)
coordinates(g2)<-list(x=c(X=0,Y=2,C=1), y=c(X=1,Y=1,C=0))
print(ggdag(g2) + theme_dag() + ggtitle("Collider C: adjustment set is EMPTY — do not condition on C"))

Adjustment set for X -> Y (collider present):

 {}

correct  Y~X    = 1.99  (true 2)

WRONG    Y~X+C  = 0.49  (conditioning on the collider CREATES bias)

3. The front-door criterion — when no back-door set exists

Mark the confounder $U$ as latent and ask dagitty for the back-door adjustment set: it returns nothing — with $U$ unobserved, no adjustment identifies the effect through the back door. But the mediator $M$ opens the front door. We verify there is no back-door solution, then compute the front-door estimate by hand — the $X\to M$ effect composed with the $M\to Y$ effect (adjusting for $X$) — and recover the true effect of 1.2 that the naive regression missed.

g3<-dagitty('dag{ U[latent] U->X U->Y X->M M->Y }')
bd<-adjustmentSets(g3, exposure="X", outcome="Y")
cat("Back-door adjustment sets with U latent:", if(length(bd)==0) "NONE (back door is unidentifiable)\n" else "\n")
set.seed(2); U<-rnorm(n); X<-U+rnorm(n); M<-0.8*X+rnorm(n); Y<-1.5*M+2*U+rnorm(n)
a<-coef(lm(M~X))["X"]; b<-coef(lm(Y~M+X))["M"]
cat(sprintf("naive     Y~X            = %.2f  (biased by unobserved U)\n", coef(lm(Y~X))["X"]))
cat(sprintf("front-door (X->M)x(M->Y|X) = %.2f x %.2f = %.2f  (true 1.2)\n", a, b, a*b))
options(repr.plot.width=6.5, repr.plot.height=4.2)
coordinates(g3)<-list(x=c(X=0,M=1,Y=2,U=1), y=c(X=0,M=0,Y=0,U=1))
print(ggdag(g3) + theme_dag() + ggtitle("Front door: U latent (no back door), identify through M"))

Back-door adjustment sets with U latent: NONE (back door is unidentifiable)

naive     Y~X            = 2.19  (biased by unobserved U)

front-door (X->M)x(M->Y|X) = 0.80 x 1.51 = 1.22  (true 1.2)

4. Causal mediation with inference — mediation::mediate

The mediation package is the reference implementation of the Imai-Keele-Tingley-Yamamoto framework. Given a mediator model ($M\sim X$) and an outcome model ($Y\sim X+M$), mediate() returns the ACME (indirect effect), the ADE (direct effect), and the total effect — each with confidence intervals from a quasi-Bayesian or bootstrap simulation, the inference the from-scratch product-of-coefficients cannot easily provide. It recovers the simulated decomposition (true direct 0.5, indirect 0.84).

suppressMessages(library(mediation))
set.seed(3); nn<-3000; Xt<-rbinom(nn,1,.5); Mm<-0.7*Xt+rnorm(nn); Yy<-0.5*Xt+1.2*Mm+rnorm(nn)
df<-data.frame(Xt,Mm,Yy)
med<-mediate(lm(Mm~Xt,df), lm(Yy~Xt+Mm,df), treat="Xt", mediator="Mm", sims=500)
cat(sprintf("ACME  (indirect) = %.3f   95%% CI [%.3f, %.3f]   (true 0.84)\n", med$d0, med$d0.ci[1], med$d0.ci[2]))
cat(sprintf("ADE   (direct)   = %.3f   95%% CI [%.3f, %.3f]   (true 0.50)\n", med$z0, med$z0.ci[1], med$z0.ci[2]))
cat(sprintf("Total effect     = %.3f\n", med$tau.coef))
cat(sprintf("Prop. mediated   = %.2f\n", med$n0))
options(repr.plot.width=7, repr.plot.height=4)
plot(med, main="mediation::mediate — ACME, ADE, and total effect with CIs")

Warning message:
"package 'mediation' was built under R version 4.6.1"

ACME  (indirect) = 0.896   95% CI [0.803, 0.989]   (true 0.84)

ADE   (direct)   = 0.452   95% CI [0.376, 0.529]   (true 0.50)

Total effect     = 1.349

Prop. mediated   = 0.66

5. Summary

dagitty turned each identification question into a one-line query: { Z } for the confounder (adjust), { } for the collider (adjust for nothing — the warning against bad controls), and no set when a latent confounder blocks the back door (forcing the front-door route). ggdag drew the graphs and mediation::mediate decomposed the effect into direct and indirect parts with valid confidence intervals — all matching the from-scratch Python results.

This is the practical case for the graphical framework: draw your causal assumptions, and the identification strategy — which variables to adjust for, whether the effect is identified at all — follows mechanically, avoiding the pervasive errors of conditioning on colliders and mediators. Combined with the potential-outcomes language of the rest of the arc (the back-door criterion is unconfoundedness; the do-operator is the potential outcome), it completes the toolkit for reasoning about identification before touching an estimator. Next: Heterogeneous Effects & Double/Debiased Machine Learning, where these identified effects are allowed to vary across units and are estimated with machine learning.