⚙️ Setup¶
import numpy as np, pandas as pd, matplotlib.pyplot as plt
from scipy import stats
import seaborn as sns # seaborn = high-level statistical plots (heatmaps, regplots, pairplots)
import statsmodels.api as sm
from statsmodels.formula.api import ols
PUR="#9333ea"; DEEP="#7e22ce"; LIGHT="#c084fc"; INK="#1a2138"; GRID="#e6e9f2"; GREEN="#059669"; RED="#ef4444"; AMBER="#d97706"; BLUE="#2563eb"
plt.rcParams.update({"figure.facecolor":"white","axes.facecolor":"white","figure.dpi":110,"font.size":11,
"axes.edgecolor":GRID,"axes.grid":True,"grid.color":GRID,"axes.axisbelow":True,"axes.spines.top":False,
"axes.spines.right":False,"axes.titlesize":12,"axes.titleweight":"bold","legend.frameon":False})
sns.set_style("whitegrid")
BASE="https://raw.githubusercontent.com/johnfisher-ai/Statistics-Data-Science-AI-Visual-Book/main/data/"
rng = np.random.default_rng(91)
print("x ~ y observed. Possible stories:")
print(" 1) x -> y (x causes y)")
print(" 2) y -> x (REVERSE causation)")
print(" 3) z -> x, z -> y (CONFOUNDING: a lurking common cause)")
print(" 4) coincidence (spurious, especially with many comparisons)")
print("\nThe correlation number is identical in all four; only DESIGN or domain knowledge decides.")
x ~ y observed. Possible stories: 1) x -> y (x causes y) 2) y -> x (REVERSE causation) 3) z -> x, z -> y (CONFOUNDING: a lurking common cause) 4) coincidence (spurious, especially with many comparisons) The correlation number is identical in all four; only DESIGN or domain knowledge decides.
This is the single most important caveat in statistics. A strong r is evidence that two variables are associated, never proof that one causes the other. The rest of the chapter is about telling the four stories apart.
z=rng.normal(0,1,400) # the hidden common cause
x=2.0*z + rng.normal(0,0.6,400) # z drives x
y=1.5*z + rng.normal(0,0.6,400) # z drives y (NOT x)
print(f"raw corr(x, y) = {stats.pearsonr(x,y)[0]:.2f} <- strong, but spurious")
def resid(v,w): return sm.OLS(v, sm.add_constant(w)).fit().resid
print(f"partial corr(x, y | z) = {stats.pearsonr(resid(x,z), resid(y,z))[0]:.2f} <- collapses once z is held constant")
fig,ax=plt.subplots(1,2,figsize=(10,3.4))
sc=ax[0].scatter(x,y,c=z,cmap="viridis",s=14); ax[0].set_title("x vs y, colored by z (the hidden driver)"); plt.colorbar(sc,ax=ax[0],label="z")
ax[1].scatter(resid(x,z),resid(y,z),s=14,color=PUR,alpha=0.5); ax[1].set_title("after removing z: no relationship left")
ax[1].set_xlabel("x residual"); ax[1].set_ylabel("y residual"); plt.tight_layout(); plt.show()
raw corr(x, y) = 0.88 <- strong, but spurious partial corr(x, y | z) = -0.02 <- collapses once z is held constant
Color the points by z and the illusion is obvious: the "relationship" is just z sliding both variables together. Strip z out (right panel) and the cloud goes round, the partial correlation near 0 is the statistical fingerprint of a confounder.
# spurious-by-chance: 1000 pairs of PURE NOISE, how strong does the luckiest correlation get?
best=max(abs(np.corrcoef(rng.normal(0,1,30), rng.normal(0,1,30))[0,1]) for _ in range(1000))
print(f"1000 unrelated 30-point pairs: the strongest |r| by pure luck = {best:.2f}")
print("With enough comparisons, a strong correlation will appear from nothing -> demand a plausible mechanism.")
1000 unrelated 30-point pairs: the strongest |r| by pure luck = 0.63 With enough comparisons, a strong correlation will appear from nothing -> demand a plausible mechanism.
Small samples and many comparisons breed impressive-looking correlations from noise. The defenses are the same throughout science: ask which direction is plausible?, look for a mechanism, and reserve causal claims for randomized experiments (Chapter 65) or careful causal inference.
Daily records show ice-cream sales and drownings rise together (correlation-vs-causation--confounding.xlsx). Does ice cream cause drownings? We explore, measure the raw correlation, then control for temperature, the lurking common cause.
try: d = pd.read_excel("../../data/correlation-vs-causation--confounding.xlsx", sheet_name="Days")
except FileNotFoundError: d = pd.read_excel(BASE+"correlation-vs-causation--confounding.xlsx", sheet_name="Days")
# EXPLORE FIRST: shape, missing, and the suspicious raw correlation
print("shape:", d.shape, "| missing:", d.isna().sum().sum())
print(d[["temperature_f","ice_cream_sales","drownings"]].describe().round(1).T[["mean","std","min","max"]])
print(f"\nraw corr(ice_cream_sales, drownings) = {d.ice_cream_sales.corr(d.drownings):.2f} (alarming, if taken at face value)")
shape: (200, 5) | missing: 0
mean std min max
temperature_f 75.4 12.4 42.0 102.0
ice_cream_sales 280.3 105.0 11.0 540.0
drownings 5.6 3.0 0.0 14.0
raw corr(ice_cream_sales, drownings) = 0.58 (alarming, if taken at face value)
# temperature is correlated with BOTH -> the prime suspect for confounding
print(f"corr(temperature, ice_cream) = {d.temperature_f.corr(d.ice_cream_sales):.2f}")
print(f"corr(temperature, drownings) = {d.temperature_f.corr(d.drownings):.2f}")
def resid(v): return sm.OLS(d[v], sm.add_constant(d.temperature_f)).fit().resid
pc=stats.pearsonr(resid("ice_cream_sales"), resid("drownings"))[0]
print(f"\nPARTIAL corr(ice_cream, drownings | temperature) = {pc:.2f} -> the link vanishes!")
fig,ax=plt.subplots(1,2,figsize=(11,3.6))
sc=ax[0].scatter(d.ice_cream_sales,d.drownings,c=d.temperature_f,cmap="coolwarm",s=22)
ax[0].set_xlabel("ice-cream sales"); ax[0].set_ylabel("drownings"); ax[0].set_title("Raw: rises together")
plt.colorbar(sc,ax=ax[0],label="temp (F)")
ax[1].scatter(resid("ice_cream_sales"),resid("drownings"),s=22,color=PUR,alpha=0.55)
ax[1].set_xlabel("ice-cream residual"); ax[1].set_ylabel("drownings residual"); ax[1].set_title("After removing temperature: no link")
plt.tight_layout(); plt.show()
corr(temperature, ice_cream) = 0.87 corr(temperature, drownings) = 0.68 PARTIAL corr(ice_cream, drownings | temperature) = -0.03 -> the link vanishes!
The raw correlation is about 0.58, but the color tells the truth: hot days (red) sit top-right, cold days (blue) bottom-left. Temperature drives both, more ice cream AND more swimming, hence more drownings. Control for it and the partial correlation collapses to about −0.03: ice cream and drownings have no direct link at all. Banning ice cream would not save a single swimmer; the cause is the heat.