{
"nbformat": 4,
"nbformat_minor": 5,
"metadata": {
"kernelspec": {
"name": "python3",
"display_name": "Python 3",
"language": "python"
},
"language_info": {
"name": "python"
},
"colab": {
"name": "stochGD.ipynb"
}
},
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": "# stochastic gradient descent (SGD) \u2014 Python demo\n\nNumerical companion to the entry [stochastic gradient descent (SGD)](https://dictionaryofml.org/terms/stochGD.html) of the [Dictionary of Applied Machine Learning](https://dictionaryofml.org/): it recomputes what the entry states and prints one line per check.\n\nOne block per claim/figure of the entry (marked [B-...], linked from the entry's paragraphs by content): each block verifies numerically what the corresponding statement asserts or generates the data behind a figure. Fixed seeds; numpy/matplotlib plus stdlib urllib for the one-time download of the weather observations behind Fig. 1 (cached in the committed CSVs \u2014 re-runs read those and need no network).\n\nRequires NumPy and Matplotlib only, and uses fixed seeds, so the printed numbers reproduce exactly. Generated from [`pythondemos/stochGD.py`](https://dictionaryofml.org/terms/stochGD.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(), \"stochGD.py\")\nos.makedirs(\"pythondemos\", exist_ok=True)"
},
{
"cell_type": "code",
"metadata": {},
"execution_count": null,
"outputs": [],
"source": "\"\"\"\nstochGD.py \u2014 numerical companion to the glossary entry\n'stochastic gradient descent (SGD)'.\n\nOne block per claim/figure of the entry (marked [B-...], linked from the\nentry's paragraphs by content): each block verifies numerically what the\ncorresponding statement asserts or generates the data behind a figure.\nFixed seeds; numpy/matplotlib plus stdlib urllib for the one-time\ndownload of the weather observations behind Fig. 1 (cached in the\ncommitted CSVs \u2014 re-runs read those and need no network).\n\nBlocks\n------\n[B-weather] Real training set behind Fig. 1 of the entry: daily\n minimum temperature (feature x) and daily maximum\n temperature (label y) observed on 2026-07-01 at 6 Austrian\n weather stations (downloaded from the GeoSphere Austria\n open data API, dataset klima-v2-1d) and 6 Finnish stations\n (downloaded from the FMI open data WFS, daily\n observations). A linear hypothesis map h(x) = w1 x + w0 is\n learned by ERM on the m = 12 data points; a fixed random\n subset of 4 data points is the batch that one SGD update\n uses in place of the full training set.\n[B-approx] SGD replaces the full gradient \u2014 a sum of per-data-point\n gradients over the entire trainset \u2014 with the sum over a\n randomly chosen batch: the batch gradient is an unbiased\n approximation (its average over many random batches\n matches the full gradient), and one SGD step touches only\n |B| of the m data points.\n[B-batchsize] The batch size trades gradient accuracy against cost: the\n approximation error of the batch gradient shrinks like\n 1/sqrt(|B|) as the batch grows (variance scaling), while\n the per-step cost grows linearly in |B|.\n[B-minibatch] Mini-batch SGD (|B| > 1) converges to the ERM solution on\n a least-squares problem while evaluating only a small\n fraction of the per-data-point gradients that full GD\n uses for the same number of passes.\n\nOutputs\n-------\nstochGD_geosphere.csv : tmin, tmax \u2014 the 6 Austrian data points.\nstochGD_fmi.csv : tmin, tmax \u2014 the 6 Finnish data points.\nstochGD_batch.csv : tmin, tmax \u2014 the 4 data points of the batch.\nstochGD_hypothesis.csv : tmin, tmax \u2014 endpoints of the learned h(x).\nstochGD_errorsegs.csv : x, y \u2014 prediction-error segments (one segment\n per data point, nan-separated for pgfplots\n 'unbounded coords=jump').\nstochGD.png : preview figure (checking only).\n\nData generated by pythondemos/stochGD.py.\n\"\"\"\n\nimport numpy as np\nimport matplotlib\n\nmatplotlib.use(\"Agg\")\nimport matplotlib.pyplot as plt\nfrom pathlib import Path\n\nOUT_DIR = Path(__file__).parent\n\nrng = np.random.default_rng(42)\nreport = []\n\n\ndef check(name, ok):\n report.append((name, bool(ok)))\n print(f\" [{'ok' if ok else 'FAIL'}] {name}\")\n\n\n# least-squares ERM objective on m data points\nm, d = 2000, 5\nX = rng.normal(size=(m, d))\nw_true = rng.normal(size=d)\ny = X @ w_true + 0.3 * rng.normal(size=m)\nfull_grad = lambda w: (2 / m) * X.T @ (X @ w - y)\ndef batch_grad(w, B):\n idx = rng.choice(m, B, replace=False)\n return (2 / B) * X[idx].T @ (X[idx] @ w - y[idx]), B"
},
{
"cell_type": "markdown",
"metadata": {},
"source": "**[B-weather]** Real training set behind Fig. 1 of the entry: daily minimum temperature (feature x) and daily maximum temperature (label y) observed on 2026-07-01 at 6 Austrian weather stations (downloaded from the GeoSphere Austria open data API, dataset klima-v2-1d) and 6 Finnish stations (downloaded from the FMI open data WFS, daily observations). A linear hypothesis map h(x) = w1 x + w0 is learned by ERM on the m = 12 data points; a fixed random subset of 4 data points is the batch that one SGD update uses in place of the full training set."
},
{
"cell_type": "code",
"metadata": {},
"execution_count": null,
"outputs": [],
"source": "print(\"[B-weather] real tmin/tmax data points behind Fig. 1 of the entry\")\n\nDATE = \"2026-07-01\"\nGEO_STATIONS = {105: \"Wien Hohe Warte\", 30: \"Graz Universitaet\",\n 39: \"Innsbruck Universitaet\", 48: \"Klagenfurt Flughafen\",\n 131: \"Salzburg Flughafen\", 5000: \"Linz Hoersching\"}\nFMI_PLACES = [\"Helsinki\", \"Tampere\", \"Turku\", \"Oulu\", \"Rovaniemi\",\n \"Kuopio\"]\n\n\ndef download_geosphere():\n \"\"\"One (tmin, tmax) pair per Austrian station for DATE.\"\"\"\n import json\n from urllib.request import urlopen\n ids = \",\".join(str(s) for s in GEO_STATIONS)\n url = (\"https://dataset.api.hub.geosphere.at/v1/station/historical/\"\n f\"klima-v2-1d?parameters=tlmin,tlmax&station_ids={ids}\"\n f\"&start={DATE}&end={DATE}\")\n feats = json.load(urlopen(url, timeout=60))[\"features\"]\n return sorted(\n (f[\"properties\"][\"parameters\"][\"tlmin\"][\"data\"][0],\n f[\"properties\"][\"parameters\"][\"tlmax\"][\"data\"][0])\n for f in feats)\n\n\ndef download_fmi():\n \"\"\"One (tmin, tmax) pair per Finnish place for DATE.\"\"\"\n import re\n from urllib.request import urlopen\n pairs = []\n for place in FMI_PLACES:\n url = (\"https://opendata.fmi.fi/wfs?service=WFS&version=2.0.0\"\n \"&request=getFeature&storedquery_id=\"\n \"fmi::observations::weather::daily::simple\"\n f\"&place={place}&starttime={DATE}T00:00:00Z\"\n f\"&endtime={DATE}T23:59:59Z¶meters=tmin,tmax\")\n xml = urlopen(url, timeout=60).read().decode()\n vals = dict(re.findall(\n r\"(\\w+)\\s*\"\n r\"([-\\d.]+)\",\n xml))\n pairs.append((float(vals[\"tmin\"]), float(vals[\"tmax\"])))\n return sorted(pairs)\n\n\ndef write_csv(name, header, rows):\n with open(OUT_DIR / name, \"w\") as f:\n f.write(header + \"\\n\")\n for row in rows:\n f.write(\",\".join(row) + \"\\n\")\n\n\ndef read_csv(name):\n lines = (OUT_DIR / name).read_text().strip().splitlines()[1:]\n return [tuple(float(v) for v in ln.split(\",\")) for ln in lines]\n\n\ntry: # one-time download ...\n geo, fmi = download_geosphere(), download_fmi()\n write_csv(\"stochGD_geosphere.csv\", \"tmin,tmax\",\n [(f\"{a:.1f}\", f\"{b:.1f}\") for a, b in geo])\n write_csv(\"stochGD_fmi.csv\", \"tmin,tmax\",\n [(f\"{a:.1f}\", f\"{b:.1f}\") for a, b in fmi])\n print(f\" downloaded {len(geo)}+{len(fmi)} stations for {DATE}\")\nexcept OSError as err: # ... offline: reuse committed CSVs\n print(f\" download failed ({err}); reading committed CSVs\")\n geo = read_csv(\"stochGD_geosphere.csv\")\n fmi = read_csv(\"stochGD_fmi.csv\")\n\npts = geo + fmi # the training set, m = 12\nx_all = np.array([p[0] for p in pts])\ny_all = np.array([p[1] for p in pts])\n\n# linear hypothesis h(x) = w1 x + w0 learned by ERM (average squared\n# prediction error) on all m data points\nfit_w1, fit_w0 = np.polyfit(x_all, y_all, 1)\nyhat = fit_w1 * x_all + fit_w0\nxg = (np.floor(x_all.min()) - 1.0, np.ceil(x_all.max()) + 1.0)\nwrite_csv(\"stochGD_hypothesis.csv\", \"tmin,tmax\",\n [(f\"{x:.1f}\", f\"{fit_w1 * x + fit_w0:.2f}\") for x in xg])\n\n# prediction-error segments, nan-separated for pgfplots\nseg_rows = []\nfor xi, yi, hi in zip(x_all, y_all, yhat):\n seg_rows += [(f\"{xi:.1f}\", f\"{yi:.1f}\"), (f\"{xi:.1f}\", f\"{hi:.2f}\"),\n (\"nan\", \"nan\")]\nwrite_csv(\"stochGD_errorsegs.csv\", \"x,y\", seg_rows)\n\n# the batch: a fixed random subset of 4 of the m = 12 data points\nbatch_idx = np.sort(np.random.default_rng(1).choice(len(pts), 4,\n replace=False))\nwrite_csv(\"stochGD_batch.csv\", \"tmin,tmax\",\n [(f\"{x_all[i]:.1f}\", f\"{y_all[i]:.1f}\") for i in batch_idx])\n\ncheck(\"6 Austrian + 6 Finnish data points (tmin, tmax) for \" + DATE,\n len(geo) == 6 and len(fmi) == 6)\ncheck(\"the batch is a strict random subset of the training set\",\n 0 < len(batch_idx) < len(pts))\ncheck(\"the batch holds data points from both servers\",\n batch_idx.min() < len(geo) <= batch_idx.max())\ncheck(\"ERM residuals of the linear hypothesis average to zero\",\n abs(np.mean(y_all - yhat)) < 1e-8)"
},
{
"cell_type": "markdown",
"metadata": {},
"source": "**[B-approx]** SGD replaces the full gradient \u2014 a sum of per-data-point gradients over the entire trainset \u2014 with the sum over a randomly chosen batch: the batch gradient is an unbiased approximation (its average over many random batches matches the full gradient), and one SGD step touches only |B| of the m data points."
},
{
"cell_type": "code",
"metadata": {},
"execution_count": null,
"outputs": [],
"source": "print(\"[B-approx] the batch gradient approximates the full-sum gradient\")\nw0 = np.zeros(d)\ng_full = full_grad(w0)\ng_avg = np.mean([batch_grad(w0, 20)[0] for _ in range(4000)], axis=0)\ncheck(\"averaging batch gradients over many draws recovers the full \"\n \"gradient (unbiasedness)\",\n np.linalg.norm(g_avg - g_full) < 0.05 * np.linalg.norm(g_full))\ncheck(\"one SGD step touches |B| = 20 of the m = 2000 data points\",\n batch_grad(w0, 20)[1] == 20 < m)"
},
{
"cell_type": "markdown",
"metadata": {},
"source": "**[B-batchsize]** The batch size trades gradient accuracy against cost: the approximation error of the batch gradient shrinks like 1/sqrt(|B|) as the batch grows (variance scaling), while the per-step cost grows linearly in |B|."
},
{
"cell_type": "code",
"metadata": {},
"execution_count": null,
"outputs": [],
"source": "print(\"[B-batchsize] batch size trades accuracy against cost\")\nerrs = []\nfor B in (10, 100, 1000):\n errs.append(np.mean([np.linalg.norm(batch_grad(w0, B)[0] - g_full)\n for _ in range(300)]))\nprint(f\" mean gradient error at |B| = 10, 100, 1000: \"\n f\"{errs[0]:.3f}, {errs[1]:.3f}, {errs[2]:.3f}\")\ncheck(\"the gradient error shrinks as the batch grows\",\n errs[0] > errs[1] > errs[2])\ncheck(\"error scaling is consistent with 1/sqrt(|B|) \"\n \"(10x batch -> ~3.2x smaller error)\",\n 2.0 < errs[0] / errs[1] < 5.0 and 2.0 < errs[1] / errs[2] < 5.0)"
},
{
"cell_type": "markdown",
"metadata": {},
"source": "**[B-minibatch]** Mini-batch SGD (|B| > 1) converges to the ERM solution on a least-squares problem while evaluating only a small fraction of the per-data-point gradients that full GD uses for the same number of passes."
},
{
"cell_type": "code",
"metadata": {},
"execution_count": null,
"outputs": [],
"source": "print(\"[B-minibatch] mini-batch SGD reaches the ERM solution cheaply\")\nw_hat = np.linalg.solve(X.T @ X, X.T @ y) # ERM solution\nw = np.zeros(d)\ngrads_evaluated = 0\nfor t in range(1, 1201):\n g, B = batch_grad(w, 20)\n w -= (0.05 / np.sqrt(t)) * g\n grads_evaluated += B\ngd_grads = 1200 * m # full GD, same steps\ncheck(\"mini-batch SGD converges near the ERM solution\",\n np.linalg.norm(w - w_hat) < 0.1)\ncheck(\"using 1% of the per-data-point gradient evaluations of full GD\",\n grads_evaluated == 0.01 * gd_grads)\n\n# ------------------------------------------------------------ preview\nfig, (ax0, ax1) = plt.subplots(1, 2, figsize=(9.6, 3.4))\ngx = np.array(xg)\nax0.plot(gx, fit_w1 * gx + fit_w0, \"-\", color=\"tab:blue\",\n label=\"hypothesis h(x)\")\nfor xi, yi, hi in zip(x_all, y_all, yhat):\n ax0.plot([xi, xi], [yi, hi], \":\", color=\"gray\", lw=1)\nax0.plot([p[0] for p in geo], [p[1] for p in geo], \"o\", color=\"black\",\n label=\"stored at Geosphere.at\")\nax0.plot([p[0] for p in fmi], [p[1] for p in fmi], \"s\", color=\"black\",\n label=\"stored at FMI.fi\")\nax0.plot(x_all[batch_idx], y_all[batch_idx], \"o\", mfc=\"none\",\n mec=\"tab:red\", ms=12, label=\"batch\")\nax0.set_xlabel(\"daily minimum temperature (deg C)\")\nax0.set_ylabel(\"daily maximum temperature (deg C)\")\nax0.set_title(f\"[B-weather] trainset, batch and errors ({DATE})\")\nax0.legend(frameon=False, fontsize=8)\nax1.loglog([10, 100, 1000], errs, \"o-\")\nax1.set_xlabel(\"batch size |B|\"); ax1.set_ylabel(\"gradient error\")\nax1.set_title(\"[B-batchsize] accuracy vs batch size\")\nfig.tight_layout()\nfig.savefig(OUT_DIR / \"stochGD.png\", dpi=110)\nprint(f\"\\n{sum(ok for _, ok in report)}/{len(report)} checks passed\")\nassert all(ok for _, ok in report)"
}
]
}