{
 "nbformat": 4,
 "nbformat_minor": 5,
 "metadata": {
  "kernelspec": {
   "name": "python3",
   "display_name": "Python 3",
   "language": "python"
  },
  "language_info": {
   "name": "python"
  },
  "colab": {
   "name": "bagging.ipynb"
  }
 },
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "# bootstrap aggregating (bagging) \u2014 Python demo\n\nNumerical companion to the entry [bootstrap aggregating (bagging)](https://dictionaryofml.org/terms/bagging.html) of the [Dictionary of Applied Machine Learning](https://dictionaryofml.org/): it recomputes what the entry states and prints one line per check.\n\nBacks the entry's claims with numbers: a decision tree grown until it reproduces every label changes its predictions noticeably when a handful of training points are replaced, while linear regression barely does; bagging reduces the variance of the tree's predictions and leaves its bias unchanged; it does little for linear regression; and adding more base learners does not make the aggregated hypothesis fit the training set more closely. Self-contained (numpy/matplotlib only), fixed seed.\n\nRequires NumPy and Matplotlib only, and uses fixed seeds, so the printed numbers reproduce exactly. Generated from [`pythondemos/bagging.py`](https://dictionaryofml.org/terms/bagging.py); CC BY 4.0."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "# Notebook shim: the script resolves output paths relative to __file__,\n# which a notebook kernel does not define; everything lands in the\n# working directory instead.\nimport os\n__file__ = os.path.join(os.getcwd(), \"bagging.py\")\nos.makedirs(\"pythondemos\", exist_ok=True)"
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "\"\"\"\nbagging.py \u2014 numerical companion to the glossary entry 'bootstrap\naggregating (bagging)'.\n\nPurpose\n-------\nBacks the entry's claims with numbers: a decision tree grown until it\nreproduces every label changes its predictions noticeably when a handful\nof training points are replaced, while linear regression barely\ndoes; bagging reduces the variance of the tree's predictions and leaves\nits bias unchanged; it does little for linear regression; and adding\nmore base learners does not make the aggregated hypothesis fit the\ntraining set more closely.  Self-contained (numpy/matplotlib only),\nfixed seed.\n\nSetup\n-----\nRegression with one feature x in [0, 1] and label\ny = sin(2 pi x) + 0.3 x + noise (noise variance 0.09).  A training set\nholds 50 points.  Base learners: a decision tree grown until every\nleaf holds a single training point (so it reproduces every label) and linear\nregression (a straight line).  Bagging draws B bootstrap resamples of\nthe training set, trains one base learner on each, and averages their\npredictions.  Variance and bias are measured on a grid of 200 feature\nvalues over 100 independent training sets.\n\nBlocks\n------\n[B-unstable] Replacing 5 of the 50 training points, repeated 50 times,\n             changes the tree's predictions by more than five times as\n             much as linear regression's, measured by the largest change\n             over the grid (the tree's change is local but large, the\n             line's global but small).\n[B-variance] Over 100 training sets, bagging with B = 25 lowers the\n             variance of the tree's predictions by a factor above 2\n             while its squared bias stays below 0.01 (as does the single\n             tree's);\n             for linear regression the variance drops by less than a\n             factor 1.3.\n[B-nofit]    The average loss of the bagged tree on its own training\n             set does not fall toward zero as B grows: at B = 100 it is\n             still above half of its value at B = 10, unlike for a\n             method that accumulates its base learners.\n\nOutputs\n-------\nbagging_variance.csv : per base learner: variance and squared bias of\n                       the single and of the bagged (B = 25) predictions,\n                       averaged over the grid.\nbagging_trainloss.csv : B and the average loss of the bagged tree on\n                       the training set.\nbagging.png           : matplotlib preview (checking only).\n\"\"\"\n\nimport numpy as np\nimport matplotlib\n\nmatplotlib.use(\"Agg\")\nimport matplotlib.pyplot as plt\n\nfrom pathlib import Path\n\nOUT_DIR = Path(__file__).parent\n\nreport = []\n\n\ndef check(name, ok):\n    report.append((name, bool(ok)))\n    print(f\"  [{'ok' if ok else 'FAIL'}] {name}\")\n\n\nrng = np.random.default_rng(0)\nm = 50\ngrid = np.linspace(0.0, 1.0, 200)\n\n\ndef truth(x):\n    return np.sin(2 * np.pi * x) + 0.3 * x\n\n\ndef draw_trainset(size=m):\n    x = np.sort(rng.random(size))\n    return x, truth(x) + 0.3 * rng.standard_normal(size)\n\n\n# ------------------------------------------------------- base learners\ndef tree_fit(x, y):\n    \"\"\"Decision tree grown until every leaf holds one training point: a list\n    of (threshold, value) leaves, i.e. a piecewise-constant map.\"\"\"\n    order = np.argsort(x); x, y = x[order], y[order]\n    cuts = []\n\n    def split(lo, hi):                       # indices lo..hi-1\n        if hi - lo <= 1:\n            return\n        best, pos = None, None\n        for k in range(lo + 1, hi):\n            if x[k] == x[k - 1]:\n                continue\n            l, r = y[lo:k], y[k:hi]\n            sse = ((l - l.mean()) ** 2).sum() + ((r - r.mean()) ** 2).sum()\n            if best is None or sse < best:\n                best, pos = sse, k\n        if pos is None:\n            return\n        cuts.append(0.5 * (x[pos] + x[pos - 1]))\n        split(lo, pos); split(pos, hi)\n    split(0, len(x))\n    cuts = np.sort(np.array(cuts))\n    bins = np.searchsorted(cuts, x, side=\"right\")\n    values = np.array([y[bins == b].mean() for b in range(len(cuts) + 1)])\n    return lambda z: values[np.searchsorted(cuts, z, side=\"right\")]\n\n\ndef linreg_fit(x, y):\n    w = np.polyfit(x, y, 1)\n    return lambda z: np.polyval(w, z)\n\n\ndef bagged(fit, x, y, B):\n    hs = []\n    for _ in range(B):\n        idx = rng.integers(0, len(x), len(x))             # bootstrap resample\n        hs.append(fit(x[idx], y[idx]))\n    return lambda z: np.mean([h(z) for h in hs], axis=0)"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "**[B-unstable]** Replacing 5 of the 50 training points, repeated 50 times, changes the tree's predictions by more than five times as much as linear regression's, measured by the largest change over the grid (the tree's change is local but large, the line's global but small)."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "chg_tree, chg_lin = [], []\nfor _ in range(50):\n    x0, y0 = draw_trainset()\n    x1, y1 = x0.copy(), y0.copy()\n    swap = rng.choice(m, 5, replace=False)                  # a handful\n    x1[swap] = rng.random(5); y1[swap] = truth(x1[swap]) + 0.3 * rng.standard_normal(5)\n    chg_tree.append(np.abs(tree_fit(x0, y0)(grid) - tree_fit(x1, y1)(grid)).max())\n    chg_lin.append(np.abs(linreg_fit(x0, y0)(grid) - linreg_fit(x1, y1)(grid)).max())\nchg_tree, chg_lin = float(np.mean(chg_tree)), float(np.mean(chg_lin))\ncheck(f\"[B-unstable] replacing 5 of 50 training points changes the tree's predictions \"\n      f\"by up to {chg_tree:.2f} (average over 50 repetitions), linear regression's \"\n      f\"by up to {chg_lin:.2f}\", chg_tree > 5 * chg_lin)"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "**[B-variance]** Over 100 training sets, bagging with B = 25 lowers the variance of the tree's predictions by a factor above 2 while its squared bias stays below 0.01 (as does the single tree's); for linear regression the variance drops by less than a factor 1.3."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "K, B = 100, 25\nres = {}\nfor name, fit in ((\"tree\", tree_fit), (\"linreg\", linreg_fit)):\n    single, bag = [], []\n    for _ in range(K):\n        x, y = draw_trainset()\n        single.append(fit(x, y)(grid)); bag.append(bagged(fit, x, y, B)(grid))\n    single, bag = np.array(single), np.array(bag)\n    res[name] = dict(\n        var_single=float(single.var(axis=0).mean()), var_bag=float(bag.var(axis=0).mean()),\n        bias_single=float(((single.mean(axis=0) - truth(grid)) ** 2).mean()),\n        bias_bag=float(((bag.mean(axis=0) - truth(grid)) ** 2).mean()))\nt, l = res[\"tree\"], res[\"linreg\"]\ncheck(f\"[B-variance] tree: variance {t['var_single']:.3f} -> {t['var_bag']:.3f} \"\n      f\"(factor {t['var_single'] / t['var_bag']:.1f}), squared bias \"\n      f\"{t['bias_single']:.3f} -> {t['bias_bag']:.3f}; linear regression: variance \"\n      f\"{l['var_single']:.4f} -> {l['var_bag']:.4f} (factor {l['var_single'] / l['var_bag']:.2f})\",\n      t[\"var_single\"] > 2 * t[\"var_bag\"]\n      and max(t[\"bias_bag\"], t[\"bias_single\"]) < 0.01\n      and l[\"var_single\"] < 1.3 * l[\"var_bag\"])"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "**[B-nofit]** The average loss of the bagged tree on its own training set does not fall toward zero as B grows: at B = 100 it is still above half of its value at B = 10, unlike for a method that accumulates its base learners."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "x, y = draw_trainset()\nBs = [1, 5, 10, 25, 50, 100]\ntrain_loss = [float(np.mean((y - bagged(tree_fit, x, y, b)(x)) ** 2)) for b in Bs]\ncheck(\"[B-nofit]    average loss of the bagged tree on its training set: \"\n      + \", \".join(f\"B={b}: {v:.3f}\" for b, v in zip(Bs, train_loss)),\n      train_loss[-1] > 0.5 * train_loss[2])\n\n# ---------------------------------------------------------------- CSV\nwith open(OUT_DIR / \"bagging_variance.csv\", \"w\") as fh:\n    fh.write(\"learner,var_single,var_bagged,bias2_single,bias2_bagged\\n\")\n    for name in (\"tree\", \"linreg\"):\n        r = res[name]\n        fh.write(f\"{name},{r['var_single']:.4f},{r['var_bag']:.4f},\"\n                 f\"{r['bias_single']:.4f},{r['bias_bag']:.4f}\\n\")\nwith open(OUT_DIR / \"bagging_trainloss.csv\", \"w\") as fh:\n    fh.write(\"B,train_loss\\n\")\n    for b, v in zip(Bs, train_loss):\n        fh.write(f\"{b},{v:.4f}\\n\")\n\n# -------------------------------------------------------------- preview\nfrom matplotlib.patches import Patch\nfig, (ax, ax2) = plt.subplots(1, 2, figsize=(9.4, 3.8))\npos = np.array([0, 1, 2.5, 3.5])\nvals = [t[\"var_single\"], t[\"var_bag\"], l[\"var_single\"], l[\"var_bag\"]]\nbars = ax.bar(pos, vals, 0.8, color=[\"0.75\", \"0.3\"] * 2, edgecolor=\"k\")\nfor b in bars[1::2]:\n    b.set_hatch(\"//\")\nax.plot(pos[:2], [t[\"bias_single\"], t[\"bias_bag\"]], \"k_\", ms=18, mew=2)\nax.plot(pos[2:], [l[\"bias_single\"], l[\"bias_bag\"]], \"k_\", ms=18, mew=2)\nax.set_xticks([0.5, 3.0]); ax.set_xticklabels([\"decision tree\", \"linear regression\"])\nax.set_xlabel(\"base learner\")\nax.set_ylabel(\"variance (bars), squared bias (dashes)\")\nax.set_title(\"bagging with B = 25 over 100 training sets\")\nax.legend(handles=[Patch(facecolor=\"0.75\", edgecolor=\"k\", label=\"single base learner\"),\n                   Patch(facecolor=\"0.3\", edgecolor=\"k\", hatch=\"//\", label=\"bagged\")],\n          frameon=False, fontsize=8)\nax2.plot(Bs, train_loss, \"k.-\", lw=1.2, ms=7, label=\"bagged decision tree\")\nax2.set_xscale(\"log\")\nax2.set_xlabel(\"number of base learners B\")\nax2.set_ylabel(\"average loss on the training set\")\nax2.set_title(\"the training-set loss does not fall with B\")\nax2.set_ylim(0, max(train_loss) * 1.2)\nax2.legend(frameon=False, fontsize=8)\nfig.tight_layout()\nfig.savefig(OUT_DIR / \"bagging.png\", dpi=110)\n\nn_ok = sum(ok for _, ok in report)\nprint(f\"\\n{n_ok}/{len(report)} checks pass\")\nprint(f\"wrote {OUT_DIR / 'bagging_variance.csv'}, {OUT_DIR / 'bagging_trainloss.csv'}, \"\n      f\"{OUT_DIR / 'bagging.png'}\")\nif n_ok != len(report):\n    raise SystemExit(1)"
  }
 ]
}