{ "nbformat": 4, "nbformat_minor": 5, "metadata": { "colab": { "provenance": [], "toc_visible": true }, "kernelspec": { "name": "python3", "display_name": "Python 3" }, "language_info": { "name": "python" } }, "cells": [ { "id": "m07-l3-0", "cell_type": "markdown", "metadata": {}, "source": [ "# Lesson 7.3 · DoubleML in Practice: PLR, IRM and a Real Health Application\n", "\n", "*Accelerated Foundation in Applied Causal Inference*\n", "\n", "**By the end of this lesson you can:**\n", "\n", "- Estimate partially linear (PLR) and interactive (IRM) models with the DoubleML library\n", "- Choose and compare nuisance learners (Lasso, random forests, LightGBM) by out-of-sample loss\n", "- Estimate the effect of quitting smoking on weight gain with the NHEFS data\n", "- Report a DML estimate with its assumptions, learners and robustness checks\n", "\n", "> Read the full explanations on the course website, then run each cell here in order (**Runtime → Run all** works on a fresh runtime)." ] }, { "id": "m07-l3-1", "cell_type": "markdown", "metadata": {}, "source": [ "## The DoubleML workflow" ] }, { "id": "m07-l3-2", "cell_type": "code", "metadata": {}, "execution_count": null, "outputs": [], "source": [ "!pip install -q DoubleML lightgbm scikit-learn causaldata\n", "import numpy as np, pandas as pd\n", "import doubleml as dml\n", "from lightgbm import LGBMRegressor, LGBMClassifier\n", "from sklearn.ensemble import RandomForestRegressor, RandomForestClassifier\n", "from sklearn.linear_model import LassoCV, LogisticRegressionCV\n", "from sklearn.pipeline import make_pipeline\n", "from sklearn.preprocessing import PolynomialFeatures, StandardScaler\n", "rng = np.random.default_rng(3)" ] }, { "id": "m07-l3-3", "cell_type": "markdown", "metadata": {}, "source": [ "## PLR on the digital-skills simulation" ] }, { "id": "m07-l3-4", "cell_type": "markdown", "metadata": {}, "source": [ "> **Question it answers:** What does DoubleML estimate for the digital-skills programme, with boosting as the nuisance learner?" ] }, { "id": "m07-l3-5", "cell_type": "code", "metadata": {}, "execution_count": null, "outputs": [], "source": [ "N, P = 3_000, 20\n", "X = rng.normal(0, 1, size=(N, P))\n", "p = 1 / (1 + np.exp(-(0.5 * X[:, 0] - 0.3 * X[:, 1] ** 2 + 0.2 * X[:, 2] * X[:, 4])))\n", "D = rng.binomial(1, p)\n", "Y = 5_000 * D + 2_000 * np.sin(X[:, 0]) + 3_000 * X[:, 3] + 1_000 * X[:, 1] ** 2 + rng.normal(0, 4_000, N)\n", "df = pd.DataFrame(X, columns=[f\"X_{i}\" for i in range(P)]).assign(D=D, Y=Y)\n", "\n", "data = dml.DoubleMLData(df, y_col=\"Y\", d_cols=\"D\")\n", "ml_l = LGBMRegressor(n_estimators=300, learning_rate=0.05, max_depth=4, verbose=-1, random_state=42)\n", "ml_m = LGBMClassifier(n_estimators=300, learning_rate=0.05, max_depth=4, verbose=-1, random_state=42)\n", "\n", "plr = dml.DoubleMLPLR(data, ml_l, ml_m, n_folds=5)\n", "plr.fit()\n", "print(plr.summary.round(2))\n", "print(\"95% CI:\\n\", plr.confint().round(0))" ] }, { "id": "m07-l3-6", "cell_type": "markdown", "metadata": {}, "source": [ "## Comparing learners" ] }, { "id": "m07-l3-7", "cell_type": "code", "metadata": {}, "execution_count": null, "outputs": [], "source": [ "learners = {\n", " \"Lasso (+ squares)\": (make_pipeline(PolynomialFeatures(2, include_bias=False), StandardScaler(), LassoCV(cv=3)),\n", " make_pipeline(StandardScaler(), LogisticRegressionCV(cv=3, max_iter=2_000))),\n", " \"Random forest\": (RandomForestRegressor(n_estimators=300, min_samples_leaf=5, n_jobs=-1, random_state=1),\n", " RandomForestClassifier(n_estimators=300, min_samples_leaf=5, n_jobs=-1, random_state=1)),\n", " \"LightGBM\": (ml_l, ml_m),\n", "}\n", "rows = []\n", "for name, (l, m) in learners.items():\n", " model = dml.DoubleMLPLR(data, l, m, n_folds=5).fit()\n", " losses = model.evaluate_learners()\n", " rows.append((name, model.coef[0], model.se[0], losses[\"ml_l\"].item(), losses[\"ml_m\"].item()))\n", "print(pd.DataFrame(rows, columns=[\"learner\", \"theta\", \"se\", \"RMSE E[Y|X]\", \"RMSE E[D|X]\"]).round(2))" ] }, { "id": "m07-l3-8", "cell_type": "markdown", "metadata": {}, "source": [ "## IRM on NHEFS: smoking cessation and weight" ] }, { "id": "m07-l3-9", "cell_type": "markdown", "metadata": {}, "source": [ "> **Question it answers:** How much weight do people gain because they quit smoking?" ] }, { "id": "m07-l3-10", "cell_type": "code", "metadata": {}, "execution_count": null, "outputs": [], "source": [ "from causaldata import nhefs_complete\n", "nh = nhefs_complete.load_pandas().data\n", "covs = [\"sex\", \"race\", \"age\", \"school\", \"smokeintensity\", \"smokeyrs\", \"wt71\", \"exercise\", \"active\"]\n", "nh = nh[covs + [\"qsmk\", \"wt82_71\"]].dropna()\n", "nh = pd.get_dummies(nh, columns=[\"exercise\", \"active\"], drop_first=True, dtype=float)\n", "\n", "naive = nh.wt82_71[nh.qsmk == 1].mean() - nh.wt82_71[nh.qsmk == 0].mean()\n", "nh_data = dml.DoubleMLData(nh, y_col=\"wt82_71\", d_cols=\"qsmk\")\n", "irm = dml.DoubleMLIRM(\n", " nh_data,\n", " ml_g=RandomForestRegressor(n_estimators=500, min_samples_leaf=10, max_features=0.5, random_state=0),\n", " ml_m=RandomForestClassifier(n_estimators=500, min_samples_leaf=10, max_features=0.5, random_state=0),\n", " n_folds=5, n_rep=3, trimming_threshold=0.01,\n", ")\n", "irm.fit()\n", "print(f\"Observations: {len(nh):,} | quitters: {int(nh.qsmk.sum())}\")\n", "print(f\"Naive difference: {naive:.2f} kg\")\n", "print(irm.summary.round(2))" ] }, { "id": "m07-l3-11", "cell_type": "markdown", "metadata": {}, "source": [ "## Reporting template" ] }, { "id": "m07-l3-12", "cell_type": "markdown", "metadata": {}, "source": [ "### ✍️ Exercise\n", "\n", "Re-run the NHEFS IRM with LightGBM learners and with the ATT score (`score=\"ATTE\"`). Put the three estimates\n", "(RF ATE, LightGBM ATE, RF ATT) in one table. Is the conclusion robust?" ] }, { "id": "m07-l3-13", "cell_type": "code", "metadata": {}, "execution_count": null, "outputs": [], "source": [ "# Your answer here" ] }, { "id": "m07-l3-14", "cell_type": "markdown", "metadata": {}, "source": [ "## Checklist and recap" ] } ] }