Econometrics: pulling cause out of observational data¶
Library-first with statsmodels. We watch omitted-variable bias inflate a pooled estimate, remove it with fixed effects, prove the fix by revealing the hidden confounder, and run a two-stage least squares instrumental-variable demo.
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
from statsmodels.stats.outliers_influence import variance_inflation_factor
from statsmodels.nonparametric.smoothers_lowess import lowess
from sklearn.linear_model import LinearRegression, Ridge, Lasso, LogisticRegression
from sklearn.preprocessing import StandardScaler, PolynomialFeatures
from sklearn.pipeline import make_pipeline
from sklearn.model_selection import train_test_split, cross_val_score
from sklearn.metrics import r2_score, mean_squared_error
ORG="#ea580c"; DEEP="#c2410c"; LIGHT="#fdba74"; INK="#1a2138"; GRID="#e6e9f2"; GREEN="#059669"; RED="#ef4444"; AMBER="#d97706"; BLUE="#2563eb"; PUR="#9333ea"; GREY="#94a3b8"; SLATE="#475569"
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/"
try:
panel = pd.read_excel('../../data/econometrics-and-panel-data--panel.xlsx', sheet_name='Panel')
except FileNotFoundError:
panel = pd.read_excel(BASE + 'econometrics-and-panel-data--panel.xlsx', sheet_name='Panel')
print(f'{panel.worker_id.nunique()} workers x {panel.year.nunique()} years = {len(panel)} rows')
print(f'corr(union, ability) = {panel[["union","ability_index"]].corr().iloc[0,1]:.2f} (the endogeneity)')
panel.head(3)
160 workers x 4 years = 640 rows corr(union, ability) = 0.37 (the endogeneity)
| worker_id | year | tenure | union | ability_index | log_wage | wage | |
|---|---|---|---|---|---|---|---|
| 0 | W000 | 2021 | 8 | 1 | 0.082 | 3.210 | 24.78 |
| 1 | W000 | 2022 | 9 | 0 | 0.082 | 3.039 | 20.89 |
| 2 | W000 | 2023 | 10 | 0 | 0.082 | 3.263 | 26.12 |
pooled = ols('log_wage ~ union + tenure', panel).fit()
b, (lo, hi) = pooled.params['union'], pooled.conf_int().loc['union']
print(f'Pooled OLS union coefficient = {b:.3f} (95% CI {lo:.3f} to {hi:.3f})')
print('True premium built into the data = 0.15 -> pooled OLS is biased UP.')
Pooled OLS union coefficient = 0.544 (95% CI 0.462 to 0.626) True premium built into the data = 0.15 -> pooled OLS is biased UP.
Bias is not variance. Notice the pooled estimate (about 0.54) comes with a tight confidence interval, it looks precise. But precise and correct are different things. The interval is narrow because we have plenty of data, yet it does not contain the true 0.15, because the error here is bias (a systematic tilt from the omitted ability), not sampling noise. More data would only shrink the interval further around the wrong number. Only a better design, not a bigger sample, can remove bias.
fe = ols('log_wage ~ union + tenure + C(worker_id)', panel).fit()
b, (lo, hi) = fe.params['union'], fe.conf_int().loc['union']
print(f'Fixed-effects union coefficient = {b:.3f} (95% CI {lo:.3f} to {hi:.3f})')
print('Recovered the true ~0.15 premium, without ever measuring ability.')
Fixed-effects union coefficient = 0.152 (95% CI 0.123 to 0.181) Recovered the true ~0.15 premium, without ever measuring ability.
# Show the correction visually: pooled vs FE vs truth
ctrl = ols('log_wage ~ union + tenure + ability_index', panel).fit()
ests = {'pooled OLS': pooled, 'control ability': ctrl, 'fixed effects': fe}
fig, ax = plt.subplots(figsize=(7.6,3.6))
for i,(name,m) in enumerate(ests.items()):
lo,hi = m.conf_int().loc['union']; c = m.params['union']
col = RED if name=='pooled OLS' else GREEN
ax.plot([lo,hi],[i,i], color=col, lw=3); ax.scatter(c,i,color=col,s=60,zorder=3)
ax.axvline(0.15, color=INK, ls='--', label='true = 0.15')
ax.set_yticks(range(len(ests))); ax.set_yticklabels(list(ests.keys()))
ax.set(title='Union premium: bias, then correction', xlabel='union coefficient (log wage)'); ax.legend()
plt.tight_layout(); plt.show()
print(f"control for ability : union = {ctrl.params['union']:.3f}")
print(f"fixed effects : union = {fe.params['union']:.3f}")
print('Nearly identical: FE absorbs exactly the ability confounder, without needing to observe it.')
# Direction of omitted-variable bias = sign(corr(union,ability)) * sign(effect of ability)
print('\nOVB direction: union-ability corr > 0 and ability raises wage > 0 -> upward bias, as seen.')
control for ability : union = 0.157 fixed effects : union = 0.152 Nearly identical: FE absorbs exactly the ability confounder, without needing to observe it. OVB direction: union-ability corr > 0 and ability raises wage > 0 -> upward bias, as seen.
rng = np.random.default_rng(7); n = 2000
u = rng.normal(0, 1, n) # unobserved confounder
z = rng.normal(0, 1, n) # instrument (independent of u)
x = 0.8*z + 0.9*u + rng.normal(0,0.5,n) # endogenous: depends on z AND u
y = 2.0*x + 1.5*u + rng.normal(0,0.5,n) # true causal effect of x on y is 2.0
d = pd.DataFrame({'y':y,'x':x,'z':z})
ols_biased = ols('y ~ x', d).fit().params['x']
stage1 = ols('x ~ z', d).fit(); d['x_hat'] = stage1.fittedvalues
iv_2sls = ols('y ~ x_hat', d).fit().params['x_hat']
print(f'true causal effect = 2.00')
print(f'naive OLS (biased by u) = {ols_biased:.2f}')
print(f'2SLS with instrument z = {iv_2sls:.2f} (recovers the truth)')
true causal effect = 2.00 naive OLS (biased by u) = 2.77 2SLS with instrument z = 2.00 (recovers the truth)
How two-stage least squares works. When the troublesome predictor X is correlated with the error (through a time-varying confounder that fixed effects cannot remove), we use an instrument Z, a variable that shifts X but affects Y only through X. Stage 1 regresses X on Z and keeps the fitted X-hat, the part of X explained by the instrument, scrubbed clean of the confounder. Stage 2 regresses Y on that X-hat. Here naive OLS was biased by the confounder u, but 2SLS with the valid instrument z recovered the true causal slope of 2.0. The catch in practice: a valid instrument must be relevant (it really moves X) and satisfy the exclusion restriction (it touches Y through no other path), which can be argued but never fully proven.
Putting it all together: designing for cause¶
Econometrics asks a harder question than the earlier chapters: not 'does the model fit?' but 'does the design identify a causal effect?' Here is the full run.
Step 1, check the threat. We wanted the causal effect of union membership on wages. Union status was correlated with unobserved ability (r about 0.37), so the exogeneity condition failed, a plain regression would be biased.
Steps 2-3, fit naively and diagnose. Pooled OLS put the union premium at 0.54, more than triple the true 0.15, with a tight interval that did NOT contain the truth. That is bias, not variance: a confident wrong answer.
Step 4, fix the design with fixed effects. Because we had panel data (the same workers over time), fixed effects added a per-worker intercept and used only within-worker variation, differencing out time-invariant ability without ever measuring it. The estimate fell to 0.15, matching what we got by controlling for ability directly, proof the method works.
When fixed effects are not enough: instrumental variables. Fixed effects only remove time-invariant confounders. For a time-varying confounder, reverse causation, or measurement error, we used two-stage least squares with an instrument: naive OLS was biased, 2SLS recovered the true slope.
Step 5, interpret. A union coefficient of 0.15 on the log scale means wages are e^0.15 = 1.16 times higher, a 16% premium.
Takeaway. Clean residuals are necessary but not sufficient. The deepest question is whether the design identifies a causal effect, and fixed effects, instruments, and difference-in-differences are how econometrics answers it, closing the correlation-is-not-causation loop opened in Chapter 91.