Causal Inference Causal Discovery from Association to Causation

# Causal Inference & Causal Discovery: From Association to Causation

## 1. Introduction & Motivation

Causal inference represents a fundamental shift from traditional machine learning's focus on prediction to understanding cause-and-effect relationships. While classical ML excels at identifying correlations (association), real-world decision-making requires understanding causality: What happens if we intervene? How would policy changes affect outcomes?

The distinction is critical: two variables may be highly correlated without one causing the other. Age and shoe size correlate in children due to common cause (growth), not causal relationship between them. Traditional supervised learning memorizes such associations but fails under intervention or distribution shift.

Causal inference enables:
1. Causal Effect Estimation: Quantifying impact of interventions from observational data
2. Counterfactual Reasoning: "What if" scenarios for policy evaluation
3. Robustness to Distribution Shift: Models generalizing across environments
4. Feature Selection: Identifying truly relevant variables
5. Algorithm Fairness: Understanding discriminatory mechanisms

The field synthesizes statistics (causal models, potential outcomes), philosophy (causality definitions), and ML (algorithms for discovery and estimation).

## 2. Core Concepts & Theory

### Causal Graphs and DAGs

Causal relationships represented as directed acyclic graphs (DAGs) where edges represent direct causal effects. Nodes are variables, edges indicate causation direction.

Example: Health outcome affected by treatment and age:

Age → Outcome
Treatment → Outcome

DAGs encode assumptions about causal structure. Two graphs differing in edge direction may imply different causal relationships.

### Confounding and Confounders

Confounder: variable affecting both treatment and outcome. Creates spurious association.

Example: Income confounds education-health relationship (both affected by family wealth):

Family Wealth → Education
Family Wealth → Health
Education → Health (true causal effect)

Ignoring confounding leads to biased causal estimates.

### D-Separation and Conditional Independence

d-separation: Criterion determining whether variable sets are independent in DAG.

Paths between variables: blocked by:
1. Collider node (V-shape) if not conditioned on
2. Any other node if conditioned on

d-separation implies conditional independence: if A ⊥ B | C in DAG, then P(A|B,C) = P(A|C).

### Potential Outcomes Framework

Neyman-Rubin causal model: each unit has potential outcomes under different treatment levels.

For binary treatment:
- Y_i(1): outcome if unit i receives treatment
- Y_i(0): outcome if control

Causal effect: tau_i = Y_i(1) - Y_i(0) (unobservable for each unit—fundamental problem of causal inference)

Observed outcome: Y_i = T_i * Y_i(1) + (1 - T_i) * Y_i(0) where T_i in {0, 1} is treatment indicator.

### Identifiability Conditions

Causal effects identifiable when:

1. Unconfoundedness (CIA): Y(t) _|_ T | X for all treatment levels t
- No unmeasured confounders
- All common causes observed

2. Overlap (Positivity): 0 < P(T=t|X) < 1 for all treatment levels t and X values
- All units have positive probability receiving each treatment level
- No deterministic treatment assignment

3. Consistency: Y = Y(T) (no hidden versions of treatment)

Under these assumptions, causal effects identifiable from observational data.

## 3. Mathematical Formulation

### Average Treatment Effect (ATE)

Most common causal estimand:
$$ au = E[Y(1) - Y(0)] = E[Y(1)] - E[Y(0)]$$

Under unconfoundedness and overlap:
$$ au = E_X[E[Y|T=1,X] - E[Y|T=0,X]]$$

Nonparametric estimator: average treatment effect estimate over observed confounders.

### Heterogeneous Treatment Effects (HTE)

Treatment effects vary across units:
$$ au(x) = E[Y(1) - Y(0)|X=x]$$

Conditional average treatment effect (CATE) at covariate values x.

### Inverse Probability Weighting (IPW)

Estimate causal effects by reweighting observations:
$$\hat{ au}_{ ext{IPW}} = \frac{1}{n}\sum_i \left(\frac{T_i}{\hat{p}(X_i)}Y_i - \frac{1-T_i}{1-\hat{p}(X_i)}Y_i ight)$$

where p_hat(x) = P(T=1|X=x) is propensity score. Reweighting creates pseudo-population where treatment independent of confounders.

### Doubly Robust Estimation

Combine regression and IPW for robustness:
$$\hat{ au}_{ ext{DR}} = \frac{1}{n}\sum_i \left[\frac{T_i}{\hat{p}(X_i)}(Y_i - \hat{m}_1(X_i)) + \hat{m}_1(X_i) - \frac{1-T_i}{1-\hat{p}(X_i)}(Y_i - \hat{m}_0(X_i)) - \hat{m}_0(X_i) ight]$$

where m_hat_t(x) estimates E[Y|T=t,X=x]. Consistent if either propensity score or outcome regression correct.

### Regression Adjustment

Estimate using conditional expectation:
$$\hat{ au} = \frac{1}{n}\sum_i[\hat{E}[Y|T=1,X_i] - \hat{E}[Y|T=0,X_i]]$$

Simple but biased if model misspecified; doubly robust methods address this.

## 4. Advanced Theory & Extensions

### Instrumental Variables (IV)

When unconfoundedness violated, use instrumental variable: variable Z affecting outcome only through treatment.

Requirements:
1. Relevance: Z correlated with T
2. Exogeneity: Z independent of unobserved confounders (U)
3. Exclusion: Z affects Y only through T

Graphically: Z → T → Y with U → T and U → Y but Z independent of U.

Two-stage least squares (2SLS) standard IV estimator:
$$T_i = \alpha + \gamma Z_i + u_i$$
$$Y_i = \beta_0 + \beta_1 \hat{T}_i + \epsilon_i$$

First stage predicts treatment using IV; second stage regresses outcome on predicted treatment.

### Sensitivity Analysis

Quantify robustness to unmeasured confounding. Rotnitzky-Robins sensitivity parameter Lambda bounds confounder strength.

If maximum confounder association bounded by Lambda, causal effect estimate still valid if inference accounts for sensitivity bounds.

### Front-Door Criterion

When IV unavailable, front-door criterion enables identification if causal path goes through intermediate mediator with complete information.

Requires: Z independent of U, all paths from T to U go through M (mediator).

## 5. Computational Considerations

### Propensity Score Estimation

Critical for IPW and matching methods. Models:
- Logistic Regression: Fast, interpretable, assumes linearity
- Random Forests: Flexible, captures interactions
- Gradient Boosting: Often best balance of flexibility and stability
- Neural Networks: Most flexible but risk overfitting with small sample sizes

Overlap assessment critical: plot propensity score distributions for treated/control; should overlap substantially. Regions without overlap (common support violation) lead to unreliable estimates.

### Scalability Issues

Standard causal methods scale to millions of observations but thousands of confounders problematic. Solutions:

1. Variable Selection: Use high-dimensional variable selection (Lasso, elastic net) to reduce confounder set
2. Representation Learning: Learn low-dimensional confounder representation via autoencoders
3. Balancing Weights: Iteratively reweight to balance confounder distributions

### Causal Discovery Complexity

Discovering DAG structure from data computationally hard (NP-hard in general). Practical algorithms:

• Constraint-based (PC algorithm): O(p^d) complexity where p is variables, d is graph degree
• Score-based (structure learning): Greedy search through exponential space
• Functional models (LiNGAM): O(p³) if linear non-Gaussian relationships

Discovery requires larger sample sizes than effect estimation (rule of thumb: n >> p²).

## 6. Practical Implementation Strategies

### Causal Inference Workflow

1. Problem Formulation: Define treatment, outcome, and confounders
2. DAG Construction: Draw causal graph encoding assumptions
3. Identifiability Check: Verify conditions met for effect estimation
4. Confounder Adjustment: Select appropriate method (matching, weighting, regression, IV)
5. Robustness: Sensitivity analysis, alternative specifications
6. Interpretation: Communicate uncertainty and assumptions

### Software Frameworks

DoWhy (Microsoft): High-level API for causal inference, supports multiple methods.

EconML (Microsoft): Focus on heterogeneous treatment effects, includes multiple learners.

CausalML: Open-source implementation of state-of-the-art HTE algorithms.

Causalimpact: Time-series causal inference for intervention analysis.

R causal packages: CausalForesight, DoWhy integration, grf (generalized random forests).

### Method Selection Guide

### Propensity Score Matching

Algorithm:
1. Fit logistic regression: P(T=1|X)
2. For each treated unit, find control unit with closest propensity score (caliper matching)
3. Estimate ATE on matched sample

Reduces bias from overt bias but may amplify hidden bias if unconfoundedness false.

## 7. Benchmark Datasets & Evaluation

### Benchmark Datasets

IHDP (Infant Health and Development Program): 747 samples, real medical intervention study. True causal effect known. Standard for HTE evaluation.

ACIC (Atlantic Causal Inference Conference) Challenge: Synthetic data with known ground truth causal DAGs of varying complexity. Multiple data-generating mechanisms.

Jobs Program Data: Real quasi-experiment (Job Training Partnership Act) with experimental gold standard available. Tests real-world causal methods.

NSW Labar Dataset: Econometric dataset comparing treatment (job training) effect estimates across methods.

### Evaluation Metrics

Bias: E[tau_hat - tau*] where tau* is true effect. Root Mean Squared Error: sqrt(E[(tau_hat - tau*)^2]).

Coverage: Do 95% confidence intervals contain true effect in 95% of simulations? Critical for valid inference.

Heterogeneous Effect Recovery: How well estimated tau(x) correlates with true tau*(x)?

Robustness to Violation: How estimates degrade under unmeasured confounding?

### Comparison Baselines

• Naive: Compare mean outcomes treated vs. control (maximum bias)
• Linear OLS: Simple regression adjustment
• Propensity Score: IPW or matching
• BART: Flexible nonparametric regression
• Causal Forests: State-of-art for HTE
• Experimental (Gold Standard): Randomized trial results

## 8. Key Challenges & Limitations

### Unconfoundedness is Unverifiable

Core assumption (no unmeasured confounding) untestable from data alone. Sensitivity analyses provide robustness bounds but cannot rule out arbitrarily strong unmeasured confounding.

Implication: Causal claims from observational data always rest on subjective judgment about completeness of confounder measurement.

### Sample Size Requirements

Discovering causal structure requires much larger n than estimating effects. Rule of thumb:
- Effect estimation: n > 10p (treating as p confounders)
- Causal discovery: n > p² or worse

Most real datasets insufficient for reliable structure learning.

### Violations of Overlap

If some treatment combinations rare (sparse regions in confounder space), estimation unreliable. Example: elderly rarely receive certain treatments.

Options: restrict to common support region (reduces generalizability), use extrapolation methods (risky), rely on model assumptions.

### Definition Ambiguity

What constitutes "treatment"? Must be:
- Well-defined: multiple versions wouldn't produce different outcomes
- Controllable: could be assigned in principle

Violations: education (versions unclear), "use social media" (can't separate mechanisms).

### Causal Claim Lifecycle

Causal conclusions erode with:
1. Different datasets → no confounding measured
2. Population shifts → confounders change
3. Mechanistic understanding → causal model evolves

Updating claims requires new data or refined assumptions.

## 9. Hyperparameter Tuning & Optimization

### Propensity Score Tuning

Model Complexity: Underfitting leaves residual confounding; overfitting wastes data. Cross-validation selects optimal complexity, but causal context different from prediction (goal: balance, not prediction accuracy).

Bandwidth Selection: For nearest-neighbor matching, bandwidth (caliper width) critical. Smaller → better matches but reduced sample. Cross-validate on balance, not prediction.

### Outcome Model Tuning

For regression adjustment and doubly robust methods, outcome model E_hat[Y|T,X] tuned similarly. Regularization (L2, L1) prevents overfitting.

Cross-Fitting: For doubly robust methods, divide data into K folds:
1. Estimate propensity/outcome models on fold k
2. Compute treatment effect estimates on other folds
3. Average across folds

Prevents overfitting-induced bias.

### BART Hyperparameters

Bayesian Additive Regression Trees popular for causal forests. Key parameters:
- Number of Trees: 200-400 typical
- Tree Depth: Shallow (2-5) to avoid overfitting
- Shrinkage: Controls regularization strength

### Causal Forest Parameters

• Number of Trees: Larger ensemble better but computationally expensive
• Minimum Node Size: Larger → smoother estimates, biased
• Splitting Criterion: Honest splitting avoids overfitting effect heterogeneity

Cross-validation guides selection, accounting for causal efficiency (not prediction MSE).

## 10. Real-World Applications & Case Studies

### Healthcare: Treatment Effect Heterogeneity

Electronic health records (EHR) enable personalized medicine. Example: which heart attack patients benefit most from aggressive antiplatelet therapy?

Causal forests on EHR data (10,000 patients, 500 covariates) identified effect modifiers:
- High benefit: elderly with diabetes
- Low/negative benefit: young without prior cardiac disease

Heterogeneity validates precision medicine: tailor treatment to patient characteristics.

### Marketing: Campaign Effectiveness

A/B test limited to small scale; observational data from past campaigns larger. Propensity score methods estimate campaign lift (incremental sales from ad exposure).

Challenge: selection bias (likely users more likely exposed). Solution: IPW reweighting by propensity score creates pseudo-experiment.

Result: 30% of observed difference in purchase rates due to campaign, 70% to selection.

### Policy: Job Training Programs

NSW study compared treatment effect estimators against experimental gold standard. Found:
- Naive estimator severely biased
- Propensity score, double robust methods close to truth
- Domain knowledge about confounders critical for method selection

Lesson: observational methods can work but require careful confounder measurement and multiple robustness checks.

### Epidemiology: COVID-19 Interventions

Identifying causal effects of lockdowns, mask mandates, vaccination with observational data complicated by:
- Simultaneous interventions (confounding)
- Time-varying treatments
- Spatial spillovers

Causal methods (especially difference-in-differences, synthetic controls) inform policy despite imperfect identifiability.

## 11. Integration with Other Methods

### Combining with Machine Learning

Causal + Predictive: Causal estimates inform features for ML models. Example: estimate patient subgroup treatment effects, use heterogeneity as feature in risk prediction.

Honest Forests: Use separate sample (honest splitting) for estimation vs. estimation to avoid overfitting treatment effect heterogeneity.

Double Machine Learning: Debiased ML approach combining ML flexibility with causal validity:
1. Fit nuisance models (propensity, outcome) via ML
2. Residualize: remove predicted confounding
3. Estimate causal parameters on residuals
4. Valid inference despite ML estimation

### Causal Discovery in Multivariate Settings

Combine constraint-based (PC algorithm) and score-based (GES, NOTEARS) approaches. Two-stage:
1. Score-based for skeleton (undirected graph)
2. Constraint-based for orientation

Improves over either approach alone.

### Measurement Error in Confounders

When confounders measured with error, causal estimates biased. Latent confounder models handle:
- Multiple imperfect measurements of true confounder
- Factor analysis estimates latent confounder
- Causal inference proceeds on factor scores

Accounting for measurement error → more honest uncertainty quantification.

## 12. Future Research Directions

### Causal Discovery from Interventional Data

Randomized experiments (interventions on variables) dramatically improve discovery. Moving beyond observational discovery to incorporate experimental intervention data.

### Dynamical Causal Systems

Extend beyond static snapshots to time-evolving systems with feedback loops. Challenges: temporal confounding, time-varying treatments, long-term effects.

### Fairness Through Causality

Formal fairness definitions (demographic parity, equalized odds) connected to causal mechanisms. Identify discriminatory causal paths and intervene at appropriate points.

### Large Language Models & Causal Reasoning

LLMs encode causal knowledge but inconsistently. Research: extract causal graphs from LLMs, ground them in data, verify against observational tests.

### Nonparametric Causal Methods

Move beyond linear/additive models to fully nonparametric causal estimation. Kernel methods, deep learning for nuisance models with valid inference.

## 13. Summary & Key Takeaways

Causal inference transforms ML from prediction to understanding mechanisms and effects of interventions. Key insights:

1. Causality Distinct from Correlation: Association doesn't imply causation; confounding creates spurious correlations that causal methods address.

2. Identifiability Requires Assumptions: Unconfoundedness, overlap, consistency must hold for valid effect estimation. Sensitivity analysis quantifies robustness.

3. Multiple Methods Available: Matching, IPW, regression, doubly robust, IV—each suitable for different contexts. No universal best method.

4. Heterogeneity is Reality: Average effects mislead; treatment effects vary by individual. Modern methods (forests, BART, CATE) capture heterogeneity.

5. Robustness Critical: Multiple specifications, sensitivity analysis, domain knowledge essential. Single analysis insufficient.

6. Integration with ML: Causal inference benefits from ML flexibility; ML benefits from causal grounding for robustness and interpretability.

7. Discovery Possible but Hard: Structure learning requires strong assumptions and large samples. Hybrid constraint/score-based methods most practical.

---

## Appendix: Practical Labs

### Lab 1: Propensity Score Matching

import numpy as np
import pandas as pd
from sklearn.linear_model import LogisticRegression
from sklearn.neighbors import NearestNeighbors
from scipy.spatial.distance import cdist

# Generate observational data with confounding
np.random.seed(42)
n = 500
X_confounder = np.random.randn(n) # Confounder affects both T and Y
T = (0.5 + 0.5*X_confounder + 0.3*np.random.randn(n) > 0).astype(int)
Y = 2*T + X_confounder + np.random.randn(n) # True causal effect is 2

# Naive estimate (biased)
naive_ate = Y[T==1].mean() - Y[T==0].mean()
print(f"Naive ATE: {naive_ate:.4f} (True: 2.0)")

# 1. Estimate propensity scores
prop_model = LogisticRegression()
prop_model.fit(X_confounder.reshape(-1,1), T)
propensity = prop_model.predict_proba(X_confounder.reshape(-1,1))[:, 1]

# 2. Find matches with similar propensity scores
treated_idx = np.where(T == 1)[0]
control_idx = np.where(T == 0)[0]
caliper = 0.1

matched_pairs = []
for t_idx in treated_idx:
 ps_treated = propensity[t_idx]
 dists = np.abs(propensity[control_idx] - ps_treated)
 close_controls = control_idx[dists < caliper]
 if len(close_controls) > 0:
 best_match = close_controls[np.argmin(dists[dists < caliper])]
 matched_pairs.append((t_idx, best_match))

# 3. Estimate ATE on matched sample
if len(matched_pairs) > 0:
 treated_Y = Y[np.array(matched_pairs)[:, 0]]
 control_Y = Y[np.array(matched_pairs)[:, 1]]
 ps_matched_ate = treated_Y.mean() - control_Y.mean()
 print(f"PS Matched ATE: {ps_matched_ate:.4f} (True: 2.0)")

### Lab 2: Inverse Probability Weighting

import numpy as np
from sklearn.linear_model import LinearRegression

# Continue with previous data
np.random.seed(42)
n = 500
X_confounder = np.random.randn(n)
T = (0.5 + 0.5*X_confounder + 0.3*np.random.randn(n) > 0).astype(int)
Y = 2*T + X_confounder + np.random.randn(n)

# Propensity scores
from sklearn.linear_model import LogisticRegression
prop_model = LogisticRegression()
prop_model.fit(X_confounder.reshape(-1,1), T)
propensity = prop_model.predict_proba(X_confounder.reshape(-1,1))[:, 1]

# IPW: reweight observations
weights = np.zeros(n)
weights[T==1] = 1 / (propensity[T==1] + 0.01) # Weight for treated
weights[T==0] = 1 / (1 - propensity[T==0] + 0.01) # Weight for control
weights /= weights.sum() / n # Normalize

# Weighted ATE
treated_mean = np.average(Y[T==1], weights=weights[T==1])
control_mean = np.average(Y[T==0], weights=weights[T==0])
ipw_ate = treated_mean - control_mean
print(f"IPW ATE: {ipw_ate:.4f} (True: 2.0)")

# Doubly robust: combine IPW and outcome regression
# Fit outcome model
outcome_model = LinearRegression()
X_full = np.column_stack([T, X_confounder])
outcome_model.fit(X_full, Y)
Y_pred_1 = outcome_model.predict(np.column_stack([np.ones(n), X_confounder]))
Y_pred_0 = outcome_model.predict(np.column_stack([np.zeros(n), X_confounder]))

# Doubly robust formula
dr_term_1 = (T * (Y - Y_pred_1) / (propensity + 0.01)) + Y_pred_1
dr_term_0 = ((1-T) * (Y - Y_pred_0) / (1 - propensity + 0.01)) + Y_pred_0
dr_ate = dr_term_1.mean() - dr_term_0.mean()
print(f"Doubly Robust ATE: {dr_ate:.4f} (True: 2.0)")

### Lab 3: Causal Discovery with PC Algorithm

import numpy as np
from itertools import combinations

def gaussian_ci_test(X, i, j, S, alpha=0.05):
 """Conditional independence test: X_i ⊥ X_j | X_S"""
 from scipy.stats import norm
 if len(S) == 0:
 # Marginal independence
 rho = np.corrcoef(X[:,[i,j]].T)[0,1]
 else:
 # Partial correlation
 from sklearn.linear_model import LinearRegression
 reg_i = LinearRegression().fit(X[:,S], X[:,i])
 reg_j = LinearRegression().fit(X[:,S], X[:,j])
 residuals_i = X[:,i] - reg_i.predict(X[:,S])
 residuals_j = X[:,j] - reg_j.predict(X[:,S])
 rho = np.corrcoef([residuals_i, residuals_j])[0,1]
 
 # Fisher Z-test
 z = 0.5 * np.log((1 + rho)/(1 - rho + 1e-6) + 1e-6)
 threshold = norm.ppf(1 - alpha/2) / np.sqrt(len(X) - len(S) - 3)
 return abs(z) < threshold

# Generate data from known DAG: X1 -> X2 -> X3
np.random.seed(42)
n = 500
X1 = np.random.randn(n)
X2 = X1 + np.random.randn(n)
X3 = X2 + np.random.randn(n)
X = np.column_stack([X1, X2, X3])

# PC Algorithm: discover skeleton
adj = np.ones((3, 3)) - np.eye(3) # Start with fully connected

# Condition on empty set (marginal independence)
for i, j in combinations(range(3), 2):
 if gaussian_ci_test(X, i, j, set()):
 adj[i,j] = adj[j,i] = 0

# Condition on single variable
for i, j in combinations(range(3), 2):
 if adj[i,j] == 0: continue
 for k in range(3):
 if k not in [i,j]:
 if gaussian_ci_test(X, i, j, {k}):
 adj[i,j] = adj[j,i] = 0
 break

print("Learned skeleton (adjacency):")
print(adj)
print("Expected: [[0,1,0], [1,0,1], [0,1,0]]")

### Lab 4: Heterogeneous Treatment Effects with Causal Forests

import numpy as np
from sklearn.ensemble import RandomForestRegressor

# Generate data with heterogeneous effects
np.random.seed(42)
n = 1000
X = np.random.uniform(-1, 1, (n, 5))
T = (X[:, 0] > 0).astype(int)
# Heterogeneous effect: stronger for X1 > 0
tau = 2 * (X[:, 0] + 1)
Y = tau * T + X[:, 1] + np.random.randn(n) * 0.5

# Simple causal forest estimation
# 1. Honest splitting: split each tree's data
from sklearn.model_selection import train_test_split
train_idx = np.arange(n)
np.random.shuffle(train_idx)
split = n // 2

# Estimate treatment propensity
prop_forest = RandomForestRegressor(max_depth=5, random_state=42)
prop_forest.fit(X[train_idx[:split]], T[train_idx[:split]])
propensity = prop_forest.predict(X)

# Estimate outcome model for each treatment
y_forest_t1 = RandomForestRegressor(max_depth=5, random_state=42)
y_forest_t1.fit(X[train_idx[:split]][T[train_idx[:split]]==1], 
 Y[train_idx[:split]][T[train_idx[:split]]==1])

y_forest_t0 = RandomForestRegressor(max_depth=5, random_state=42)
y_forest_t0.fit(X[train_idx[:split]][T[train_idx[:split]]==0],
 Y[train_idx[:split]][T[train_idx[:split]]==0])

# Estimate heterogeneous treatment effects on held-out data
test_idx = train_idx[split:]
cate_estimates = y_forest_t1.predict(X[test_idx]) - y_forest_t0.predict(X[test_idx])
true_tau_test = tau[test_idx]

# Evaluate: correlation between estimates and truth
from scipy.stats import pearsonr
corr, _ = pearsonr(cate_estimates, true_tau_test)
print(f"Estimated vs True HTE Correlation: {corr:.4f}")
print(f"Mean Estimated HTE: {cate_estimates.mean():.4f} (True Mean: {true_tau_test.mean():.4f})")