{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "150c154f",
   "metadata": {},
   "source": [
    "<link rel=\"stylesheet\" href=\"berkeley.css\">\n",
    "\n",
    "<h1 class=\"cal cal-h1\">Lecture 08: Bias-Variance and Logistic Regression (1) &ndash; CS 189, Fall 2026</h1>\n",
    "\n",
    "*Reference: Bishop & Bishop,* Deep Learning: Foundations and Concepts, *&sect;4.3 (Part 0), Ch. 5 &sect;5.1 and &sect;5.4.1-5.4.3 (Parts 1-6).*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "71cc02d4",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:29.041361Z",
     "iopub.status.busy": "2026-09-22T21:10:29.041168Z",
     "iopub.status.idle": "2026-09-22T21:10:32.673953Z",
     "shell.execute_reply": "2026-09-22T21:10:32.671908Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "# Part 0 (bias-variance)\n",
    "import matplotlib.pyplot as plt\n",
    "from numpy.polynomial import legendre as L\n",
    "\n",
    "# Parts 1-6 (logistic regression)\n",
    "import plotly.express as px\n",
    "import plotly.graph_objects as go\n",
    "from sklearn.datasets import load_breast_cancer\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.linear_model import LinearRegression, LogisticRegression\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from scipy.optimize import minimize\n",
    "\n",
    "import warnings\n",
    "from sklearn.exceptions import ConvergenceWarning\n",
    "warnings.filterwarnings(\"ignore\", category=ConvergenceWarning)\n",
    "warnings.filterwarnings(\"ignore\", category=RuntimeWarning)\n",
    "warnings.filterwarnings(\"ignore\", message=\"Setting penalty=None\")\n",
    "\n",
    "import plotly.io as pio\n",
    "pio.renderers.default = \"notebook_connected\"\n",
    "\n",
    "C_DATA, C_FIT, C_ALT, C_RESID, C_SPAN = \"#003262\", \"#FDB515\", \"#C4820E\", \"#D55E00\", \"#00553A\"\n",
    "plt.rcParams.update({\"figure.figsize\": (8, 5), \"axes.grid\": True, \"grid.alpha\": 0.3,\n",
    "                     \"font.size\": 13, \"axes.titlesize\": 15, \"axes.labelsize\": 14})"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a41c3eb5",
   "metadata": {},
   "source": [
    "### Part 0. **The Bias-Variance Decomposition**  *(carried over from Lecture 07)*\n",
    "\n",
    "Take the expected squared error at an input $x$, averaged over training sets $\\mathcal{D}$:\n",
    "\n",
    "$$\n",
    "\\mathbb{E}_{\\mathcal{D}}\\big[(t - y(x;\\mathcal{D}))^2\\big]\n",
    "= \\underbrace{\\sigma^2}_{\\text{noise}}\n",
    "+ \\underbrace{\\big(\\mathbb{E}_{\\mathcal{D}}[y(x;\\mathcal{D})] - f(x)\\big)^2}_{\\text{bias}^2}\n",
    "+ \\underbrace{\\mathbb{E}_{\\mathcal{D}}\\big[(y(x;\\mathcal{D}) - \\mathbb{E}_{\\mathcal{D}}[y(x;\\mathcal{D})])^2\\big]}_{\\text{variance}}\n",
    "$$\n",
    "\n",
    "- **Noise** is irreducible. No model removes it.\n",
    "- **Bias** is how far the *average* fitted function sits from the truth: underfitting.\n",
    "- **Variance** is how much the fitted function moves when the training set changes: overfitting.\n",
    "\n",
    "Every term is an expectation over training sets, so we can estimate all three by simulation: draw many training sets, fit each, and look at the spread of the fitted functions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "70592e7d",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:32.676475Z",
     "iopub.status.busy": "2026-09-22T21:10:32.676247Z",
     "iopub.status.idle": "2026-09-22T21:10:33.155375Z",
     "shell.execute_reply": "2026-09-22T21:10:33.154338Z"
    }
   },
   "outputs": [],
   "source": [
    "SIGMA_BV, NSETS = 0.30, 400\n",
    "f_true = lambda x: np.sin(5 * x)\n",
    "xg_bv = np.linspace(-0.9, 0.9, 120)\n",
    "\n",
    "def Phi_leg(x, deg):\n",
    "    \"Legendre basis: well conditioned, so what we measure is variance and not rounding error.\"\n",
    "    return L.legvander(np.asarray(x).ravel(), deg)\n",
    "\n",
    "def decompose(n_train, degrees, seed=189):\n",
    "    rng = np.random.default_rng(seed)\n",
    "    out = []\n",
    "    for d in degrees:\n",
    "        preds = np.empty((NSETS, len(xg_bv)))\n",
    "        for s in range(NSETS):\n",
    "            xt = rng.uniform(-0.9, 0.9, n_train)\n",
    "            tt = f_true(xt) + SIGMA_BV * rng.standard_normal(n_train)\n",
    "            w = np.linalg.lstsq(Phi_leg(xt, d), tt, rcond=None)[0]\n",
    "            preds[s] = Phi_leg(xg_bv, d) @ w\n",
    "        bias2 = float(np.mean((preds.mean(0) - f_true(xg_bv)) ** 2))\n",
    "        var = float(np.mean(preds.var(0)))\n",
    "        out.append((d, bias2, var, SIGMA_BV ** 2, bias2 + var + SIGMA_BV ** 2, preds))\n",
    "    return out\n",
    "\n",
    "degrees = list(range(0, 11))\n",
    "res100 = decompose(100, degrees)\n",
    "\n",
    "print(f\"N_train = 100,  {NSETS} training sets per degree\\n\")\n",
    "print(f\"{'deg':>4} {'bias^2':>9} {'variance':>9} {'noise':>8} {'total':>9}\")\n",
    "for d, b2, v, nz, tot, _ in res100:\n",
    "    print(f\"{d:>4} {b2:>9.4f} {v:>9.4f} {nz:>8.4f} {tot:>9.4f}\")\n",
    "best = min(res100, key=lambda r: r[4])\n",
    "print(f\"\\nlowest expected error at degree {best[0]}  (total {best[4]:.4f})\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ef037752",
   "metadata": {},
   "source": [
    "#### 0.1 Seeing bias and variance directly\n",
    "\n",
    "Plot the fitted function from 40 of the training sets, plus their average, against the truth. Bias is the gap between the average (dashed) and the truth (solid). Variance is the width of the grey band."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6735ed83",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:33.158169Z",
     "iopub.status.busy": "2026-09-22T21:10:33.157954Z",
     "iopub.status.idle": "2026-09-22T21:10:33.752238Z",
     "shell.execute_reply": "2026-09-22T21:10:33.751115Z"
    }
   },
   "outputs": [],
   "source": [
    "show_degrees = [1, 3, 5, 9]\n",
    "fig, axs = plt.subplots(1, 4, figsize=(17, 4.2), sharey=True)\n",
    "for ax, d in zip(axs, show_degrees):\n",
    "    rec = [r for r in res100 if r[0] == d][0]\n",
    "    preds = rec[5]\n",
    "    for p in preds[:40]:\n",
    "        ax.plot(xg_bv, p, color=\"0.55\", lw=0.8, alpha=0.55)\n",
    "    ax.plot(xg_bv, f_true(xg_bv), lw=3, color=C_DATA, label=\"truth $f(x)$\")\n",
    "    ax.plot(xg_bv, preds.mean(0), lw=3, ls=\"--\", color=C_ALT, label=\"average fit\")\n",
    "    ax.set_title(f\"degree {d}\\nbias$^2$ = {rec[1]:.3f}\\nvariance = {rec[2]:.3f}\", fontsize=13)\n",
    "    ax.set(xlabel=\"$x$\", ylim=(-2.0, 2.0))\n",
    "axs[0].set_ylabel(\"$y$\"); axs[0].legend(fontsize=10, loc=\"lower left\")\n",
    "plt.tight_layout(); plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "55e55b2c",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:33.754328Z",
     "iopub.status.busy": "2026-09-22T21:10:33.754082Z",
     "iopub.status.idle": "2026-09-22T21:10:34.027169Z",
     "shell.execute_reply": "2026-09-22T21:10:34.026082Z"
    }
   },
   "outputs": [],
   "source": [
    "b2 = np.array([r[1] for r in res100])\n",
    "vr = np.array([r[2] for r in res100])\n",
    "nz = np.array([r[3] for r in res100])\n",
    "tot = b2 + vr + nz\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(8.5, 5))\n",
    "ax.plot(degrees, b2, \"o-\", lw=3, ms=8, color=C_DATA, label=\"bias$^2$\")\n",
    "ax.plot(degrees, vr, \"s-\", lw=3, ms=8, color=C_RESID, label=\"variance\")\n",
    "ax.plot(degrees, nz, \":\", lw=2.5, color=\"0.4\", label=\"noise $\\\\sigma^2$\")\n",
    "ax.plot(degrees, tot, \"^-\", lw=3, ms=8, color=C_SPAN, label=\"total\")\n",
    "ax.axvline(degrees[int(tot.argmin())], color=C_ALT, ls=\"--\", lw=2.5,\n",
    "           label=f\"minimum at degree {degrees[int(tot.argmin())]}\")\n",
    "ax.set(xlabel=\"polynomial degree\", ylabel=\"expected squared error\", yscale=\"log\",\n",
    "       xticks=degrees)\n",
    "ax.legend(fontsize=11); plt.tight_layout(); plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c2dd23f4",
   "metadata": {},
   "source": [
    "#### 0.2 More data supports more complexity\n",
    "\n",
    "Nothing about the tradeoff is a property of the model alone. Shrink the training set and the variance term grows at every degree, so the minimum moves left: with less data the best model is *simpler* than the truth."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "98ba74f8",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:34.029188Z",
     "iopub.status.busy": "2026-09-22T21:10:34.028987Z",
     "iopub.status.idle": "2026-09-22T21:10:34.783417Z",
     "shell.execute_reply": "2026-09-22T21:10:34.782512Z"
    }
   },
   "outputs": [],
   "source": [
    "res25 = decompose(25, degrees)\n",
    "tot25 = np.array([r[4] for r in res25])\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(8.5, 5))\n",
    "for res, tot_i, col, name in [(res25, tot25, C_RESID, \"N = 25\"),\n",
    "                              (res100, tot, C_DATA, \"N = 100\")]:\n",
    "    ax.plot(degrees, tot_i, \"o-\", lw=3, ms=8, color=col, label=f\"{name}, total error\")\n",
    "    k = int(np.argmin(tot_i))\n",
    "    ax.plot(degrees[k], tot_i[k], \"*\", ms=22, color=C_ALT, mec=\"k\", zorder=5)\n",
    "ax.set(xlabel=\"polynomial degree\", ylabel=\"expected squared error\", yscale=\"log\",\n",
    "       xticks=degrees)\n",
    "ax.legend(fontsize=12); plt.tight_layout(); plt.show()\n",
    "\n",
    "print(f\"N = 25  -> minimum at degree {degrees[int(tot25.argmin())]} (total {tot25.min():.4f})\")\n",
    "print(f\"N = 100 -> minimum at degree {degrees[int(tot.argmin())]} (total {tot.min():.4f})\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8a8edd95",
   "metadata": {},
   "source": [
    "#### 0.3 Regularization moves along the same curve\n",
    "\n",
    "Model complexity is not only the number of basis functions. Hold the degree fixed at 9 and turn $\\lambda$ instead: large $\\lambda$ buys bias and sells variance, exactly as reducing the degree does."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "adcfe88a",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:34.785933Z",
     "iopub.status.busy": "2026-09-22T21:10:34.785664Z",
     "iopub.status.idle": "2026-09-22T21:10:36.342544Z",
     "shell.execute_reply": "2026-09-22T21:10:36.341292Z"
    }
   },
   "outputs": [],
   "source": [
    "def decompose_ridge(n_train, deg, lambdas, seed=189):\n",
    "    rng = np.random.default_rng(seed)\n",
    "    out = []\n",
    "    for lam in lambdas:\n",
    "        preds = np.empty((NSETS, len(xg_bv)))\n",
    "        for s in range(NSETS):\n",
    "            xt = rng.uniform(-0.9, 0.9, n_train)\n",
    "            tt = f_true(xt) + SIGMA_BV * rng.standard_normal(n_train)\n",
    "            A = Phi_leg(xt, deg)\n",
    "            P = np.eye(deg + 1); P[0, 0] = 0.0           # intercept left unpenalized\n",
    "            w = np.linalg.solve(A.T @ A + lam * P, A.T @ tt)\n",
    "            preds[s] = Phi_leg(xg_bv, deg) @ w\n",
    "        bias2 = float(np.mean((preds.mean(0) - f_true(xg_bv)) ** 2))\n",
    "        var = float(np.mean(preds.var(0)))\n",
    "        out.append((lam, bias2, var, bias2 + var + SIGMA_BV ** 2))\n",
    "    return out\n",
    "\n",
    "lams = np.logspace(-4, 3, 22)\n",
    "rr = decompose_ridge(100, 9, lams)\n",
    "b2r = np.array([r[1] for r in rr]); vrr = np.array([r[2] for r in rr])\n",
    "totr = np.array([r[3] for r in rr])\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(8.5, 5))\n",
    "ax.plot(lams, b2r, \"o-\", lw=3, ms=7, color=C_DATA, label=\"bias$^2$\")\n",
    "ax.plot(lams, vrr, \"s-\", lw=3, ms=7, color=C_RESID, label=\"variance\")\n",
    "ax.plot(lams, totr, \"^-\", lw=3, ms=7, color=C_SPAN, label=\"total\")\n",
    "ax.axvline(lams[int(totr.argmin())], color=C_ALT, ls=\"--\", lw=2.5,\n",
    "           label=rf\"best $\\lambda$ = {lams[int(totr.argmin())]:.3g}\")\n",
    "ax.set(xscale=\"log\", yscale=\"log\", xlabel=r\"$\\lambda$  (degree fixed at 9)\",\n",
    "       ylabel=\"expected squared error\")\n",
    "ax.legend(fontsize=11); plt.tight_layout(); plt.show()\n",
    "\n",
    "print(f\"lambda -> 0     : bias^2 {b2r[0]:.4f}, variance {vrr[0]:.4f}\")\n",
    "print(f\"best lambda     : bias^2 {b2r[totr.argmin()]:.4f}, variance {vrr[totr.argmin()]:.4f}\"\n",
    "      f\"  (lambda = {lams[int(totr.argmin())]:.3g})\")\n",
    "print(f\"lambda -> large : bias^2 {b2r[-1]:.4f}, variance {vrr[-1]:.4f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a43e1293",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "Everything above assumed squared error and a real-valued target. Today's second half changes the target\n",
    "to a label in $\\{0, 1\\}$. The first thing to check is whether squared error still makes sense, and the\n",
    "answer will be no: not for the model, and not for the error function."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a95ee264",
   "metadata": {},
   "source": [
    "## The Data: Wisconsin Breast Cancer\n",
    "\n",
    "569 tumor biopsies, 30 measurements each. In scikit-learn, `target = 0` means **malignant**, so we flip it\n",
    "into a `malignant` column where 1 = malignant, 0 = benign."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f53fe48f",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:36.344799Z",
     "iopub.status.busy": "2026-09-22T21:10:36.344599Z",
     "iopub.status.idle": "2026-09-22T21:10:36.369091Z",
     "shell.execute_reply": "2026-09-22T21:10:36.368061Z"
    }
   },
   "outputs": [],
   "source": [
    "data_dict = load_breast_cancer()\n",
    "data = pd.DataFrame(data_dict[\"data\"], columns=data_dict[\"feature_names\"])\n",
    "data[\"malignant\"] = (data_dict[\"target\"] == 0).astype(int)\n",
    "\n",
    "print(\"Fraction malignant:\", data[\"malignant\"].mean().round(3))\n",
    "data[[\"mean radius\", \"mean texture\", \"malignant\"]].head()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "aa5c031f",
   "metadata": {},
   "source": [
    "### One split for the whole unit\n",
    "\n",
    "We make **one** train / validation / test split (60 / 20 / 20) and use it in both logistic regression lectures.\n",
    "Today we only touch train and test. The validation set is for choosing thresholds next lecture."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9568a179",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:36.371407Z",
     "iopub.status.busy": "2026-09-22T21:10:36.371205Z",
     "iopub.status.idle": "2026-09-22T21:10:36.383265Z",
     "shell.execute_reply": "2026-09-22T21:10:36.382149Z"
    }
   },
   "outputs": [],
   "source": [
    "data_tr, data_te = train_test_split(data, test_size=0.2, random_state=42, stratify=data[\"malignant\"])\n",
    "data_tr, data_va = train_test_split(data_tr, test_size=0.25, random_state=42, stratify=data_tr[\"malignant\"])\n",
    "\n",
    "X_train, y_train = data_tr[[\"mean radius\"]].to_numpy(), data_tr[\"malignant\"].to_numpy()\n",
    "X_test,  y_test  = data_te[[\"mean radius\"]].to_numpy(), data_te[\"malignant\"].to_numpy()\n",
    "print(len(data_tr), \"train /\", len(data_va), \"validation /\", len(data_te), \"test\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "404a8cd5",
   "metadata": {},
   "source": [
    "### Plotting helpers\n",
    "\n",
    "Every plot below is \"points colored by class plus some curves\". We write that once.\n",
    "Points are **jittered** vertically so the 0s and 1s do not sit on top of each other."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "226fd6ff",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:36.385263Z",
     "iopub.status.busy": "2026-09-22T21:10:36.385089Z",
     "iopub.status.idle": "2026-09-22T21:10:37.388742Z",
     "shell.execute_reply": "2026-09-22T21:10:37.387893Z"
    }
   },
   "outputs": [],
   "source": [
    "rng = np.random.default_rng(0)\n",
    "\n",
    "def jitter(y, amt=0.1):\n",
    "    return y + amt * rng.random(len(y)) - amt / 2\n",
    "\n",
    "def class_scatter(X, y, title=\"\"):\n",
    "    fig = go.Figure()\n",
    "    for label, name, color in [(0, \"Benign\", \"blue\"), (1, \"Malignant\", \"red\")]:\n",
    "        m = (y == label)\n",
    "        fig.add_scatter(x=X[m].ravel(), y=jitter(y[m]), mode=\"markers\", name=name,\n",
    "                        marker=dict(color=color, opacity=0.5))\n",
    "    fig.update_layout(title=title, template=\"plotly_white\",\n",
    "                      xaxis_title=\"Mean Radius\", yaxis_title=\"Malignant\")\n",
    "    return fig\n",
    "\n",
    "x_grid = np.linspace(5, 30, 200).reshape(-1, 1)\n",
    "class_scatter(X_train, y_train, \"Mean Radius vs Malignant (jittered)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5f0cf57f",
   "metadata": {},
   "source": [
    "**Question:** Looking at this plot, what rule would you use to predict malignancy?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5c65ef25",
   "metadata": {},
   "source": [
    "---\n",
    "## Step 1: Try Least Squares Anyway\n",
    "\n",
    "_\"I suppose it is tempting, if the only tool you have is a hammer, to treat everything as if it were a nail.\"_ - Abraham Maslow\n",
    "\n",
    "The labels are 0 and 1. Those are numbers, so nothing stops us from fitting ordinary least squares.\n",
    "No regularization here: one feature, 341 points."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0877afb1",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.390944Z",
     "iopub.status.busy": "2026-09-22T21:10:37.390730Z",
     "iopub.status.idle": "2026-09-22T21:10:37.438564Z",
     "shell.execute_reply": "2026-09-22T21:10:37.437505Z"
    }
   },
   "outputs": [],
   "source": [
    "lin_reg = LinearRegression().fit(X_train, y_train)\n",
    "\n",
    "fig = class_scatter(X_train, y_train, \"Least Squares Fit to a 0/1 Target\")\n",
    "fig.add_scatter(x=x_grid.ravel(), y=lin_reg.predict(x_grid), mode=\"lines\",\n",
    "                name=\"Least squares\", line=dict(color=\"orange\", width=3))\n",
    "fig"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3207f10e",
   "metadata": {},
   "source": [
    "To classify, we need a **decision rule**: predict malignant if $f(x) > 0.5$.\n",
    "For a line $f(x) = w_0 + w_1 x$ the boundary is where $f(x) = 0.5$, i.e. $x^* = (0.5 - w_0)/w_1$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5ce519a4",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.440643Z",
     "iopub.status.busy": "2026-09-22T21:10:37.440459Z",
     "iopub.status.idle": "2026-09-22T21:10:37.445880Z",
     "shell.execute_reply": "2026-09-22T21:10:37.445038Z"
    }
   },
   "outputs": [],
   "source": [
    "def ls_boundary(model):\n",
    "    return (0.5 - model.intercept_) / model.coef_[0]\n",
    "\n",
    "ls_x_star = ls_boundary(lin_reg)\n",
    "ls_train_acc = np.mean((lin_reg.predict(X_train) > 0.5) == y_train)\n",
    "print(f\"Least-squares boundary: x* = {ls_x_star:.2f}\")\n",
    "print(f\"Training accuracy:      {ls_train_acc:.3f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "005597f5",
   "metadata": {},
   "source": [
    "Accuracy looks fine. But what does a prediction of 1.3, or -0.2, mean? Clipping to [0, 1] hides the problem;\n",
    "it does not make the output a probability."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1b9e3d62",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.447902Z",
     "iopub.status.busy": "2026-09-22T21:10:37.447680Z",
     "iopub.status.idle": "2026-09-22T21:10:37.553236Z",
     "shell.execute_reply": "2026-09-22T21:10:37.552124Z"
    }
   },
   "outputs": [],
   "source": [
    "fig = class_scatter(X_train, y_train, \"Clipping the Line to [0, 1]\")\n",
    "fig.add_scatter(x=x_grid.ravel(), y=lin_reg.predict(x_grid), mode=\"lines\",\n",
    "                name=\"Least squares\", line=dict(color=\"orange\", dash=\"dot\"))\n",
    "fig.add_scatter(x=x_grid.ravel(), y=np.clip(lin_reg.predict(x_grid), 0, 1), mode=\"lines\",\n",
    "                name=\"Clipped\", line=dict(color=\"green\", width=5))\n",
    "fig"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9a034bfb",
   "metadata": {},
   "source": [
    "### The failure: one very large tumor\n",
    "\n",
    "Suppose we observe a malignant tumor with mean radius 100. It is *obviously* malignant, and it is on the\n",
    "correct side of the boundary already. Watch what least squares does with it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e2ba9887",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.555131Z",
     "iopub.status.busy": "2026-09-22T21:10:37.554787Z",
     "iopub.status.idle": "2026-09-22T21:10:37.560826Z",
     "shell.execute_reply": "2026-09-22T21:10:37.559614Z"
    }
   },
   "outputs": [],
   "source": [
    "X_ex = np.vstack([X_train, [[100.0]]])\n",
    "y_ex = np.append(y_train, 1)\n",
    "\n",
    "lin_reg_ex = LinearRegression().fit(X_ex, y_ex)\n",
    "ls_x_star_ex = ls_boundary(lin_reg_ex)\n",
    "ls_acc_ex = np.mean((lin_reg_ex.predict(X_train) > 0.5) == y_train)\n",
    "\n",
    "print(f\"Boundary moved from {ls_x_star:.2f} to {ls_x_star_ex:.2f}\")\n",
    "print(f\"Training accuracy (original points) fell from {ls_train_acc:.3f} to {ls_acc_ex:.3f}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "57f06ad9",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.563049Z",
     "iopub.status.busy": "2026-09-22T21:10:37.562866Z",
     "iopub.status.idle": "2026-09-22T21:10:37.615490Z",
     "shell.execute_reply": "2026-09-22T21:10:37.614668Z"
    }
   },
   "outputs": [],
   "source": [
    "x_wide = np.linspace(5, 105, 300).reshape(-1, 1)\n",
    "fig = class_scatter(X_ex, y_ex, \"One Easy Point Moves the Boundary\")\n",
    "fig.add_scatter(x=x_wide.ravel(), y=lin_reg.predict(x_wide), mode=\"lines\",\n",
    "                name=\"Before\", line=dict(color=\"gray\"))\n",
    "fig.add_scatter(x=x_wide.ravel(), y=lin_reg_ex.predict(x_wide), mode=\"lines\",\n",
    "                name=\"After adding x = 100\", line=dict(color=\"orange\"))\n",
    "fig.add_vline(x=ls_x_star, line_dash=\"dot\", line_color=\"gray\")\n",
    "fig.add_vline(x=ls_x_star_ex, line_dash=\"dash\", line_color=\"black\")\n",
    "fig.update_yaxes(range=[-0.5, 1.5])\n",
    "fig"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "75b5115a",
   "metadata": {},
   "source": [
    "Least squares penalizes the new point for being \"too far above the line\", so it tilts the line toward it.\n",
    "That is the Gaussian-noise assumption behind least squares at work: it treats a correct, confident point\n",
    "as a large residual. We need a different model **and** a different objective."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "814c847f",
   "metadata": {},
   "source": [
    "---\n",
    "## Step 2: The Logistic Sigmoid\n",
    "\n",
    "$$\\sigma(a) = \\frac{1}{1 + e^{-a}}$$\n",
    "\n",
    "We derive it in lecture from the assumption that the **log-odds are linear** in $x$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f58fd744",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.617656Z",
     "iopub.status.busy": "2026-09-22T21:10:37.617469Z",
     "iopub.status.idle": "2026-09-22T21:10:37.707330Z",
     "shell.execute_reply": "2026-09-22T21:10:37.706160Z"
    }
   },
   "outputs": [],
   "source": [
    "def sigmoid(a):\n",
    "    return 1 / (1 + np.exp(-a))\n",
    "\n",
    "a = np.linspace(-8, 8, 200)\n",
    "px.line(x=a, y=sigmoid(a), labels={\"x\": \"a\", \"y\": \"σ(a)\"},\n",
    "        title=\"The Logistic Sigmoid\", template=\"plotly_white\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "94a3ebe3",
   "metadata": {},
   "source": [
    "### Checking the properties numerically\n",
    "\n",
    "Each property from the slide, verified on a grid of values instead of taken on faith."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9f1d4299",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.709762Z",
     "iopub.status.busy": "2026-09-22T21:10:37.709550Z",
     "iopub.status.idle": "2026-09-22T21:10:37.716167Z",
     "shell.execute_reply": "2026-09-22T21:10:37.715265Z"
    }
   },
   "outputs": [],
   "source": [
    "a = np.linspace(-6, 6, 13)\n",
    "h = 1e-5\n",
    "\n",
    "print(\"Symmetry  1 - σ(a) = σ(-a):         \", np.allclose(1 - sigmoid(a), sigmoid(-a)))\n",
    "numeric_deriv = (sigmoid(a + h) - sigmoid(a - h)) / (2 * h)\n",
    "print(\"Derivative σ'(a) = σ(a)(1 - σ(a)):  \", np.allclose(numeric_deriv, sigmoid(a) * (1 - sigmoid(a))))\n",
    "p = sigmoid(a)\n",
    "print(\"Inverse   logit(σ(a)) = a:          \", np.allclose(np.log(p / (1 - p)), a))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e30b239e",
   "metadata": {},
   "source": [
    "### Shifting and scaling: $\\sigma(Ax + B)$\n",
    "\n",
    "Bigger $A$ makes the curve steeper. Changing $B$ shifts it; the midpoint sits at $x = -B/A$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eedd8c13",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.718661Z",
     "iopub.status.busy": "2026-09-22T21:10:37.718472Z",
     "iopub.status.idle": "2026-09-22T21:10:37.760708Z",
     "shell.execute_reply": "2026-09-22T21:10:37.759299Z"
    }
   },
   "outputs": [],
   "source": [
    "x = np.linspace(-6, 6, 200)\n",
    "fig = go.Figure()\n",
    "for A, B in [(1, 0), (3, 0), (1, 2)]:\n",
    "    fig.add_scatter(x=x, y=sigmoid(A * x + B), mode=\"lines\", name=f\"A={A}, B={B}\")\n",
    "fig.update_layout(title=\"σ(Ax + B)\", template=\"plotly_white\")\n",
    "fig"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b0714841",
   "metadata": {},
   "source": [
    "---\n",
    "## Step 3: Fit Logistic Regression\n",
    "\n",
    "**Heads up:** scikit-learn's `LogisticRegression` adds L2 regularization by default (`C=1.0`).\n",
    "We have not introduced regularization yet, so we turn it off with `C=np.inf`. We come back to `C` next lecture."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eb31da02",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.763095Z",
     "iopub.status.busy": "2026-09-22T21:10:37.762872Z",
     "iopub.status.idle": "2026-09-22T21:10:37.776183Z",
     "shell.execute_reply": "2026-09-22T21:10:37.775215Z"
    }
   },
   "outputs": [],
   "source": [
    "log_reg = LogisticRegression(C=np.inf).fit(X_train, y_train)\n",
    "\n",
    "w0, w1 = log_reg.intercept_[0], log_reg.coef_[0, 0]\n",
    "print(f\"w0 (intercept) = {w0:.3f}\")\n",
    "print(f\"w1 (slope)     = {w1:.3f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1fdced72",
   "metadata": {},
   "source": [
    "`predict_proba` returns one column per class, in the order of `log_reg.classes_`. Column 1 is $P(\\text{malignant} \\mid x)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9c96fad2",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.779283Z",
     "iopub.status.busy": "2026-09-22T21:10:37.779057Z",
     "iopub.status.idle": "2026-09-22T21:10:37.790750Z",
     "shell.execute_reply": "2026-09-22T21:10:37.786261Z"
    }
   },
   "outputs": [],
   "source": [
    "print(log_reg.classes_)\n",
    "log_reg.predict_proba(X_train[:5]).round(3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bf0886c3",
   "metadata": {},
   "source": [
    "### The real decision boundary\n",
    "\n",
    "$\\sigma(w_0 + w_1 x) = 0.5$ exactly when $w_0 + w_1 x = 0$, so $x^* = -w_0 / w_1$.\n",
    "(The original notebook shaded the regions at the *mean* radius. That is not the boundary.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0ee84bec",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.793704Z",
     "iopub.status.busy": "2026-09-22T21:10:37.793466Z",
     "iopub.status.idle": "2026-09-22T21:10:37.874948Z",
     "shell.execute_reply": "2026-09-22T21:10:37.873554Z"
    }
   },
   "outputs": [],
   "source": [
    "lr_x_star = -w0 / w1\n",
    "print(f\"Logistic regression boundary: x* = {lr_x_star:.2f}\")\n",
    "\n",
    "fig = class_scatter(X_train, y_train, \"Logistic Regression vs Least Squares\")\n",
    "fig.add_scatter(x=x_grid.ravel(), y=log_reg.predict_proba(x_grid)[:, 1], mode=\"lines\",\n",
    "                name=\"Logistic regression  P(y=1|x)\", line=dict(color=\"black\", width=3))\n",
    "fig.add_scatter(x=x_grid.ravel(), y=lin_reg.predict(x_grid), mode=\"lines\",\n",
    "                name=\"Least squares\", line=dict(color=\"orange\", dash=\"dot\"))\n",
    "fig.add_vrect(x0=5, x1=lr_x_star, fillcolor=\"blue\", opacity=0.07, line_width=0)\n",
    "fig.add_vrect(x0=lr_x_star, x1=30, fillcolor=\"red\", opacity=0.07, line_width=0)\n",
    "fig.add_vline(x=lr_x_star, line_dash=\"dash\", annotation_text=\"decision boundary\")\n",
    "fig"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e3707385",
   "metadata": {},
   "source": [
    "### Does the large tumor still break things?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6d445877",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.877124Z",
     "iopub.status.busy": "2026-09-22T21:10:37.876903Z",
     "iopub.status.idle": "2026-09-22T21:10:37.886668Z",
     "shell.execute_reply": "2026-09-22T21:10:37.885159Z"
    }
   },
   "outputs": [],
   "source": [
    "log_reg_ex = LogisticRegression(C=np.inf).fit(X_ex, y_ex)\n",
    "lr_x_star_ex = -log_reg_ex.intercept_[0] / log_reg_ex.coef_[0, 0]\n",
    "\n",
    "print(f\"Least squares boundary:        {ls_x_star:.2f} -> {ls_x_star_ex:.2f}\")\n",
    "print(f\"Logistic regression boundary:  {lr_x_star:.2f} -> {lr_x_star_ex:.2f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2af42d35",
   "metadata": {},
   "source": [
    "The point at radius 100 already has $\\sigma(\\cdot) \\approx 1$, so it contributes almost nothing to the\n",
    "cross-entropy gradient. Logistic regression ignores points it already gets right with confidence."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "69008dbf",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.889305Z",
     "iopub.status.busy": "2026-09-22T21:10:37.889095Z",
     "iopub.status.idle": "2026-09-22T21:10:37.896466Z",
     "shell.execute_reply": "2026-09-22T21:10:37.895101Z"
    }
   },
   "outputs": [],
   "source": [
    "print(f\"Test accuracy, least squares:       {np.mean((lin_reg.predict(X_test) > 0.5) == y_test):.3f}\")\n",
    "print(f\"Test accuracy, logistic regression: {log_reg.score(X_test, y_test):.3f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5ed49e92",
   "metadata": {},
   "source": [
    "Same test accuracy on this split. The win is not accuracy on clean data: it is outputs we can read as\n",
    "probabilities and a boundary that does not move when we see an easy point."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "aa8d87d7",
   "metadata": {},
   "source": [
    "### Reading the weights: odds\n",
    "\n",
    "The model says $\\log \\frac{p}{1-p} = w_0 + w_1 x$. Increasing $x$ by 1 adds $w_1$ to the log-odds,\n",
    "which **multiplies the odds by $e^{w_1}$**."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5369f465",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.898875Z",
     "iopub.status.busy": "2026-09-22T21:10:37.898625Z",
     "iopub.status.idle": "2026-09-22T21:10:37.905133Z",
     "shell.execute_reply": "2026-09-22T21:10:37.903868Z"
    }
   },
   "outputs": [],
   "source": [
    "print(f\"Each extra unit of mean radius multiplies the odds of malignancy by e^w1 = {np.exp(w1):.2f}\")\n",
    "\n",
    "for r in [12, 13]:\n",
    "    p = log_reg.predict_proba([[r]])[0, 1]\n",
    "    print(f\"radius {r}: p = {p:.3f}, odds = {p / (1 - p):.3f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e3bfcdce",
   "metadata": {},
   "source": [
    "---\n",
    "## Step 4: Two Features, Linear Boundary\n",
    "\n",
    "With $D$ features the boundary $w^\\top x + w_0 = 0$ is a $(D-1)$-dimensional hyperplane, and $w$ is\n",
    "**orthogonal** to it. Here $D = 2$: the boundary is a line.\n",
    "We standardize the features so the two axes are on the same scale and the orthogonality is visible."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "10be6742",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.907671Z",
     "iopub.status.busy": "2026-09-22T21:10:37.907499Z",
     "iopub.status.idle": "2026-09-22T21:10:37.918446Z",
     "shell.execute_reply": "2026-09-22T21:10:37.917873Z"
    }
   },
   "outputs": [],
   "source": [
    "feats2 = [\"mean radius\", \"mean texture\"]\n",
    "scaler2 = StandardScaler().fit(data_tr[feats2])\n",
    "X2 = scaler2.transform(data_tr[feats2])\n",
    "\n",
    "log_reg2 = LogisticRegression(C=np.inf).fit(X2, y_train)\n",
    "w, b = log_reg2.coef_[0], log_reg2.intercept_[0]\n",
    "print(\"w =\", w.round(3), \" w0 =\", round(b, 3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3862245f",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.921680Z",
     "iopub.status.busy": "2026-09-22T21:10:37.921028Z",
     "iopub.status.idle": "2026-09-22T21:10:37.988078Z",
     "shell.execute_reply": "2026-09-22T21:10:37.987310Z"
    }
   },
   "outputs": [],
   "source": [
    "x1 = np.linspace(-2.5, 4, 50)\n",
    "x2_boundary = -(b + w[0] * x1) / w[1]\n",
    "\n",
    "fig = px.scatter(x=X2[:, 0], y=X2[:, 1], color=y_train.astype(str),\n",
    "                 color_discrete_map={\"0\": \"blue\", \"1\": \"red\"}, opacity=0.5,\n",
    "                 labels={\"x\": \"mean radius (std.)\", \"y\": \"mean texture (std.)\", \"color\": \"malignant\"},\n",
    "                 template=\"plotly_white\", title=\"2-D Decision Boundary and the Weight Vector\")\n",
    "fig.add_scatter(x=x1, y=x2_boundary, mode=\"lines\", name=\"w·x + w0 = 0\", line=dict(color=\"black\"))\n",
    "origin = np.array([x1[20], x2_boundary[20]])\n",
    "fig.add_annotation(x=origin[0] + w[0] / 3, y=origin[1] + w[1] / 3, ax=origin[0], ay=origin[1],\n",
    "                   xref=\"x\", yref=\"y\", axref=\"x\", ayref=\"y\", showarrow=True, arrowhead=3,\n",
    "                   arrowwidth=3, text=\"w\")\n",
    "fig.update_yaxes(range=[-2.5, 4], scaleanchor=\"x\", scaleratio=1)\n",
    "fig"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e89e9205",
   "metadata": {},
   "source": [
    "---\n",
    "## Step 5: Which Objective? Squared Error vs Cross-Entropy\n",
    "\n",
    "A small toy dataset, **no intercept**, so the model is $p = \\sigma(w x)$ with a single scalar $w$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c0129406",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:37.991276Z",
     "iopub.status.busy": "2026-09-22T21:10:37.990969Z",
     "iopub.status.idle": "2026-09-22T21:10:37.999516Z",
     "shell.execute_reply": "2026-09-22T21:10:37.998552Z"
    }
   },
   "outputs": [],
   "source": [
    "toy_df = pd.DataFrame({\"x\": [-4, -2, -0.5, 1, 3, 5],\n",
    "                       \"y\": [0, 0, 1, 0, 1, 1]})\n",
    "toy_df"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7729b6bf",
   "metadata": {},
   "source": [
    "### Squared error on top of the sigmoid is non-convex"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b02ca511",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:38.001122Z",
     "iopub.status.busy": "2026-09-22T21:10:38.000973Z",
     "iopub.status.idle": "2026-09-22T21:10:38.211034Z",
     "shell.execute_reply": "2026-09-22T21:10:38.209520Z"
    }
   },
   "outputs": [],
   "source": [
    "def ss_error_on_toy_data(w):\n",
    "    p_hat = sigmoid(toy_df[\"x\"] * w)\n",
    "    return 0.5 * np.sum((toy_df[\"y\"] - p_hat) ** 2)\n",
    "\n",
    "w_grid = np.linspace(-10, 10, 400)\n",
    "w_error = pd.DataFrame({\"w\": w_grid, \"Squared Error\": [ss_error_on_toy_data(w) for w in w_grid]})\n",
    "px.line(w_error, x=\"w\", y=\"Squared Error\", title=\"Squared Error on Toy Data\", template=\"plotly_white\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5f7da4c6",
   "metadata": {},
   "source": [
    "Run the optimizer from two starting points. Same data, same objective, different answers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1b4a0d60",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:38.213435Z",
     "iopub.status.busy": "2026-09-22T21:10:38.213215Z",
     "iopub.status.idle": "2026-09-22T21:10:38.243873Z",
     "shell.execute_reply": "2026-09-22T21:10:38.242928Z"
    }
   },
   "outputs": [],
   "source": [
    "best_w = minimize(ss_error_on_toy_data, x0=0)[\"x\"][0]\n",
    "stuck_w = minimize(ss_error_on_toy_data, x0=-5)[\"x\"][0]\n",
    "print(f\"start at w =  0  ->  w = {best_w:.3f}\")\n",
    "print(f\"start at w = -5  ->  w = {stuck_w:.3f}   (stuck on the flat plateau)\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f006b2ab",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:38.246748Z",
     "iopub.status.busy": "2026-09-22T21:10:38.246541Z",
     "iopub.status.idle": "2026-09-22T21:10:38.296387Z",
     "shell.execute_reply": "2026-09-22T21:10:38.295546Z"
    }
   },
   "outputs": [],
   "source": [
    "fig = px.scatter(toy_df, x=\"x\", y=\"y\", template=\"plotly_white\")\n",
    "fig.update_traces(marker_size=14)\n",
    "xs = np.linspace(-6, 6, 200)\n",
    "fig.add_scatter(x=xs, y=sigmoid(xs * best_w), mode=\"lines\", name=f\"w = {best_w:.2f}\")\n",
    "fig.add_scatter(x=xs, y=sigmoid(xs * stuck_w), mode=\"lines\", name=f\"w = {stuck_w:.2f}\")\n",
    "fig"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4b4aeaf6",
   "metadata": {},
   "source": [
    "### Squared error is also bounded\n",
    "\n",
    "For one point with $y = 1$: the worst squared error (with the ½) is 0.5, no matter how confidently wrong\n",
    "the model is. Negative log-probability goes to infinity."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "03a18a03",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:38.299267Z",
     "iopub.status.busy": "2026-09-22T21:10:38.299068Z",
     "iopub.status.idle": "2026-09-22T21:10:38.356191Z",
     "shell.execute_reply": "2026-09-22T21:10:38.355032Z"
    }
   },
   "outputs": [],
   "source": [
    "p = np.linspace(0.001, 0.999, 200)\n",
    "loss_y1 = pd.DataFrame({\"p\": p, \"Squared Error\": 0.5 * (1 - p) ** 2, \"Neg Log\": -np.log(p)})\n",
    "px.line(loss_y1, x=\"p\", y=[\"Squared Error\", \"Neg Log\"], template=\"plotly_white\",\n",
    "        title=\"Loss for One Point with y = 1\", labels={\"value\": \"loss\", \"p\": \"predicted P(y=1)\"})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "73b994b8",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:38.358361Z",
     "iopub.status.busy": "2026-09-22T21:10:38.358125Z",
     "iopub.status.idle": "2026-09-22T21:10:38.414615Z",
     "shell.execute_reply": "2026-09-22T21:10:38.413608Z"
    }
   },
   "outputs": [],
   "source": [
    "loss_y0 = pd.DataFrame({\"p\": p, \"Squared Error\": 0.5 * p ** 2, \"Neg Log\": -np.log(1 - p)})\n",
    "px.line(loss_y0, x=\"p\", y=[\"Squared Error\", \"Neg Log\"], template=\"plotly_white\",\n",
    "        title=\"Loss for One Point with y = 0\", labels={\"value\": \"loss\", \"p\": \"predicted P(y=1)\"})"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "079dd149",
   "metadata": {},
   "source": [
    "### Cross-entropy on the same toy data\n",
    "\n",
    "Computing $\\log \\sigma(a)$ directly returns `-inf` once $\\sigma(a)$ rounds to 0 or 1.\n",
    "We use the identity $-\\log \\sigma(a) = \\log(1 + e^{-a})$, which `np.logaddexp(0, -a)` computes stably."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "074a94ce",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:38.417055Z",
     "iopub.status.busy": "2026-09-22T21:10:38.416851Z",
     "iopub.status.idle": "2026-09-22T21:10:38.649525Z",
     "shell.execute_reply": "2026-09-22T21:10:38.647974Z"
    }
   },
   "outputs": [],
   "source": [
    "def cross_entropy_from_logits(t, a):\n",
    "    # t*(-log σ(a)) + (1-t)*(-log(1-σ(a))), written stably\n",
    "    return t * np.logaddexp(0, -a) + (1 - t) * np.logaddexp(0, a)\n",
    "\n",
    "def mean_ce_on_toy_data(w):\n",
    "    return np.mean(cross_entropy_from_logits(toy_df[\"y\"], toy_df[\"x\"] * w))\n",
    "\n",
    "w_error[\"Cross-Entropy\"] = [mean_ce_on_toy_data(w) for w in w_grid]\n",
    "px.line(w_error, x=\"w\", y=\"Cross-Entropy\", title=\"Cross-Entropy on Toy Data (convex)\",\n",
    "        template=\"plotly_white\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e8f67274",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-09-22T21:10:38.652377Z",
     "iopub.status.busy": "2026-09-22T21:10:38.652160Z",
     "iopub.status.idle": "2026-09-22T21:10:38.682383Z",
     "shell.execute_reply": "2026-09-22T21:10:38.680885Z"
    }
   },
   "outputs": [],
   "source": [
    "ce_w_from_0 = minimize(mean_ce_on_toy_data, x0=0)[\"x\"][0]\n",
    "ce_w_from_neg5 = minimize(mean_ce_on_toy_data, x0=-5)[\"x\"][0]\n",
    "print(f\"start at  0 -> w = {ce_w_from_0:.3f}\")\n",
    "print(f\"start at -5 -> w = {ce_w_from_neg5:.3f}   (same answer)\")"
   ]
  }
 ],
 "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.9.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
