{
 "cells": [
  {
   "cell_type": "raw",
   "id": "3da81c08",
   "metadata": {},
   "source": [
    "---\n",
    "title: \"Lab 3 solutions: selecting and explaining a tabular model\"\n",
    "chapter: 3\n",
    "lab-solutions: true\n",
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ac9404bf",
   "metadata": {},
   "source": [
    "# Lab 3: select, explain, then test\n",
    "\n",
    "**Goal.** Compare ridge, random forests, gradient boosting and support vector regression on the same solubility problem. Evaluate the entire family-and-hyperparameter selection procedure with nested cross-validation, interpret development behaviour, then freeze a candidate before opening the final test.\n",
    "\n",
    "Run top to bottom. In the student notebook, complete each `# TODO` using the requested object names. The checks verify useful invariants, not a particular winning algorithm. The search uses a modest CPU budget; no GPU is needed. Record the search space before viewing results.\n",
    "\n",
    "**Deliverables:** outer-fold results against a baseline, a development importance plot, one final test comparison and parity plot, and a qualified scientific recommendation. Budget 90–120 minutes for the tasks and reflection."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "66d9eb36",
   "metadata": {
    "tags": [
     "setup"
    ]
   },
   "outputs": [],
   "source": [
    "# Setup: detect Colab and install missing core packages there.\n",
    "import importlib.util\n",
    "import subprocess\n",
    "import sys\n",
    "try:\n",
    "    import google.colab\n",
    "    IN_COLAB = True\n",
    "except ImportError:\n",
    "    IN_COLAB = False\n",
    "packages = {'numpy': 'numpy', 'pandas': 'pandas', 'matplotlib': 'matplotlib',\n",
    "            'sklearn': 'scikit-learn', 'threadpoolctl': 'threadpoolctl'}\n",
    "missing = [package for module, package in packages.items()\n",
    "           if importlib.util.find_spec(module) is None]\n",
    "if missing and IN_COLAB:\n",
    "    subprocess.check_call([sys.executable, '-m', 'pip', 'install', '-q', *missing])\n",
    "if missing and not IN_COLAB:\n",
    "    raise RuntimeError('Install the course requirements before running this notebook.')\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import sklearn\n",
    "from threadpoolctl import threadpool_limits\n",
    "print(f'Python {sys.version.split()[0]} | sklearn {sklearn.__version__} | Colab: {IN_COLAB}')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "56db0a97",
   "metadata": {},
   "source": [
    "## Load the data\n",
    "\n",
    "The ESOL file is downloaded rather than redistributed. This loader first checks local paths, then tries the canonical download and its DeepChem GitHub mirror and verifies the same SHA-256 as `data/README.md`. The five descriptors and target units match Lab 2. Keeping the representation fixed makes the algorithm comparison interpretable."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a12a62b0",
   "metadata": {
    "tags": [
     "setup"
    ]
   },
   "outputs": [],
   "source": [
    "# Load a verified ESOL copy locally, then try the canonical and mirrored downloads.\n",
    "from pathlib import Path\n",
    "from urllib.request import urlopen\n",
    "from urllib.error import URLError\n",
    "from io import BytesIO\n",
    "import hashlib\n",
    "candidates = [Path('data/esol.csv'), Path('../../../data/esol.csv')]\n",
    "local = next((p for p in candidates if p.exists()), None)\n",
    "urls = ['https://deepchemdata.s3-us-west-1.amazonaws.com/datasets/delaney-processed.csv',\n",
    "        'https://raw.githubusercontent.com/deepchem/deepchem/master/datasets/delaney-processed.csv']\n",
    "raw = local.read_bytes() if local is not None else None\n",
    "if raw is None:\n",
    "    for url in urls:\n",
    "        try:\n",
    "            with urlopen(url, timeout=30) as response:\n",
    "                raw = response.read()\n",
    "            break\n",
    "        except (URLError, TimeoutError):\n",
    "            continue\n",
    "if raw is None:\n",
    "    raise RuntimeError('Could not download ESOL; check your connection or run scripts/fetch_data.py.')\n",
    "expected = '8c06a76f0c6487d29ab0f903e6a7a7139f189ab3c1178f159c8be8964602f189'\n",
    "assert hashlib.sha256(raw).hexdigest() == expected, 'ESOL checksum differs; inspect the source.'\n",
    "esol = pd.read_csv(BytesIO(raw))\n",
    "print(f'Loaded {len(esol)} molecules.')\n",
    "esol.head(3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "443b9e53",
   "metadata": {},
   "source": [
    "## Task 1 — Declare the experiment and reserve the test\n",
    "\n",
    "Create `features`, `X` and `y` with the five descriptors used in Chapter 2. Make an 80/20 split with `random_state=42`, naming the parts `X_dev`, `X_test`, `y_dev`, `y_test`. The metric is RMSE in log10(mol/L). The final test remains closed until Task 7.\n",
    "\n",
    "Write down the limitation now: a random molecular split allows related chemistry across partitions and does not estimate performance on genuinely new scaffolds. The repeated teaching split is familiar to course readers; a research claim needs fresh independent evidence."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "74d35dde",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 1: define the representation and reserve the final test.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f3c583b7",
   "metadata": {
    "tags": [
     "check"
    ]
   },
   "outputs": [],
   "source": [
    "assert list(X.columns) == features and len(features) == 5, 'Use the declared five descriptors.'\n",
    "assert set(X_dev.index).isdisjoint(X_test.index), 'Partitions must not overlap.'\n",
    "assert len(X_dev) + len(X_test) == len(X), 'Account for every row.'\n",
    "assert y_dev.index.equals(X_dev.index), 'Keep feature and target rows aligned.'\n",
    "print('Check passed: representation and split.')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "99bc6fde",
   "metadata": {},
   "source": [
    "## Task 2 — Define one search across families and settings\n",
    "\n",
    "Build `pipeline` with named steps `scale` and `model`. Build a list of four dictionaries called `grids`. Each dictionary replaces the model with one family and supplies only its valid parameters. Use scaling for ridge and SVR, and `\"passthrough\"` for tree models.\n",
    "\n",
    "Use two ridge alphas (1, 100), two forest leaf sizes (2, 8; 60 trees), two boosting learning rates (0.05, 0.1; 80 iterations), and four SVR combinations (C 1 or 10; gamma 0.03 or 0.3). This is ten configurations. `inner` and `outer` are shuffled three-fold KFold objects with fixed seeds; `search` is a GridSearchCV with negative RMSE and `n_jobs=1`. Set boosting early stopping off so iteration count is explicit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9d76f9c0",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 2a: define the pipeline and modest candidate grid.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "98853996",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 2b: define the nested evaluation before looking at results.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "459b4b08",
   "metadata": {
    "tags": [
     "check"
    ]
   },
   "outputs": [],
   "source": [
    "assert len(ParameterGrid(grids)) == 10, 'Keep the declared ten-configuration budget.'\n",
    "assert search.scoring == 'neg_root_mean_squared_error', 'Use the shared RMSE criterion.'\n",
    "assert set(pipeline.named_steps) == {'scale', 'model'}, 'Both steps belong inside the search.'\n",
    "print('Check passed: search budget and pipeline.')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7cc95cdd",
   "metadata": {},
   "source": [
    "## Task 3 — Evaluate the complete choice\n",
    "\n",
    "Call `cross_validate(search, ...)` on development data using `outer_splits`, RMSE scoring and `return_estimator=True`. Store the result as `nested`. Each outer training fold runs a fresh inner search over every family. Each outer validation fold is then predicted by the winner refitted on that outer training fold.\n",
    "\n",
    "Make `outer_results` containing fold number, outer RMSE, selected inner CV RMSE and selected family. The inner selected score helped choose a model, while the outer score assesses that choice on excluded rows. Different outer folds may choose different families; the evaluated object is the selection procedure, not one fixed family."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "68b50a51",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 3: run nested CV and expose what was selected in each fold.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "89c7dcc2",
   "metadata": {
    "tags": [
     "check"
    ]
   },
   "outputs": [],
   "source": [
    "assert len(nested['estimator']) == 3, 'There should be one full search per outer fold.'\n",
    "assert np.isfinite(outer_results['outer RMSE']).all(), 'All folds need finite results.'\n",
    "assert all(hasattr(s, 'best_estimator_') for s in nested['estimator']), 'Evaluate a search, not an untuned estimator.'\n",
    "print('Check passed: each outer fold fitted its own search.')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c598d3b6",
   "metadata": {},
   "source": [
    "## Task 4 — Compare with a baseline on identical outer rows\n",
    "\n",
    "Evaluate `DummyRegressor` using the same `outer_splits`. Append its fold RMSEs to `outer_results` and report mean and sample standard deviation of both procedures. Each baseline must learn its mean from that fold's training targets.\n",
    "\n",
    "The standard deviation describes variation across these folds. The folds share training observations, so mean ± SD is not a confidence interval. A small sample of folds cannot settle whether one method is universally better. Also, do not select a family from this table and claim the same outer score independently evaluates that extra selection step."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "022759c6",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 4: evaluate the mean baseline under the same outer split design.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4fce7725",
   "metadata": {
    "tags": [
     "check"
    ]
   },
   "outputs": [],
   "source": [
    "assert summary.shape == (2, 2), 'Report mean and SD for both procedures.'\n",
    "assert (outer_results['baseline RMSE'] > 0).all(), 'Baseline errors must be computed per fold.'\n",
    "print('Check passed: comparable baseline results.')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0473cd1d",
   "metadata": {},
   "source": [
    "## Task 5 — Interpret one development-fold model\n",
    "\n",
    "Use the fitted search from the first outer fold and that fold's validation indices. Compute permutation importance with ten repeats and negative RMSE. Plot the increase in RMSE by descriptor. Store the result as `importance`.\n",
    "\n",
    "Do not use `X_test`. This is a descriptive inspection of the evaluated procedure. If it motivates a new feature set or grid, record a new experiment; do not silently retain the old outer assessment as an independent estimate of the revised procedure. Correlated descriptors may substitute for each other, and shuffle variability is not an uncertainty interval over new datasets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7c35e8d1",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 5a: compute importance on rows held out from this particular search.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b2ab8eb3",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 5b: plot the model's development-fold reliance on each descriptor.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b018dd31",
   "metadata": {
    "tags": [
     "check"
    ]
   },
   "outputs": [],
   "source": [
    "assert importance.importances.shape == (len(features), 10), 'Repeat every descriptor shuffle ten times.'\n",
    "assert set(X_dev.iloc[validation_rows].index).isdisjoint(X_test.index), 'Interpret development data only.'\n",
    "print('Check passed: importance stays within development data.')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "04a0fe87",
   "metadata": {},
   "source": [
    "## Task 6 — Refit the unchanged search and freeze its winner\n",
    "\n",
    "Clone `search`, fit it on all development rows and call it `final_search`. Assign `best_model = final_search.best_estimator_`. Record the selected class name and parameters. The final choice comes from the inner CV on all development data, not from comparing test outcomes or selecting the best outer-fold model.\n",
    "\n",
    "Before continuing, write your selection rule: “Choose the maximum negative-RMSE score among the ten declared configurations, refit on all development rows, then evaluate once.” No tuning step follows Task 7."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8211396d",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 6: choose using development CV and freeze the complete pipeline.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3870e628",
   "metadata": {
    "tags": [
     "check"
    ]
   },
   "outputs": [],
   "source": [
    "assert best_model is final_search.best_estimator_, 'Take the winner from development CV.'\n",
    "assert hasattr(best_model.named_steps['model'], 'n_features_in_'), 'Refit the selected model.'\n",
    "print('Check passed: the chosen pipeline is fitted before the final test.')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c3a70e00",
   "metadata": {},
   "source": [
    "## Task 7 — Open the final test once\n",
    "\n",
    "Fit the declared mean baseline on all development rows. Make `pred_test` from the frozen pipeline and create `test_results` with test MAE, RMSE and R² for the candidate and baseline. Plot the candidate's measured versus predicted logS.\n",
    "\n",
    "A disappointing result is a valid result. Report it without switching to a model that looks better on these test labels. Error inspection may generate hypotheses for a future experiment; it does not grant another independent use of the same test set."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9d808b8f",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 7a: evaluate the frozen candidate and declared baseline once.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c9f8cbbf",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Task 7b: plot the candidate selected before the test was opened.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "424223f2",
   "metadata": {
    "tags": [
     "check"
    ]
   },
   "outputs": [],
   "source": [
    "assert len(pred_test) == len(y_test), 'Predict each held-out row.'\n",
    "assert list(test_results.index) == ['frozen candidate', 'mean baseline'], 'Evaluate the declared comparison.'\n",
    "assert np.isfinite(test_results.to_numpy()).all(), 'Report finite metrics.'\n",
    "print('Check passed: final test report is complete.')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cfedded4",
   "metadata": {},
   "source": [
    "## Task 8 — Write a qualified recommendation\n",
    "\n",
    "Complete this paragraph with your outputs:\n",
    "\n",
    "> On five ESOL descriptors, the declared family-and-setting selection procedure achieved a mean outer CV RMSE of [fill in] log10(mol/L), with fold SD [fill in] (a descriptive spread, not a confidence interval), compared with [fill in] for the mean baseline. Refitting the same search on all development data selected [fill in]. Its final test RMSE was [fill in], MAE [fill in] and R² [fill in]. The importance analysis suggests this fold model relied on [fill in], but does not establish causation. I would / would not recommend this procedure for [fill in] because [fill in]. A random split can share related chemistry, so the next independent test should [fill in] .\n",
    "\n",
    "Also record the library versions, features, seeds, ten candidate configurations and fold design. The model is the complete fitted pipeline, including preprocessing.\n",
    "\n",
    "**Reflection:** Did outer folds choose the same family? How large was the advantage over the baseline relative to fold-to-fold variation? Would fewer descriptors or a cheaper algorithm be acceptable for the intended use? Which claim needs a scaffold or prospective split?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b93e2103",
   "metadata": {},
   "source": [
    "## Optional extension — PCA plus logistic regression versus a forest\n",
    "\n",
    "Repeat nested selection on the breast-cancer teaching dataset. Use a stratified 80/20 split, three stratified inner and outer folds, and balanced accuracy. Include PCA dimensions 5 and 15 inside a scaled logistic pipeline and forest leaf sizes 2 and 8 in an unscaled alternative. Keep the final test closed until the search is refitted on all development rows.\n",
    "\n",
    "This demonstrates a second workflow, not a clinical validation. The label mapping is 0 malignant, 1 benign. Balanced accuracy gives equal weight to class recalls but does not encode the relative costs of mistakes, nor establish probability calibration."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2f1f315b",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Optional extension: define a complete classification search.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1c5b7a72",
   "metadata": {
    "tags": [
     "solution"
    ]
   },
   "outputs": [],
   "source": [
    "# Optional extension: assess selection, refit it, then evaluate its frozen winner.\n",
    "# TODO: your code here"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2a05f1fe",
   "metadata": {},
   "source": [
    "## Closing reflection\n",
    "\n",
    "1. Why does putting the family choice in the inner loop matter?\n",
    "2. Which information would leak if PCA were fitted before cross-validation?\n",
    "3. Why are a good classification score, calibrated probabilities and a useful decision threshold separate requirements?\n",
    "4. If you now change a grid after seeing the final results, what evidence would independently evaluate that new procedure?\n",
    "\n",
    "**Further reading:** [scikit-learn's nested CV example](https://scikit-learn.org/stable/auto_examples/model_selection/plot_nested_cross_validation_iris.html) and [permutation-importance guide](https://scikit-learn.org/stable/modules/permutation_importance.html)."
   ]
  }
 ],
 "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.12.14"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}