{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "1f1b832b",
   "metadata": {},
   "source": [
    "# Causal Inference for Lending\n",
    "\n",
    "Companion notebook to [Data Science Module 8](https://rodoslay.com/data-science/08-causal-inference-in-lending) at **rodoslay.com**.\n",
    "\n",
    "Four estimators on data where the truth is known: the **rate-confounding parable** (naive vs.\n",
    "regression-adjusted vs. propensity-matched vs. randomized), an **unobserved-confounder sensitivity\n",
    "sweep**, a **regression discontinuity** at a score cutoff, and **difference-in-differences** on a\n",
    "staged limit-increase rollout."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f27b4191",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.linear_model import LogisticRegression, LinearRegression\n",
    "\n",
    "rng = np.random.default_rng(42)\n",
    "plt.rcParams[\"figure.figsize\"] = (10, 4.5)\n",
    "plt.rcParams[\"axes.grid\"] = True\n",
    "plt.rcParams[\"grid.alpha\"] = 0.3\n",
    "sigmoid = lambda x: 1 / (1 + np.exp(-x))\n",
    "logit = lambda p: np.log(p / (1 - p))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "391b8ef4",
   "metadata": {},
   "source": [
    "## 1. The rate-confounding parable\n",
    "\n",
    "True effect of a high rate: **+0.30 log-odds** (odds ratio ≈ 1.35). Risky borrowers get high rates\n",
    "far more often — and risk is only *partially* observed through a noisy score."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ad78c8a3",
   "metadata": {},
   "outputs": [],
   "source": [
    "N = 100_000\n",
    "TAU = 0.30                                          # the truth we'll try to recover\n",
    "risk = rng.standard_normal(N)                       # true riskiness (never fully observed)\n",
    "score = risk + rng.normal(0, 0.6, N)                # the lender's noisy risk measure\n",
    "\n",
    "p_high = sigmoid(1.4 * risk)                        # risk-based pricing: risk -> rate\n",
    "high_rate = rng.random(N) < p_high\n",
    "p_default = sigmoid(-3.0 + 1.1 * risk + TAU * high_rate)\n",
    "default = rng.random(N) < p_default\n",
    "\n",
    "def log_odds_effect(d, t):\n",
    "    p1, p0 = d[t].mean(), d[~t].mean()\n",
    "    return logit(p1) - logit(p0)\n",
    "\n",
    "naive = log_odds_effect(default, high_rate)\n",
    "\n",
    "# regression adjustment on the OBSERVED score\n",
    "Xr = np.column_stack([high_rate.astype(float), score])\n",
    "adj = LogisticRegression(max_iter=1000).fit(Xr, default).coef_[0][0]\n",
    "\n",
    "# propensity-score stratification (20 strata on P(high rate | score))\n",
    "ps = LogisticRegression(max_iter=1000).fit(score.reshape(-1, 1), high_rate).predict_proba(score.reshape(-1, 1))[:, 1]\n",
    "strata = np.digitize(ps, np.quantile(ps, np.linspace(0, 1, 21)[1:-1]))\n",
    "effs, ws = [], []\n",
    "for s in np.unique(strata):\n",
    "    m = strata == s\n",
    "    if default[m & high_rate].size > 50 and default[m & ~high_rate].size > 50:\n",
    "        if 0 < default[m & high_rate].mean() < 1 and 0 < default[m & ~high_rate].mean() < 1:\n",
    "            effs.append(log_odds_effect(default[m], high_rate[m])); ws.append(m.sum())\n",
    "psm = np.average(effs, weights=ws)\n",
    "\n",
    "# the randomized benchmark: same borrowers, coin-flip rates\n",
    "high_rct = rng.random(N) < 0.5\n",
    "default_rct = rng.random(N) < sigmoid(-3.0 + 1.1 * risk + TAU * high_rct)\n",
    "rct = log_odds_effect(default_rct, high_rct)\n",
    "\n",
    "pd.DataFrame({\n",
    "    \"estimator\": [\"naive comparison\", \"regression-adjusted (noisy score)\",\n",
    "                  \"propensity-stratified (noisy score)\", \"randomized experiment\", \"TRUTH\"],\n",
    "    \"log-odds effect\": [naive, adj, psm, rct, TAU],\n",
    "    \"odds ratio\": np.exp([naive, adj, psm, rct, TAU]),\n",
    "}).round(3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "67c3ff63",
   "metadata": {},
   "source": [
    "The pattern to memorize: the naive estimate is wildly inflated; adjustment and propensity methods\n",
    "close *most* of the gap but not all — because the score measures risk with noise, and **you can only\n",
    "adjust for the confounding you can see**. Only randomization nails the truth (up to sampling error).\n",
    "\n",
    "## 2. How bad is unobserved confounding? The sensitivity sweep"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3f234907",
   "metadata": {},
   "outputs": [],
   "source": [
    "noises = [0.0, 0.3, 0.6, 1.0, 1.5, 2.5]\n",
    "rows = []\n",
    "for nz in noises:\n",
    "    sc = risk + rng.normal(0, nz, N) if nz > 0 else risk.copy()\n",
    "    a = LogisticRegression(max_iter=1000).fit(np.column_stack([high_rate, sc]), default).coef_[0][0]\n",
    "    rows.append([nz, a, np.exp(a)])\n",
    "sens = pd.DataFrame(rows, columns=[\"score noise σ\", \"adjusted estimate\", \"odds ratio\"])\n",
    "plt.plot(sens[\"score noise σ\"], sens[\"adjusted estimate\"], \"o-\", color=\"#1f4d3a\")\n",
    "plt.axhline(TAU, color=\"k\", ls=\"--\", label=f\"truth ({TAU})\")\n",
    "plt.axhline(naive, color=\"#b3341e\", ls=\":\", label=f\"naive ({naive:.2f})\")\n",
    "plt.xlabel(\"noise in the observed risk measure\"); plt.ylabel(\"adjusted log-odds estimate\")\n",
    "plt.legend(); plt.title(\"Adjustment quality degrades exactly as your confounder measurement does\")\n",
    "plt.show()\n",
    "sens.round(3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "90ec79a9",
   "metadata": {},
   "source": [
    "## 3. Regression discontinuity: the natural experiment you already own\n",
    "\n",
    "Applicants just above and just below a score cutoff are statistically identical — comparing their\n",
    "outcomes estimates the causal effect of approval (here: of receiving credit on a downstream outcome)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9a3897f1",
   "metadata": {},
   "outputs": [],
   "source": [
    "N2 = 60_000\n",
    "app_score = rng.normal(650, 40, N2)\n",
    "CUT = 640\n",
    "approved = app_score >= CUT\n",
    "TRUE_RDD = -0.04                                       # approval improves the outcome by 4pp\n",
    "\n",
    "# outcome depends smoothly on score AND on the treatment\n",
    "p_bad_outcome = np.clip(0.30 - 0.0022 * (app_score - 650) + TRUE_RDD * approved, 0.01, 0.9)\n",
    "bad = rng.random(N2) < p_bad_outcome\n",
    "\n",
    "# local linear fit on each side of the cutoff, bandwidth 25 points\n",
    "bw = 25\n",
    "m = np.abs(app_score - CUT) < bw\n",
    "left = m & ~approved; right = m & approved\n",
    "fl = LinearRegression().fit(app_score[left].reshape(-1, 1), bad[left])\n",
    "fr = LinearRegression().fit(app_score[right].reshape(-1, 1), bad[right])\n",
    "rdd_est = fr.predict([[CUT]])[0] - fl.predict([[CUT]])[0]\n",
    "print(f\"RDD estimate at the cutoff: {rdd_est:+.3f}   (truth {TRUE_RDD:+.3f})\")\n",
    "\n",
    "# the classic picture: binned means with the discontinuity\n",
    "bins = np.arange(560, 740, 5)\n",
    "centers, means = [], []\n",
    "for lo in bins:\n",
    "    mm = (app_score >= lo) & (app_score < lo + 5)\n",
    "    if mm.sum() > 100:\n",
    "        centers.append(lo + 2.5); means.append(bad[mm].mean())\n",
    "plt.scatter(centers, means, s=18, color=\"#1f4d3a\")\n",
    "plt.axvline(CUT, color=\"#b3341e\", ls=\"--\")\n",
    "plt.xlabel(\"application score\"); plt.ylabel(\"bad outcome rate\")\n",
    "plt.title(\"The jump AT the cutoff is the causal effect — the slope is just risk\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "eba89b64",
   "metadata": {},
   "source": [
    "## 4. Difference-in-differences: the staged rollout\n",
    "\n",
    "Region A gets a limit-increase program at month 12; region B doesn't. Regions differ in level\n",
    "(A was always riskier) and share a common macro trend — both would poison a naive comparison."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fcb27dab",
   "metadata": {},
   "outputs": [],
   "source": [
    "months = np.arange(24)\n",
    "TRUE_DID = 0.010                                       # program adds +1pp to default rate\n",
    "macro = 0.002 * months                                 # shared worsening trend\n",
    "rate_A = 0.05 + macro + TRUE_DID * (months >= 12) + rng.normal(0, 0.0015, 24)\n",
    "rate_B = 0.03 + macro + rng.normal(0, 0.0015, 24)\n",
    "\n",
    "pre_A, post_A = rate_A[:12].mean(), rate_A[12:].mean()\n",
    "pre_B, post_B = rate_B[:12].mean(), rate_B[12:].mean()\n",
    "did = (post_A - pre_A) - (post_B - pre_B)\n",
    "naive_cross = post_A - post_B                          # 'treated region is worse' — confounded by level\n",
    "naive_time = post_A - pre_A                            # 'defaults rose after launch' — confounded by trend\n",
    "\n",
    "print(f\"naive cross-section (A vs B, post):   {naive_cross:+.4f}\")\n",
    "print(f\"naive before-after (A only):          {naive_time:+.4f}\")\n",
    "print(f\"difference-in-differences:            {did:+.4f}   (truth {TRUE_DID:+.4f})\")\n",
    "\n",
    "plt.plot(months, rate_A, \"o-\", color=\"#b3341e\", label=\"region A (treated at m12)\")\n",
    "plt.plot(months, rate_B, \"o-\", color=\"#1f4d3a\", label=\"region B (control)\")\n",
    "plt.axvline(11.5, color=\"k\", ls=\"--\", lw=0.8)\n",
    "plt.legend(); plt.title(\"Parallel trends, then a wedge — DiD measures the wedge\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dd4877ba",
   "metadata": {},
   "source": [
    "Both naive numbers are wrong in different directions; DiD cancels the level difference *and*\n",
    "the shared trend, and lands on the truth. The identifying assumption — parallel trends absent\n",
    "treatment — is visible in the pre-period: check it before you trust it.\n",
    "\n",
    "## Exercise: the credit-line-increase memo\n",
    "\n",
    "You've now estimated the program's effect four ways. Write the one-page recommendation: the\n",
    "estimate you trust most, the assumption it leans on, the analysis you'd run next quarter to check\n",
    "it, and — per the module — what you would randomize if the committee lets you. State clearly which\n",
    "question is causal and which is predictive; that sentence is the one the regulator reads twice."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "version": "3.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}