{ "cells": [ { "cell_type": "markdown", "id": "a1l6-title", "metadata": {}, "source": [ "# A1 · L6 — Which features earn their place\n", "\n", "R² cannot go down when a column is added. This notebook proves it with columns of\n", "random numbers, then measures what a column has to be worth before adjusted R² and\n", "AIC will let it in — and where newspaper lands against that bar.\n", "\n", "Every figure and number this notebook produces is what the lesson prints." ] }, { "cell_type": "code", "execution_count": 1, "id": "a1l6-header", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:09.569141Z", "iopub.status.busy": "2026-08-19T09:12:09.569029Z", "iopub.status.idle": "2026-08-19T09:12:09.846948Z", "shell.execute_reply": "2026-08-19T09:12:09.815790Z" } }, "outputs": [], "source": [ "# hide — plumbing so this file runs both in the repo and in Colab.\n", "try:\n", " from _figkit import save_fig, record\n", "except ImportError: # Colab — no repo, no problem\n", " def save_fig(name, plot, **kw):\n", " import matplotlib.pyplot as plt\n", " fig, ax = plt.subplots(figsize=kw.get('figsize', (7, 4.2))); plot(ax); plt.show()\n", " def record(key, value):\n", " print(f'{key} = {value}'); return value\n", "\n", "from pathlib import Path\n", "REL = 'public/datasets/a1-regression/advertising-media-mix.csv'\n", "here = Path.cwd()\n", "LOCAL = next((p / REL for p in [here, *here.parents] if (p / REL).exists()), None)\n", "CSV = LOCAL or 'https://raw.githubusercontent.com/Dr-Shashank-S-Sharma/expedify-ai-courses/main/datasets/a1-regression/advertising-media-mix.csv'" ] }, { "cell_type": "code", "execution_count": 2, "id": "a1l6-fit", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:09.848722Z", "iopub.status.busy": "2026-08-19T09:12:09.848575Z", "iopub.status.idle": "2026-08-19T09:12:10.048398Z", "shell.execute_reply": "2026-08-19T09:12:10.048132Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "model R2 adj R2 AIC\n", "TV 0.61188 0.60991 476.52\n", "TV + radio 0.89719 0.89615 212.82\n", "TV + radio + newspaper 0.89721 0.89564 214.79\n" ] } ], "source": [ "#| caption: Three scores for the same fit — R², adjusted R² and AIC\n", "import numpy as np\n", "import pandas as pd\n", "\n", "df = pd.read_csv(CSV)\n", "y = df['sales'].to_numpy()\n", "n = len(df)\n", "\n", "def score(columns):\n", " \"\"\"Least squares on any list of columns (names or arrays). Three scores back.\n", "\n", " R2 share of variation explained\n", " adj R2 charged a fee per column: 1 - (1-R2)(n-1)/(n-p-1)\n", " AIC n*ln(RSS/n) + 2k, k = p + 2 (the slopes, the intercept, the error variance).\n", " Only DIFFERENCES in AIC mean anything, so the constant is left out.\n", " \"\"\"\n", " cols = [df[c].to_numpy() if isinstance(c, str) else c for c in columns]\n", " X = np.column_stack([np.ones(n)] + cols)\n", " coef, *_ = np.linalg.lstsq(X, y, rcond=None)\n", " rss = float(((y - X @ coef) ** 2).sum())\n", " p = len(cols)\n", " r2 = 1 - rss / ((y - y.mean()) ** 2).sum()\n", " adj = 1 - (1 - r2) * (n - 1) / (n - p - 1)\n", " aic = n * np.log(rss / n) + 2 * (p + 2)\n", " return dict(r2=r2, adj=adj, aic=aic, p=p)\n", "\n", "LADDER = {\n", " 'TV': ['tv_spend'],\n", " 'TV + radio': ['tv_spend', 'radio_spend'],\n", " 'TV + radio + newspaper': ['tv_spend', 'radio_spend', 'newspaper_spend'],\n", "}\n", "\n", "print(f\"{'model':24s} {'R2':>8s} {'adj R2':>9s} {'AIC':>8s}\")\n", "for label, columns in LADDER.items():\n", " s = score(columns)\n", " print(f\"{label:24s} {s['r2']:8.5f} {s['adj']:9.5f} {s['aic']:8.2f}\")" ] }, { "cell_type": "code", "execution_count": 3, "id": "a1l6-ladder-metrics", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:10.049659Z", "iopub.status.busy": "2026-08-19T09:12:10.049541Z", "iopub.status.idle": "2026-08-19T09:12:10.054744Z", "shell.execute_reply": "2026-08-19T09:12:10.054495Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n_markets = 200\n", "r2_tv = 0.61188\n", "r2_tv_radio = 0.89719\n", "r2_full = 0.89721\n", "adj_tv_radio = 0.89615\n", "adj_full = 0.89564\n", "aic_tv_radio = 212.8\n", "aic_full = 214.8\n", "adj_drop_newspaper = 0.00051\n", "aic_rise_newspaper = 2.0\n" ] }, { "data": { "text/plain": [ "np.float64(2.0)" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# hide — metrics for the ladder\n", "two = score(LADDER['TV + radio'])\n", "three = score(LADDER['TV + radio + newspaper'])\n", "one = score(LADDER['TV'])\n", "\n", "record('n_markets', n)\n", "record('r2_tv', round(one['r2'], 5))\n", "record('r2_tv_radio', round(two['r2'], 5))\n", "record('r2_full', round(three['r2'], 5))\n", "record('adj_tv_radio', round(two['adj'], 5))\n", "record('adj_full', round(three['adj'], 5))\n", "record('aic_tv_radio', round(two['aic'], 1))\n", "record('aic_full', round(three['aic'], 1))\n", "record('adj_drop_newspaper', round(two['adj'] - three['adj'], 5))\n", "record('aic_rise_newspaper', round(three['aic'] - two['aic'], 1))" ] }, { "cell_type": "markdown", "id": "a1l6-md-noise", "metadata": {}, "source": [ "## 1. Columns of pure noise\n", "\n", "Add columns of random numbers — numbers that cannot possibly explain sales, because\n", "they were made up after the fact — and watch each score react. Averaged over many\n", "draws, so what you see is the behaviour of the score, not the luck of one draw." ] }, { "cell_type": "code", "execution_count": 4, "id": "a1l6-noise", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:10.055993Z", "iopub.status.busy": "2026-08-19T09:12:10.055898Z", "iopub.status.idle": "2026-08-19T09:12:10.206827Z", "shell.execute_reply": "2026-08-19T09:12:10.206476Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "junk cols R2 adj R2 AIC\n", " 0 0.89721 0.89564 214.79\n", " 1 0.89759 0.89549 216.05\n", " 5 0.89966 0.89546 219.94\n", " 10 0.90215 0.89531 224.90\n", " 15 0.90491 0.89545 229.13\n" ] } ], "source": [ "#| caption: What junk columns do to each score\n", "DRAWS, MAX_NOISE = 200, 15\n", "rng = np.random.default_rng(3)\n", "noise = rng.standard_normal((DRAWS, MAX_NOISE, n))\n", "\n", "CHANNELS = ['tv_spend', 'radio_spend', 'newspaper_spend']\n", "curve = {'r2': [], 'adj': [], 'aic': []}\n", "for j in range(MAX_NOISE + 1):\n", " runs = [score(CHANNELS + [noise[d, k] for k in range(j)]) for d in range(DRAWS)]\n", " for key in curve:\n", " curve[key].append(float(np.mean([r[key] for r in runs])))\n", "\n", "print(f\"{'junk cols':>9s} {'R2':>9s} {'adj R2':>9s} {'AIC':>8s}\")\n", "for j in (0, 1, 5, 10, 15):\n", " print(f\"{j:9d} {curve['r2'][j]:9.5f} {curve['adj'][j]:9.5f} {curve['aic'][j]:8.2f}\")" ] }, { "cell_type": "code", "execution_count": 5, "id": "a1l6-noise-figs", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:10.208027Z", "iopub.status.busy": "2026-08-19T09:12:10.207945Z", "iopub.status.idle": "2026-08-19T09:12:10.341063Z", "shell.execute_reply": "2026-08-19T09:12:10.340793Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n_draws = 200\n", "max_noise = 15\n", "r2_noise_15 = 0.90491\n", "r2_rise_noise = 0.0077\n", "adj_noise_15 = 0.89545\n", "adj_move_noise = -0.00018\n", "aic_noise_15 = 229.1\n", "aic_rise_noise = 14.3\n", "r2_rising_steps = 15 of 15\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "figure noise-climb -> which-features-earn-their-place.noise-climb.{light,dark}.svg\n", "figure adjusted -> which-features-earn-their-place.adjusted.{light,dark}.svg\n" ] } ], "source": [ "# hide — metrics + the two noise figures\n", "record('n_draws', DRAWS)\n", "record('max_noise', MAX_NOISE)\n", "record('r2_noise_15', round(curve['r2'][MAX_NOISE], 5))\n", "record('r2_rise_noise', round(curve['r2'][MAX_NOISE] - curve['r2'][0], 5))\n", "record('adj_noise_15', round(curve['adj'][MAX_NOISE], 5))\n", "record('adj_move_noise', round(curve['adj'][MAX_NOISE] - curve['adj'][0], 5))\n", "record('aic_noise_15', round(curve['aic'][MAX_NOISE], 1))\n", "record('aic_rise_noise', round(curve['aic'][MAX_NOISE] - curve['aic'][0], 1))\n", "\n", "# The lesson says the line rises at every step. That is a claim about a PICTURE,\n", "# which is the one thing the truth rule cannot check — so count the steps here\n", "# and let the prose quote the count instead of describing the shape.\n", "rising = sum(curve['r2'][j] < curve['r2'][j + 1] for j in range(MAX_NOISE))\n", "record('r2_rising_steps', f'{rising} of {MAX_NOISE}')\n", "\n", "steps = list(range(MAX_NOISE + 1))\n", "\n", "def plot_r2(ax):\n", " ax.plot(steps, curve['r2'], marker='o', markersize=3.5)\n", " ax.set_xlabel('Columns of random numbers added')\n", " ax.set_ylabel('R² — share of variation explained')\n", " ax.annotate(f\"{curve['r2'][0]:.4f}\", xy=(0, curve['r2'][0]), xytext=(8, -12),\n", " textcoords='offset points', fontsize=9)\n", " ax.annotate(f\"{curve['r2'][MAX_NOISE]:.4f}\", xy=(MAX_NOISE, curve['r2'][MAX_NOISE]),\n", " xytext=(-48, -14), textcoords='offset points', fontsize=9)\n", "\n", "save_fig('noise-climb', plot_r2, figsize=(7, 3.8))\n", "\n", "def plot_both(ax):\n", " ax.plot(steps, curve['r2'], marker='o', markersize=3.5, label='R²')\n", " ax.plot(steps, curve['adj'], marker='s', markersize=3.5, label='adjusted R²')\n", " ax.set_xlabel('Columns of random numbers added')\n", " ax.set_ylabel('Score')\n", " ax.legend(frameon=False, fontsize=9)\n", "\n", "save_fig('adjusted', plot_both, figsize=(7, 3.8))" ] }, { "cell_type": "markdown", "id": "a1l6-md-random", "metadata": {}, "source": [ "## 2. Is newspaper better than a random column?\n", "\n", "Newspaper adds something to R². So does a column of random numbers. The question\n", "worth asking is which of the two adds more — so draw 500 random columns, add each\n", "to the TV + radio model, and see where newspaper's gain falls in that distribution." ] }, { "cell_type": "code", "execution_count": 6, "id": "a1l6-random", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:10.342501Z", "iopub.status.busy": "2026-08-19T09:12:10.342405Z", "iopub.status.idle": "2026-08-19T09:12:10.358921Z", "shell.execute_reply": "2026-08-19T09:12:10.358696Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "newspaper adds 0.000016 to R²\n", "a random column adds 0.000258 at the median, 0.000520 on average, up to 0.005101\n", "theory says a junk column is worth 0.000525\n", "random columns that beat newspaper: 85.4%\n" ] } ], "source": [ "#| caption: Newspaper's R² gain against 500 random columns\n", "N_RANDOM = 500\n", "rng = np.random.default_rng(11)\n", "base = ['tv_spend', 'radio_spend']\n", "r2_base = score(base)['r2']\n", "\n", "gains = np.array([score(base + [rng.standard_normal(n)])['r2'] - r2_base\n", " for _ in range(N_RANDOM)])\n", "gain_newspaper = score(base + ['newspaper_spend'])['r2'] - r2_base\n", "\n", "# A useless column is not worth zero — it is worth the variation it mops up by\n", "# luck, and that has a closed form: the variation still unexplained, shared out\n", "# over the degrees of freedom left. Print both to see they agree.\n", "expected_junk = (1 - r2_base) / (n - 3 - 1)\n", "\n", "print(f'newspaper adds {gain_newspaper:.6f} to R²')\n", "print(f'a random column adds {np.median(gains):.6f} at the median, '\n", " f'{gains.mean():.6f} on average, up to {gains.max():.6f}')\n", "print(f'theory says a junk column is worth {expected_junk:.6f}')\n", "print(f'random columns that beat newspaper: {(gains > gain_newspaper).mean():.1%}')" ] }, { "cell_type": "code", "execution_count": 7, "id": "a1l6-random-fig", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:10.360061Z", "iopub.status.busy": "2026-08-19T09:12:10.359963Z", "iopub.status.idle": "2026-08-19T09:12:10.457079Z", "shell.execute_reply": "2026-08-19T09:12:10.456831Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n_random = 500\n", "gain_newspaper = 0.000016\n", "gain_random_median = 0.000258\n", "gain_random_mean = 0.000520\n", "gain_random_max = 0.005101\n", "gain_random_min = 0.000000\n", "expected_junk_gain = 0.000525\n", "share_random_beating_newspaper = 85%\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "figure random-vs-newspaper -> which-features-earn-their-place.random-vs-newspaper.{light,dark}.svg\n" ] } ], "source": [ "# hide — metrics + figure\n", "record('n_random', N_RANDOM)\n", "record('gain_newspaper', f'{gain_newspaper:.6f}')\n", "record('gain_random_median', f'{np.median(gains):.6f}')\n", "record('gain_random_mean', f'{gains.mean():.6f}')\n", "record('gain_random_max', f'{gains.max():.6f}')\n", "record('gain_random_min', f'{gains.min():.6f}')\n", "record('expected_junk_gain', f'{expected_junk:.6f}')\n", "record('share_random_beating_newspaper', f'{(gains > gain_newspaper).mean():.0%}')\n", "\n", "def plot(ax):\n", " ax.hist(gains, bins=40, zorder=1)\n", " top = ax.get_ylim()[1]\n", " ax.axvline(gain_newspaper, linestyle='--', linewidth=1.8, color='#e2574c', zorder=4)\n", " ax.axvline(np.median(gains), linestyle='-', linewidth=1.8, color='#2f9e6e', zorder=4)\n", " ax.annotate('newspaper', xy=(gain_newspaper, top * 0.92), xytext=(46, 0),\n", " textcoords='offset points', fontsize=9, color='#e2574c', va='center',\n", " arrowprops=dict(arrowstyle='->', color='#e2574c', linewidth=1.2))\n", " ax.annotate('median random column', xy=(np.median(gains), top * 0.62), xytext=(52, 0),\n", " textcoords='offset points', fontsize=9, color='#2f9e6e', va='center',\n", " arrowprops=dict(arrowstyle='->', color='#2f9e6e', linewidth=1.2))\n", " ax.set_xlabel('R² gained by adding this column to TV + radio')\n", " ax.set_ylabel(f'Random columns (of {N_RANDOM})')\n", "\n", "save_fig('random-vs-newspaper', plot, figsize=(7, 3.8))" ] }, { "cell_type": "markdown", "id": "a1l6-md-subsets", "metadata": {}, "source": [ "## 3. Score every model you could have built\n", "\n", "Three channels make eight possible models, including the one with nothing in it.\n", "Score all eight and let AIC rank them." ] }, { "cell_type": "code", "execution_count": 8, "id": "a1l6-subsets", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:10.458314Z", "iopub.status.busy": "2026-08-19T09:12:10.458216Z", "iopub.status.idle": "2026-08-19T09:12:10.462382Z", "shell.execute_reply": "2026-08-19T09:12:10.462161Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "model R2 adj R2 AIC\n", "TV + radio 0.8972 0.8962 212.8\n", "TV + radio + newspaper 0.8972 0.8956 214.8\n", "TV + newspaper 0.6458 0.6422 460.2\n", "TV 0.6119 0.6099 476.5\n", "radio 0.3320 0.3287 585.1\n", "radio + newspaper 0.3327 0.3259 586.9\n", "newspaper 0.0521 0.0473 655.1\n", "(nothing) 0.0000 0.0000 663.8\n" ] } ], "source": [ "#| caption: All eight models, ranked by AIC (lower is better)\n", "from itertools import combinations\n", "\n", "subsets = [c for r in range(4) for c in combinations(CHANNELS, r)]\n", "ranked = sorted(((score(list(c)), c) for c in subsets), key=lambda pair: pair[0]['aic'])\n", "\n", "PRETTY = {'tv_spend': 'TV', 'radio_spend': 'radio', 'newspaper_spend': 'newspaper'}\n", "name = lambda combo: ' + '.join(PRETTY[c] for c in combo) or '(nothing)'\n", "\n", "print(f\"{'model':32s} {'R2':>8s} {'adj R2':>9s} {'AIC':>8s}\")\n", "for s, combo in ranked:\n", " print(f\"{name(combo):32s} {s['r2']:8.4f} {s['adj']:9.4f} {s['aic']:8.1f}\")" ] }, { "cell_type": "code", "execution_count": 9, "id": "a1l6-subsets-fig", "metadata": { "execution": { "iopub.execute_input": "2026-08-19T09:12:10.463518Z", "iopub.status.busy": "2026-08-19T09:12:10.463378Z", "iopub.status.idle": "2026-08-19T09:12:10.526725Z", "shell.execute_reply": "2026-08-19T09:12:10.526484Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n_subsets = 8\n", "best_model = TV + radio\n", "best_aic = 212.8\n", "runner_up = TV + radio + newspaper\n", "best_by_adj = TV + radio\n", "best_by_r2 = TV + radio + newspaper\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "figure all-subsets -> which-features-earn-their-place.all-subsets.{light,dark}.svg\n" ] } ], "source": [ "# hide — metrics + figure\n", "record('n_subsets', len(subsets))\n", "record('best_model', name(ranked[0][1]))\n", "record('best_aic', round(ranked[0][0]['aic'], 1))\n", "record('runner_up', name(ranked[1][1]))\n", "record('best_by_adj', name(max(ranked, key=lambda pair: pair[0]['adj'])[1]))\n", "record('best_by_r2', name(max(ranked, key=lambda pair: pair[0]['r2'])[1]))\n", "\n", "def plot(ax):\n", " # A lollipop, not bars: AIC has no meaningful zero, so a bar drawn from one\n", " # would be an invented baseline and a bar drawn from anywhere else is a lie\n", " # about proportion. Position is the only honest encoding here.\n", " labels = [name(combo) for _, combo in ranked][::-1]\n", " values = [s['aic'] for s, _ in ranked][::-1]\n", " best = min(values)\n", " ax.hlines(labels, best, values, linewidth=1.4, color='#8a8a94', zorder=1)\n", " ax.scatter(values, labels, s=52,\n", " color=['#2f9e6e' if v == best else '#3b6fd4' for v in values], zorder=3)\n", " ax.set_xlabel('AIC — lower is better')\n", " ax.set_xlim(best - 30, max(values) + 45)\n", " ax.grid(axis='y', visible=False)\n", " for i, v in enumerate(values):\n", " ax.annotate(f'{v:.0f}', xy=(v, i), xytext=(8, 0), textcoords='offset points',\n", " va='center', fontsize=9)\n", "\n", "save_fig('all-subsets', plot, figsize=(7, 4.2))" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.5" } }, "nbformat": 4, "nbformat_minor": 5 }