nb3.4 — 2SLS with `linearmodels`; the First-Stage F

Companion to Ch 3.4, “Instrumental Variables.”

The chapter made a promise that sounds like magic: even when an unobserved confounder poisons your regressor — so that ordinary least squares (OLS) returns a biased number no amount of control variables can fix — a valid instrument Z can still recover the true causal effect. This notebook turns that promise into running code and watches it succeed, then watches it fail.

We do five things, in order.

Build a world where OLS is guaranteed to lie. We simulate a treatment D that is endogenous: an unmeasured “motivation” variable drives both D and the outcome Y. Because motivation sits in the error term, OLS confuses the effect of the course with the effect of being the kind of person who takes the course.
Run 2SLS three ways and prove they agree. By hand (first stage → fitted \hat D → second stage), with the cov-ratio Frisch–Waugh–Lovell (FWL) formula from Ch 2.3, and with the production routine linearmodels.iv.IV2SLS. All three return the same point estimate — but only IV2SLS gives honest standard errors.
Read the first-stage F and hold it against the “F>10” rule of thumb and the Stock–Yogo logic.
Sabotage the instrument — shrink \operatorname{Cov}(Z,D) to near zero — and watch 2SLS slide back toward the OLS bias while its standard error explodes and the F collapses to nothing.
Plot the whole story: the 2SLS estimate, its standard error, and the first-stage F as functions of instrument strength, tracing the boundary between a trustworthy IV and a dangerous one.

The chapter’s headline numbers — OLS \approx -3.62 (biased), 2SLS \approx -1.92 against a true -2.0, first-stage F \approx 286 — are reproduced exactly in §2. Everything here previews Chapter 3.5, where the weak-instrument disaster gets its own lab and its own cure (Anderson–Rubin inference).

Setup

We fix one seeded generator so the notebook reproduces bit-for-bit (the CONVENTIONS.md rule), and force matplotlib’s non-interactive Agg backend so the notebook runs headless — laptop, server, or CI — with no display attached.

The one import to notice is the last one. In linearmodels 7.0 the 2SLS estimator lives at linearmodels.iv.IV2SLS. That is the import path we use throughout.

import matplotlib
matplotlib.use("Agg")  # headless, non-interactive backend

import numpy as np
import pandas as pd
import statsmodels.api as sm
import matplotlib.pyplot as plt
from linearmodels.iv import IV2SLS  # linearmodels 7.0 import path

rng = np.random.default_rng(42)  # one seeded generator for the whole notebook

print("numpy", np.__version__, "| pandas", pd.__version__)
print("statsmodels", sm.__version__)
import linearmodels
print("linearmodels", linearmodels.__version__)

1. A world built to break OLS

We reproduce the chapter’s simulation exactly. Maya wants the effect of completing a financial-literacy course (D) on the change in credit-card balance (Y, coded so lower is better). The catch is financial motivation: an unmeasured trait that makes a person both more likely to take the course and better at managing debt on their own. Motivation lives in the error term, so it is the textbook unobserved confounder.

Here is the data-generating process, with every channel labeled. The true causal effect is \beta_1 = -2.0: completing the course lowers the balance change by 2 units.

Confounder motivation \sim \mathcal N(0,1) — unobserved; we simulate it but never let any estimator see it.
Instrument Z — a random mailer, a coin flip. Because it is assigned by chance, it is independent of motivation: the exclusion restriction’s hard half holds by construction.
Treatment D — crosses an enrollment threshold driven by the mailer (0.8 * Z, the relevance channel) and by motivation (the endogeneity). It is binary: enrolled or not.
Outcome Y — gets the true effect beta1_true * D, but motivation also enters Y directly (-1.5 * motivation), and that term is invisible to any regression of Y on D, so it contaminates the error.

The identifying assumption, stated in the discipline of CONVENTIONS §4: outcome = balance change; treatment = course completion D; instrument = random mailer Z; controls = none (just a constant) in this first pass; identifying assumption = the mailer affects debt outcomes only by raising course completion, and (being randomized) is independent of motivation.

# --- The chapter's exact simulation (its own seed, to match the printed numbers) ---
sim = np.random.default_rng(20260528)
N = 5_000
beta1_true = -2.0

# Unobserved confounder: financial motivation (lowers balance change, raises take-up)
motivation = sim.normal(0.0, 1.0, size=N)

# Instrument: random mailer (coin flip -> independent of motivation; exclusion holds by design)
Z = sim.integers(0, 2, size=N).astype(float)

# Treatment: pushed by the mailer (relevance) AND by motivation (the endogeneity)
latent = 0.8 * Z + 1.0 * motivation + sim.normal(0.0, 1.0, size=N)
D = (latent > 0.7).astype(float)

# Outcome: true effect of D is beta1_true; motivation ALSO enters Y -> lands in the error
eps = sim.normal(0.0, 1.0, size=N)
Y = 1.0 + beta1_true * D - 1.5 * motivation + eps

df = pd.DataFrame({"Y": Y, "D": D, "Z": Z, "const": 1.0})

print(f"N = {N},  true beta1 = {beta1_true}")
print(f"share treated (D=1): {D.mean():.3f}")
print(f"share mailed  (Z=1): {Z.mean():.3f}")
# The smoking gun: motivation correlates with D (and is hidden in Y's error)
print(f"corr(motivation, D) = {np.corrcoef(motivation, D)[0, 1]:+.3f}  <- the endogeneity")

2. OLS is biased — and we can see exactly why

Run the regression Maya is tempted to report: Y on a constant and D. With no endogenous regressor declared, IV2SLS simply runs OLS, so we get a clean apples-to-apples comparison with the 2SLS call that follows.

The bias has a sign you can predict. Motivated people both take the course (D=1) and have lower balance changes on their own (the -1.5 * motivation term). So treated people look better than the course alone can explain, and OLS hands the course credit it did not earn — pushing \hat\beta_1 more negative than the truth. We expect roughly -3.62 against a true -2.0.

# Naive OLS via IV2SLS: dependent, exog=[const, D], no endog, no instruments.
ols = IV2SLS(df["Y"], df[["const", "D"]], None, None).fit()
beta_ols = ols.params["D"]

print(f"OLS  beta_hat(D) = {beta_ols:.3f}   (true = {beta1_true})")
print(f"OLS bias         = {beta_ols - beta1_true:+.3f}   "
      f"(more negative than truth, exactly as predicted)")

Now the instrument earns its keep. We declare D endogenous and hand IV2SLS the mailer Z as its instrument. The call signature is

IV2SLS(dependent, exog=[const + controls], endog=[D], instruments=[Z])

We fit with cov_type="robust" (heteroskedasticity-robust standard errors — the Ch 2.4 default for real data). The point estimate should land near the true -2.0, and the chapter’s printed value is -1.92.

# 2SLS: instrument the endogenous D with Z. linearmodels returns CORRECT standard errors.
iv = IV2SLS(df["Y"], df[["const"]], df[["D"]], df[["Z"]]).fit(cov_type="robust")
beta_2sls = iv.params["D"]

print(f"OLS  beta_hat(D) = {beta_ols:.3f}   (biased)")
print(f"2SLS beta_hat(D) = {beta_2sls:.3f}   (true = {beta1_true})   SE = {iv.std_errors['D']:.3f}")
print()
print("The OLS-vs-2SLS gap IS the endogeneity bias, made visible:")
print(f"   gap = {beta_ols - beta_2sls:+.3f}")

3. 2SLS by hand, and why the hand-rolled standard errors are a trap

The chapter insists you understand the two stages mechanically before trusting the black box. So we build \hat\beta_1^{\text{2SLS}} from scratch:

First stage. Regress the endogenous D on the instrument and the exogenous columns (here just the constant): D = \pi_0 + \pi_1 Z + \nu. Keep the fitted values \hat D — the slice of D explained by exogenous variation, laundered of motivation.
Second stage. Regress Y on \hat D (and the same exogenous columns). The slope is the 2SLS point estimate.

We compute both stages with plain matrix algebra and compare to IV2SLS. They must match to many decimals.

Then the warning the chapter shouts in bold: the second-stage regression thinks \hat D is real data, but \hat D is an estimate carrying first-stage noise. So the hand-rolled standard error is wrong — usually too small. We print both to see the lie.

# --- First stage: D on [const, Z], keep fitted D-hat ---
X_first = df[["const", "Z"]].to_numpy()
d = df["D"].to_numpy()
pi_hat, *_ = np.linalg.lstsq(X_first, d, rcond=None)
D_hat = X_first @ pi_hat
print(f"first-stage coef on Z (pi_1) = {pi_hat[1]:.4f}  (>0: the mailer raises take-up)")

# --- Second stage: Y on [const, D-hat] ---
X_second = np.column_stack([np.ones(N), D_hat])
y = df["Y"].to_numpy()
b_hand, *_ = np.linalg.lstsq(X_second, y, rcond=None)
beta_hand = b_hand[1]

# Naive (WRONG) second-stage SE: treats D-hat as if it were data
resid_second = y - X_second @ b_hand
sigma2 = (resid_second @ resid_second) / (N - 2)
XtX_inv = np.linalg.inv(X_second.T @ X_second)
se_hand_wrong = np.sqrt(sigma2 * XtX_inv[1, 1])

print(f"\nby-hand 2SLS beta(D) = {beta_hand:.6f}")
print(f"IV2SLS  2SLS beta(D) = {beta_2sls:.6f}")
print(f"agree to 6 dp? {np.isclose(beta_hand, beta_2sls, atol=1e-6)}")
print()
print(f"hand-rolled SE (WRONG, ignores first-stage noise) = {se_hand_wrong:.4f}")
print(f"IV2SLS SE      (CORRECT)                           = {iv.std_errors['D']:.4f}")
print("-> the hand-rolled SE understates the true uncertainty; always report IV2SLS's.")

4. The first-stage F: how strong is the lever?

Relevance is the testable condition. The standard measure of how hard the instrument shoves the treatment is the first-stage F-statistic: the F for the joint significance of the excluded instrument(s) — here, just Z — in the first-stage regression. With a single instrument it is literally the squared t-statistic on Z.

linearmodels hands you this in iv.first_stage.diagnostics. We read off f.stat (the partial F) and f.pval. Because we built a strong first stage (0.8 * Z), the F should land far above 10 — the chapter’s value is \approx 286.

How to read it. The famous rule of thumb is F>10, from Stock–Yogo (2005): clearing it gives a formal guarantee that the 2SLS bias is at most about 10% of the OLS bias. But treat 10 as a smoke alarm, not a certificate — it assumes classical (homoskedastic, non-clustered) errors and bounds bias, not the honesty of your t-statistic. With clustered or heteroskedastic finance data you report the Olea–Pflueger effective F instead (Ch 3.5). At F \approx 286 we are comfortably, unambiguously strong, so the distinction does not bite here.

diag = iv.first_stage.diagnostics
print("First-stage diagnostics (linearmodels):")
print(diag[["rsquared", "f.stat", "f.pval"]])

F_strong = float(diag["f.stat"].iloc[0])
print(f"\nfirst-stage partial F = {F_strong:.1f}")

# With ONE instrument, F = (t-stat on Z in the first stage)^2. Confirm it.
fs = iv.first_stage.individual["D"]          # the first-stage regression of D
t_on_Z = fs.tstats["Z"]
print(f"t-stat on Z in first stage = {t_on_Z:.3f}   ->   t^2 = {t_on_Z**2:.1f}")
print(f"matches F? {np.isclose(t_on_Z**2, F_strong, rtol=1e-3)}")

verdict = "STRONG (clears F>10 with enormous margin)" if F_strong > 10 else "WEAK"
print(f"\nVerdict: {verdict}")

5. 2SLS as Frisch–Waugh–Lovell, with controls

Real specifications carry control variables. The chapter showed (and Ch 2.3 proved) that 2SLS is FWL with an instrument: partial the controls out of Y, D, and Z alike, then the 2SLS slope is the cov-ratio

\hat\beta_1^{\text{2SLS}} = \frac{\operatorname{Cov}(\tilde Z,\tilde Y)}{\operatorname{Cov}(\tilde Z,\tilde D)},

where the tildes mean “residualized on the controls.” It is the Wald ratio of Ch 3.4 §4, living in the residualized world of Ch 2.3.

To show it cleanly we switch to a continuous treatment with two observed controls — age and income — exactly the “scale up” the chapter’s Your Turn asks for. The confounder motivation is still unobserved and still in both D and Y’s error. We then estimate \beta_1=-2.0 four ways and confirm all four agree:

IV2SLS with controls in exog,
by-hand two stages (controls in both stages),
the FWL cov-ratio after residualizing on [const, age, income],
(for contrast) OLS, which stays biased.

# --- Continuous-treatment world with two controls; motivation still unobserved ---
Nc = 4_000
age    = rng.normal(0.0, 1.0, Nc)     # observed control
income = rng.normal(0.0, 1.0, Nc)     # observed control
mot    = rng.normal(0.0, 1.0, Nc)     # UNOBSERVED confounder
Zc     = rng.normal(0.0, 1.0, Nc)     # continuous instrument (strong: 0.9 loading)

D_c = 0.5 + 0.9 * Zc + 0.4 * age - 0.3 * income + 1.0 * mot + rng.normal(0, 1, Nc)
eps_c = rng.normal(0.0, 1.0, Nc)
Y_c = 1.0 + beta1_true * D_c + 0.5 * age + 0.7 * income - 1.5 * mot + eps_c

dc = pd.DataFrame({"Y": Y_c, "D": D_c, "Z": Zc, "age": age, "income": income, "const": 1.0})

exog_cols = ["const", "age", "income"]

# (4) OLS with controls -- still biased, because motivation is uncontrolled
ols_c = IV2SLS(dc["Y"], dc[exog_cols + ["D"]], None, None).fit()

# (1) IV2SLS with controls
iv_c = IV2SLS(dc["Y"], dc[exog_cols], dc[["D"]], dc[["Z"]]).fit(cov_type="robust")

print(f"OLS    beta(D) = {ols_c.params['D']:+.4f}   (biased; controls cannot fix unobserved motivation)")
print(f"IV2SLS beta(D) = {iv_c.params['D']:+.4f}   (true = {beta1_true})")
print(f"first-stage F  = {float(iv_c.first_stage.diagnostics['f.stat'].iloc[0]):.1f}")

# (2) By-hand 2SLS WITH controls in both stages
Xf = dc[["const", "Z", "age", "income"]].to_numpy()
coef_f, *_ = np.linalg.lstsq(Xf, dc["D"].to_numpy(), rcond=None)
Dhat_c = Xf @ coef_f
Xs = np.column_stack([np.ones(Nc), Dhat_c, dc["age"].to_numpy(), dc["income"].to_numpy()])
coef_s, *_ = np.linalg.lstsq(Xs, dc["Y"].to_numpy(), rcond=None)
beta_hand_c = coef_s[1]

# (3) FWL cov-ratio: residualize Y, D, Z on [const, age, income], then ratio of covariances
C = dc[exog_cols].to_numpy()
def residualize(v):
    coef, *_ = np.linalg.lstsq(C, v, rcond=None)
    return v - C @ coef
Yt = residualize(dc["Y"].to_numpy())
Dt = residualize(dc["D"].to_numpy())
Zt = residualize(dc["Z"].to_numpy())
beta_fwl = np.cov(Zt, Yt, bias=True)[0, 1] / np.cov(Zt, Dt, bias=True)[0, 1]

print(f"IV2SLS   beta(D) = {iv_c.params['D']:.6f}")
print(f"by-hand  beta(D) = {beta_hand_c:.6f}")
print(f"FWL ratio beta(D)= {beta_fwl:.6f}")
print()
print(f"all three agree to 6 dp? "
      f"{np.allclose([beta_hand_c, beta_fwl], iv_c.params['D'], atol=1e-6)}")
print("-> 2SLS IS FWL with an instrument: partial out controls from Y, D, Z, then take the ratio.")

6. Sabotage: a weak instrument drags 2SLS back toward OLS

Now the cautionary tale. The whole IV machine divides by the first stage — “how much D moves with Z.” Shrink that denominator toward zero and three pathologies (Ch 3.4 §8) appear at once:

Bias toward OLS. With Z barely moving D, the fitted \hat D is mostly first-stage noise, and that noise carries the confounder. So 2SLS drifts back to the very OLS bias it was meant to escape.
Exploding standard error. Dividing by a near-zero denominator inflates the sampling variance, and the usual SE/t-stat stop being trustworthy.
Amplified exclusion violations (not shown here): any tiny forbidden Z\to Y path, divided by a near-zero first stage, blows up.

We rebuild the chapter’s binary world but shrink the relevance loading from 0.8 * Z to 0.05 * Z. Watch 2SLS abandon the truth and the F collapse toward zero.

# Rebuild the binary world, but with a NEARLY DEAD instrument (0.05 instead of 0.8)
weak = np.random.default_rng(20260528)   # same seed/structure as the strong run
mot_w = weak.normal(0, 1, N)
Z_w   = weak.integers(0, 2, N).astype(float)
latent_w = 0.05 * Z_w + 1.0 * mot_w + weak.normal(0, 1, N)   # <-- weak lever
D_w = (latent_w > 0.7).astype(float)
eps_w = weak.normal(0, 1, N)
Y_w = 1.0 + beta1_true * D_w - 1.5 * mot_w + eps_w
dw = pd.DataFrame({"Y": Y_w, "D": D_w, "Z": Z_w, "const": 1.0})

ols_w  = IV2SLS(dw["Y"], dw[["const", "D"]], None, None).fit()
iv_w   = IV2SLS(dw["Y"], dw[["const"]], dw[["D"]], dw[["Z"]]).fit(cov_type="robust")
F_weak = float(iv_w.first_stage.diagnostics["f.stat"].iloc[0])

print("                     beta(D)      SE        first-stage F")
print(f"true                 {beta1_true:+.3f}")
print(f"OLS (weak world)     {ols_w.params['D']:+.3f}")
print(f"2SLS STRONG (0.8*Z)  {beta_2sls:+.3f}     {iv.std_errors['D']:6.3f}    {F_strong:8.1f}")
print(f"2SLS WEAK   (0.05*Z) {iv_w.params['D']:+.3f}   {iv_w.std_errors['D']:8.3f}    {F_weak:8.2f}")
print()
print("The weak F is ~0 (far below 10). 2SLS no longer recovers -2.0, and its SE has exploded:")
print("IV done with a dead lever is WORSE than the OLS it was meant to fix. This is the Ch 3.5 disaster.")

7. The picture: estimate, precision, and F vs. instrument strength

One sabotage run is anecdote; a sweep is the lesson. We vary the instrument’s relevance loading \pi from near-zero (dead lever) to strong, and at each \pi run a small Monte Carlo of fresh datasets, recording the 2SLS estimate, its standard error, and the first-stage F.

Two summary choices matter and are themselves teaching points:

We plot the median 2SLS estimate across draws, not the mean. Weak-IV estimates have such fat tails that the mean is dominated by a few wild draws — a symptom of the pathology, not a number to trust.
We anchor two horizontal lines: the true effect -2.0 and the OLS bias (a badly biased \approx -3, from stronger confounding in this sweep’s DGP). The story is the median 2SLS line starting glued to the OLS bias when F is tiny, then peeling off and converging to the truth as F climbs past 10.

# Sweep instrument strength; small Monte Carlo at each level.
Ns = 2_000
B1 = -2.0
def one_dataset(pi, gen):
    a  = gen.normal(0, 1, Ns); inc = gen.normal(0, 1, Ns)
    m  = gen.normal(0, 1, Ns); z   = gen.normal(0, 1, Ns)
    d  = 0.5 + pi * z + 0.4 * a - 0.3 * inc + 2.0 * m + gen.normal(0, 1, Ns)  # strong confounding
    e  = gen.normal(0, 1, Ns)
    yy = 1.0 + B1 * d + 0.5 * a + 0.7 * inc - 3.0 * m + e
    frame = pd.DataFrame({"Y": yy, "D": d, "Z": z, "age": a, "income": inc, "const": 1.0})
    fit = IV2SLS(frame["Y"], frame[["const", "age", "income"]],
                 frame[["D"]], frame[["Z"]]).fit(cov_type="robust")
    Ff = float(fit.first_stage.diagnostics["f.stat"].iloc[0])
    return fit.params["D"], fit.std_errors["D"], Ff

sweep_rng = np.random.default_rng(7)
pis = [0.0, 0.02, 0.05, 0.10, 0.15, 0.20, 0.30, 0.50, 0.90]
n_mc = 120
med_beta, med_se, mean_F = [], [], []
for pi in pis:
    bs, ses, Fs = [], [], []
    for _ in range(n_mc):
        b, s, Ff = one_dataset(pi, sweep_rng)
        bs.append(b); ses.append(s); Fs.append(Ff)
    med_beta.append(np.median(bs)); med_se.append(np.median(ses)); mean_F.append(np.mean(Fs))
med_beta, med_se, mean_F = map(np.array, (med_beta, med_se, mean_F))

# OLS bias anchor in this DGP (run one big strong-instrument dataset, read its OLS coef)
big = np.random.default_rng(99)
a = big.normal(0,1,Ns); inc = big.normal(0,1,Ns); m = big.normal(0,1,Ns); z = big.normal(0,1,Ns)
d = 0.5 + 0.9*z + 0.4*a - 0.3*inc + 2.0*m + big.normal(0,1,Ns)
yy = 1.0 + B1*d + 0.5*a + 0.7*inc - 3.0*m + big.normal(0,1,Ns)
fr = pd.DataFrame({"Y":yy,"D":d,"age":a,"income":inc,"const":1.0})
ols_anchor = IV2SLS(fr["Y"], fr[["const","D","age","income"]], None, None).fit().params["D"]

print("pi    medianF   median 2SLS   median SE")
for pi, F, b, s in zip(pis, mean_F, med_beta, med_se):
    print(f"{pi:4.2f}  {F:7.1f}   {b:+10.3f}   {s:8.3f}")
print(f"\nOLS bias anchor (truth {B1}): {ols_anchor:+.3f}")

fig, axes = plt.subplots(1, 3, figsize=(15, 4.2))

# Panel A: median 2SLS estimate vs first-stage F
ax = axes[0]
ax.plot(mean_F, med_beta, "o-", color="#1f77b4", label="median 2SLS")
ax.axhline(B1, color="green", ls="--", label=f"true effect = {B1}")
ax.axhline(ols_anchor, color="red", ls=":", label=f"OLS bias = {ols_anchor:.2f}")
ax.axvline(10, color="gray", ls="-.", lw=1, label="F = 10 rule")
ax.set_xscale("log")
ax.set_xlabel("first-stage F (log scale)"); ax.set_ylabel("2SLS estimate of beta1")
ax.set_title("A. Weak IV hugs the OLS bias")
ax.legend(fontsize=8, loc="lower right")

# Panel B: median SE vs F -- precision collapse
ax = axes[1]
ax.plot(mean_F, med_se, "s-", color="#d62728")
ax.axvline(10, color="gray", ls="-.", lw=1)
ax.set_xscale("log"); ax.set_yscale("log")
ax.set_xlabel("first-stage F (log scale)"); ax.set_ylabel("median 2SLS SE (log scale)")
ax.set_title("B. Weak IV explodes the standard error")

# Panel C: first-stage F vs the relevance loading pi
ax = axes[2]
ax.plot(pis, mean_F, "^-", color="#2ca02c")
ax.axhline(10, color="gray", ls="-.", lw=1, label="F = 10 rule")
ax.set_yscale("log")
ax.set_xlabel("instrument loading  pi  (Cov(Z,D) strength)"); ax.set_ylabel("first-stage F (log)")
ax.set_title("C. Strength of the lever sets the F")
ax.legend(fontsize=8)

fig.suptitle("2SLS vs. instrument strength: the boundary between trustworthy and dangerous IV",
             fontsize=12)
fig.tight_layout(rect=[0, 0, 1, 0.96])
fig.savefig("nb3.4-instrument-strength.png", dpi=110, bbox_inches="tight")
plt.close(fig)
print("saved figure: nb3.4-instrument-strength.png")
print("Panel A: median 2SLS starts near the OLS bias at tiny F, converges to the truth past F=10.")
print("Panel B: the standard error balloons as F -> 0 (precision collapse).")

What to carry forward

Four things from this notebook will keep paying off.

OLS-vs-2SLS is the bias, visible. When an unobserved confounder makes D endogenous, OLS and a valid 2SLS disagree, and the gap is the endogeneity bias itself. Here OLS read -3.62 and 2SLS read -1.92 against a true -2.0; the 1.7-unit gap was motivation’s fingerprint.

Three roads to the same point estimate. By-hand two stages, the FWL cov-ratio, and IV2SLS all return the identical 2SLS coefficient — because 2SLS is FWL with the instrument standing in for the contaminated regressor. But only IV2SLS gives honest standard errors; the hand-rolled second-stage SE is wrong because it forgets that \hat D is an estimate. Hand-roll to understand; report with IV2SLS.

The first-stage F is your strength gauge. Report it always. Clear F>10 comfortably (Stock–Yogo), but remember it assumes classical errors and bounds bias, not inference — with clustered finance data, reach for the Olea–Pflueger effective F in Ch 3.5.

A weak instrument is worse than no instrument. As \operatorname{Cov}(Z,D)\to 0, 2SLS slides back to the OLS bias and its standard error explodes — you get a wrong answer dressed up as a precise one. The sweep made the boundary visible: trust the IV only when the lever is genuinely strong. Chapter 3.5 reproduces this disaster in full and arms you with Anderson–Rubin inference that stays valid even when the first stage is weak.

Your Turn — vary instrument strength and trace the bias yourself

You have seen the sweep we built; now run your own and reason about what you find.

Re-trace the bias. Re-run the §7 sweep but change the confounding strength: in one_dataset, flip the sign of motivation’s effect on Y from -3.0 * m to +3.0 * m. Now OLS is biased in the opposite direction. Re-plot Panel A. Does weak 2SLS still hug the (new) OLS bias? Which direction does it converge from as F grows?
Find the danger zone. Add finer \pi values between 0.05 and 0.20 to pis. For each, record the median F and the fraction of Monte Carlo draws whose 95% CI excludes the truth -2.0 (a draw’s CI is param ± 1.96 * SE). At roughly what F does that fraction fall to a tolerable ~5%? Compare your answer to the F>10 rule — is 10 enough here?
Mean vs. median. Re-plot Panel A using the mean 2SLS instead of the median. Explain in two sentences why the mean line is so much wilder at low F, and why the median is the honest summary. (Hint: think about what dividing by a near-zero, sometimes wrong-signed denominator does to the tails.)
Strengthen, don’t sabotage. Suppose you could add a second valid instrument. Add a second column Z2 (another rng.normal) with its own loading into D, pass both to IV2SLS as instruments=frame[["Z", "Z2"]], and watch the first-stage F. Why does adding a second relevant instrument tend to raise the F and stabilize 2SLS — and what new worry (over-identification) does it introduce that a single just-identified instrument did not have?