{"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.11.0","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":14774,"sourceType":"competition"}],"isInternetEnabled":true,"language":"python","sourceType":"notebook"},"accelerator":"GPU"},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Diabetic Retinopathy Stage Detection — Computer Vision Coursework\n\n**Module:** Computer Vision &nbsp;·&nbsp; **Batch:** BSCCOMP24.2P &nbsp;·&nbsp; **Assessment:** CW1 (100 marks)\n\n---\n\n## Abstract\n\nDiabetic retinopathy (DR) is the leading cause of preventable blindness in working-age adults, yet\nscreening depends on scarce ophthalmologist time. This project builds an automated DR **stage**\nclassifier (5 ordinal classes) from retinal fundus photographs, developed as **four progressively\nrefined versions** so that the contribution of each design decision can be measured in isolation:\n\n| Version | What changes | What it isolates |\n|---|---|---|\n| **V1 — Baseline** | CNN written from scratch, minimal preprocessing, no augmentation | Lower bound; does a hand-built CNN learn anything? |\n| **V2 — Data-centric** | *Same architecture* + fundus-specific preprocessing, augmentation, class balancing | The contribution of **data quality**, holding the model fixed |\n| **V3 — Transfer learning** | *Same data* + ImageNet-pretrained EfficientNet-B3, two-stage fine-tuning | The contribution of **learned representations**, holding data fixed |\n| **V4 — Optimised** | *Same backbone* + MixUp/CutMix, label smoothing, cosine warm-up, TTA | The contribution of **training regularisation and inference strategy** |\n\nBecause exactly one factor changes between consecutive versions, the version-to-version deltas form\na controlled **ablation study** rather than a sequence of unrelated experiments.\n\n---\n\n## How this notebook maps to the marking rubric\n\n| # | Rubric criterion | Marks | Where it is addressed |\n|---|---|---|---|\n| 1 | Problem Understanding & Dataset Justification | 10 | §1 — clinical background, EDA, stratified **train/val/test** split, ethics & limitations |\n| 2 | Data Preprocessing Techniques | 10 | §2 — circle crop, Ben Graham illumination correction, CLAHE, denoising, edge enhancement, reproducible cache |\n| 3 | Data Augmentation & Dataset Balancing | 10 | §3 — geometric/photometric augmentation, class weights, `WeightedRandomSampler` |\n| 4 | CNN Architecture & Transfer Learning | **20** | §5 (scratch CNN), §7 (EfficientNet-B3 two-stage fine-tuning, discriminative LRs, **empirical backbone study** + fine-tuning strategy study) |\n| 5 | Training Strategy & Experimental Design | 10 | §4 — early stopping, checkpointing, LR scheduling, AMP, gradient clipping, controlled ablation design |\n| 6 | Model Evaluation & Performance Analysis | **15** | §9 — accuracy, macro/weighted precision·recall·F1, confusion matrices, curves, **quadratic weighted kappa**, ordinal error analysis |\n| 7 | Code Quality & Documentation | 10 | Typed config, docstringed functions, single source of truth, deterministic seeding, no copy-paste between versions |\n| 8 | Report Quality & Presentation | 10 | Every figure is saved to `figures/` at 150 dpi and bundled for the report; §12 exports the archive |\n| 9 | Innovation, Practical Impact & Critical Discussion | 5 | §8 Grad-CAM explainability, §11 Gradio prototype, §12 deployment feasibility, ethics, future work |\n\n---\n\n> **Before running:** see §0.1 for the two Kaggle settings that must be enabled (GPU and Internet).","metadata":{}},{"cell_type":"markdown","source":"## 0. Environment setup\n\n### 0.1 Required Kaggle settings — read this first\n\nThis notebook is designed for a **Kaggle Notebook** session, and also runs unmodified in Google Colab.\n\nIn the Kaggle editor, open the **right-hand sidebar** and set:\n\n1. **Add Data** → search `aptos2019-blindness-detection` → add the competition dataset.\n   (You must first visit the competition page and click *Join Competition* to accept the rules.)\n2. **Accelerator** → `GPU T4 x2` or `GPU P100`.\n   Without a GPU the notebook still runs, but training will take many hours instead of ~1 hour.\n3. **Internet** → `On`.\n   This is required **once**, to download the ImageNet-pretrained weights used in V3/V4.\n   Enabling internet requires a phone-verified Kaggle account.\n   If you cannot enable it, see the fallback note in §7.1 — the notebook degrades gracefully\n   rather than crashing.\n\nThen: **Run All**. Total runtime is roughly **60–90 minutes** on a P100.\n\n> **Tip for a fast smoke test:** set `CFG.fast_dev_run = True` in §0.3. The whole notebook then\n> completes in a few minutes on a tiny subset so you can verify everything works end-to-end\n> before committing to a full run.","metadata":{}},{"cell_type":"code","source":"# ─────────────────────────────────────────────────────────────────────────────\n# Imports. Everything used here ships with the standard Kaggle/Colab image,\n# so no `pip install` is required for the core pipeline.\n# ─────────────────────────────────────────────────────────────────────────────\nimport os, sys, json, time, math, random, shutil, warnings\nfrom pathlib import Path\nfrom dataclasses import dataclass, asdict\nfrom concurrent.futures import ThreadPoolExecutor\n\nimport numpy as np\nimport pandas as pd\nimport cv2\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nfrom PIL import Image\n\nimport torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.data import Dataset, DataLoader, WeightedRandomSampler\n\nimport torchvision\nfrom torchvision import transforms as T\n\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.metrics import (\n    accuracy_score, precision_recall_fscore_support, classification_report,\n    confusion_matrix, cohen_kappa_score,\n)\n\nwarnings.filterwarnings(\"ignore\")\n\nprint(f\"Python       : {sys.version.split()[0]}\")\nprint(f\"NumPy        : {np.__version__}\")\nprint(f\"OpenCV       : {cv2.__version__}\")\nprint(f\"PyTorch      : {torch.__version__}\")\nprint(f\"torchvision  : {torchvision.__version__}\")\nprint(f\"CUDA available: {torch.cuda.is_available()}\")\nif torch.cuda.is_available():\n    print(f\"GPU          : {torch.cuda.get_device_name(0)}\")\n    print(f\"GPU memory   : {torch.cuda.get_device_properties(0).total_memory / 1e9:.1f} GB\")\nelse:\n    print(\"WARNING: no GPU detected. Enable the accelerator in the Kaggle sidebar.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### 0.2 Reproducibility\n\nEvery stochastic component — weight initialisation, data shuffling, augmentation sampling,\ndropout — is seeded from a single constant. `cudnn.deterministic` is left **off** by default\nbecause it roughly halves throughput; set `CFG.deterministic = True` if exact bitwise\nreproducibility matters more than speed.\n\nThis matters for the ablation: if V2 beats V1, we need confidence that the gap comes from the\npreprocessing change and not from a lucky random initialisation. All versions use the same seed.","metadata":{}},{"cell_type":"code","source":"def seed_everything(seed: int = 42, deterministic: bool = False) -> None:\n    \"\"\"Seed every RNG the pipeline touches so runs are repeatable.\n\n    Args:\n        seed: The master seed applied to Python, NumPy and PyTorch (CPU + CUDA).\n        deterministic: If True, force cuDNN into deterministic mode. This removes\n            non-determinism from convolution algorithm selection at a ~2x throughput\n            cost, so it is off by default.\n    \"\"\"\n    random.seed(seed)\n    os.environ[\"PYTHONHASHSEED\"] = str(seed)\n    np.random.seed(seed)\n    torch.manual_seed(seed)\n    torch.cuda.manual_seed_all(seed)\n    if deterministic:\n        torch.backends.cudnn.deterministic = True\n        torch.backends.cudnn.benchmark = False\n    else:\n        torch.backends.cudnn.benchmark = True   # autotune conv algorithms for fixed input sizes\n\n\ndef worker_init_fn(worker_id: int) -> None:\n    \"\"\"Give every DataLoader worker a distinct but reproducible seed.\n\n    Without this, forked workers inherit the parent's NumPy state and can emit\n    identical \"random\" augmentations within a batch.\n    \"\"\"\n    np.random.seed(torch.initial_seed() % 2**32 + worker_id)","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### 0.3 Central configuration\n\nEvery tunable lives in one frozen dataclass. Nothing downstream hard-codes a magic number, so a\nsingle edit here changes the whole experiment — which is what makes the four versions genuinely\ncomparable, and what the rubric means by *reproducible*.\n\n**Note on image size:** 300×300 is used for *all four* versions. It is the native training\nresolution of EfficientNet-B3, and keeping it fixed across versions means the V2→V3 comparison\nmeasures the architecture change alone rather than a resolution change confounded with it.","metadata":{}},{"cell_type":"code","source":"@dataclass\nclass Config:\n    \"\"\"Single source of truth for every hyperparameter in the project.\"\"\"\n\n    # ── Reproducibility ──────────────────────────────────────────────────────\n    seed: int = 42\n    deterministic: bool = False\n\n    # ── Data ─────────────────────────────────────────────────────────────────\n    img_size: int = 300           # native resolution of EfficientNet-B3\n    cache_size: int = 380         # images cached slightly larger, so random crops have headroom\n    num_classes: int = 5\n    val_frac: float = 0.15        # of the full dataset\n    test_frac: float = 0.15       # of the full dataset -> 70/15/15 train/val/test\n    num_workers: int = 2\n\n    # ── Training ─────────────────────────────────────────────────────────────\n    batch_size: int = 32\n    epochs_v1: int = 25\n    epochs_v2: int = 30\n    epochs_v3_head: int = 5       # stage A: frozen backbone, train new head only\n    epochs_v3_full: int = 25      # stage B: unfrozen, discriminative learning rates\n    epochs_v4: int = 30\n    early_stop_patience: int = 8\n    grad_clip: float = 1.0\n    use_amp: bool = True          # mixed precision: ~2x faster, ~40% less VRAM\n\n    # ── Learning rates ───────────────────────────────────────────────────────\n    lr_scratch: float = 1e-3      # training from random init\n    lr_head: float = 1e-3         # new classifier head on a frozen backbone\n    lr_backbone: float = 1e-5     # pretrained features: tiny LR to avoid catastrophic forgetting\n    weight_decay: float = 1e-4\n\n    # ── V4 regularisation ────────────────────────────────────────────────────\n    label_smoothing: float = 0.1\n    mixup_alpha: float = 0.4\n    cutmix_alpha: float = 1.0\n    mix_prob: float = 0.5         # probability a given batch is mixed at all\n    warmup_epochs: int = 3\n    tta_rounds: int = 5           # test-time augmentation passes\n\n    # ── Model selection ──────────────────────────────────────────────────────\n    monitor: str = \"qwk\"          # checkpoint on quadratic weighted kappa (see §4.2)\n\n    # ── Debug ────────────────────────────────────────────────────────────────\n    fast_dev_run: bool = False    # True -> tiny subset + 2 epochs, for smoke-testing\n\n    def __post_init__(self):\n        if self.fast_dev_run:\n            self.epochs_v1 = self.epochs_v2 = self.epochs_v4 = 2\n            self.epochs_v3_head, self.epochs_v3_full = 1, 2\n            self.early_stop_patience = 99\n            self.tta_rounds = 2\n\n\nCFG = Config()\nseed_everything(CFG.seed, CFG.deterministic)\n\nDEVICE = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nprint(f\"Device: {DEVICE}\")\nprint(\"\\nConfiguration:\")\nfor k, v in asdict(CFG).items():\n    print(f\"  {k:24s} = {v}\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### 0.4 Paths, plotting style and output directories\n\nThe dataset root is **discovered** rather than hard-coded, so the same notebook runs on Kaggle\n(`/kaggle/input/...`), Colab (`/content/...`) or a local machine with no edits.\n\nA single matplotlib style and a fixed colour-blind-safe palette are set here. Every figure in the\nnotebook inherits them, so the plots that go into the report read as one consistent visual system\ninstead of nine different default styles.","metadata":{}},{"cell_type":"code","source":"def find_dataset_root() -> Path:\n    \"\"\"Locate the APTOS 2019 dataset directory across Kaggle / Colab / local.\n\n    Searches common roots for a directory containing both `train.csv` and a\n    `train_images` folder.\n\n    Returns:\n        Path to the dataset root.\n\n    Raises:\n        FileNotFoundError: if no matching directory is found, with guidance on\n            how to attach the dataset.\n    \"\"\"\n    search_roots = [\"/kaggle/input\", \"/content\", \"./dataset\", \".\", \"..\"]\n    for base in search_roots:\n        base_p = Path(base)\n        if not base_p.is_dir():\n            continue\n        for root, dirs, files in os.walk(base_p):\n            if \"train.csv\" in files and \"train_images\" in dirs:\n                return Path(root)\n            # keep the walk shallow; the datasets we want are never deeply nested\n            if len(Path(root).parts) - len(base_p.parts) >= 3:\n                dirs[:] = []\n    raise FileNotFoundError(\n        \"Could not find the APTOS 2019 dataset.\\n\"\n        \"On Kaggle: sidebar -> Add Data -> search 'aptos2019-blindness-detection' -> Add.\\n\"\n        \"Elsewhere: place train.csv and train_images/ under ./dataset/.\"\n    )\n\n\nDATA_ROOT  = find_dataset_root()\nTRAIN_IMG  = DATA_ROOT / \"train_images\"\nTRAIN_CSV  = DATA_ROOT / \"train.csv\"\n\n# Writable working directory (Kaggle gives /kaggle/working, 20 GB)\nWORK       = Path(\"/kaggle/working\") if Path(\"/kaggle/working\").is_dir() else Path(\"./work\")\nFIG_DIR    = WORK / \"figures\"        # every figure is saved here for the report\nCACHE_RAW  = WORK / \"cache_raw\"      # minimally processed images  -> V1\nCACHE_PROC = WORK / \"cache_proc\"     # fully preprocessed images   -> V2, V3, V4\nCKPT_DIR   = WORK / \"checkpoints\"\nfor d in (FIG_DIR, CACHE_RAW, CACHE_PROC, CKPT_DIR):\n    d.mkdir(parents=True, exist_ok=True)\n\nprint(f\"Dataset root : {DATA_ROOT}\")\nprint(f\"Train images : {TRAIN_IMG}  ({len(list(TRAIN_IMG.glob('*'))):,} files)\")\nprint(f\"Working dir  : {WORK}\")\n\n# ── Consistent visual identity for every figure in the report ────────────────\nCLASS_NAMES = [\"No DR\", \"Mild\", \"Moderate\", \"Severe\", \"Proliferative\"]\nCLASS_FULL  = {\n    0: \"No DR\",\n    1: \"Mild NPDR\",\n    2: \"Moderate NPDR\",\n    3: \"Severe NPDR\",\n    4: \"Proliferative DR\",\n}\n# Sequential severity ramp: green (healthy) -> red/purple (sight-threatening).\n# Ordered and perceptually increasing, because the classes are ordinal, not categorical.\nSEVERITY_COLORS = [\"#2E7D52\", \"#8CA632\", \"#E0A32E\", \"#D96D28\", \"#A6372E\"]\n# Distinct hues for comparing the four versions (categorical, colour-blind safe).\nVERSION_COLORS  = {\"V1\": \"#5B7FA6\", \"V2\": \"#4E9E7E\", \"V3\": \"#D99B3C\", \"V4\": \"#9B5FA6\"}\n\nplt.rcParams.update({\n    \"figure.dpi\": 110, \"savefig.dpi\": 150, \"savefig.bbox\": \"tight\",\n    \"font.size\": 10, \"axes.titlesize\": 12, \"axes.titleweight\": \"600\",\n    \"axes.labelsize\": 10, \"axes.spines.top\": False, \"axes.spines.right\": False,\n    \"axes.grid\": True, \"grid.alpha\": 0.25, \"grid.linewidth\": 0.6,\n    \"legend.frameon\": False, \"figure.facecolor\": \"white\", \"axes.facecolor\": \"white\",\n})\n\n\ndef savefig(name: str) -> None:\n    \"\"\"Save the current figure into FIG_DIR under a report-ready filename.\"\"\"\n    plt.savefig(FIG_DIR / f\"{name}.png\", dpi=150, bbox_inches=\"tight\", facecolor=\"white\")\n\n\nprint(f\"\\nFigures will be saved to: {FIG_DIR}\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"---\n\n# 1. Problem understanding and dataset justification\n\n> **Rubric criterion 1 — 10 marks.** Clinical significance · dataset relevance · class description ·\n> image distribution · train/validation/test splits · ethical concerns · dataset limitations.\n\n## 1.1 Clinical background\n\nDiabetic retinopathy is a microvascular complication of chronic hyperglycaemia. Persistently\nelevated blood glucose damages the endothelium of the retinal capillaries, which become leaky and\neventually occlude. The resulting ischaemia drives the retina to release VEGF, which in the late\nstage provokes the growth of fragile new vessels that bleed into the vitreous and pull the retina\naway from the choroid.\n\nIt is the leading cause of blindness in working-age adults, and critically, **it is asymptomatic\nuntil it is advanced**. A patient's vision is typically normal through stages 1–3; by the time they\nnotice a problem, irreversible damage has usually occurred. This is precisely why population-level\n*screening* — not symptom-driven presentation — is the standard of care, and why an automated\ntriage tool is clinically valuable: the bottleneck is not treatment, it is the ophthalmologist\nhours needed to read millions of screening photographs.\n\n## 1.2 The five stages and their visual signatures\n\nThe task is graded on the International Clinical Diabetic Retinopathy severity scale. Understanding\nwhat the network must actually *see* explains most of the design decisions later in this notebook:\n\n| Label | Stage | Visual signature on a fundus photograph |\n|---|---|---|\n| **0** | No DR | Healthy retina. No lesions. |\n| **1** | Mild NPDR | **Microaneurysms only** — tiny red dots, often 15–60 µm, a handful of pixels wide. |\n| **2** | Moderate NPDR | Microaneurysms **plus** dot/blot haemorrhages, hard exudates (yellow lipid deposits), cotton-wool spots. |\n| **3** | Severe NPDR | The *4-2-1 rule*: >20 intraretinal haemorrhages in each of 4 quadrants, **or** venous beading in ≥2 quadrants, **or** IRMA in ≥1 quadrant. |\n| **4** | Proliferative DR | **Neovascularisation** and/or vitreous/pre-retinal haemorrhage. Sight-threatening. |\n\nThree consequences follow directly, and each one justifies a design choice made later:\n\n1. **The discriminative lesions are tiny.** A microaneurysm separating class 0 from class 1 may be\n   5–10 pixels across in the original image. Aggressive downsampling destroys the very signal the\n   task depends on — this is the argument for 300×300 rather than the 224×224 that is conventional\n   for ImageNet models, and for contrast enhancement in §2.\n2. **The classes are ordinal, not categorical.** Confusing *Severe* with *Proliferative* is a minor\n   error; confusing *No DR* with *Proliferative* is a catastrophic one. Plain accuracy treats these\n   identically, so §4.2 adopts **quadratic weighted kappa** as the primary metric and §9 analyses\n   adjacent-versus-distant errors separately.\n3. **The boundary between stages is a judgement call, even for experts.** Inter-grader agreement on\n   this scale is itself only moderate, which puts a real ceiling on achievable accuracy and is an\n   important caveat for the report's discussion.\n\n## 1.3 Dataset: APTOS 2019 Blindness Detection\n\nCollected by the **Aravind Eye Hospital** in Tamil Nadu, India. Technicians travelled to rural\nvillages and captured fundus photographs which were later graded by clinicians on the 0–4 scale.\n\n**Why this dataset is appropriate for the brief:** it is sourced from Kaggle as required; it is\nlabelled with the *stage*, not merely presence/absence, which is what the brief asks the model to\npredict; and, unlike the cleaner research datasets, it was captured in real screening conditions\nwith varying cameras, lighting and focus — so a model that works here is closer to one that would\nwork in deployment.","metadata":{}},{"cell_type":"code","source":"# ── Load labels and take a first look ────────────────────────────────────────\ndf = pd.read_csv(TRAIN_CSV)\ndf[\"path\"] = df[\"id_code\"].apply(lambda c: str(TRAIN_IMG / f\"{c}.png\"))\n\n# Sanity check: every referenced file must actually exist before we go further.\nmissing = [p for p in df[\"path\"] if not os.path.exists(p)]\nassert not missing, f\"{len(missing)} referenced images are missing, e.g. {missing[:3]}\"\n\nprint(f\"Total labelled images : {len(df):,}\")\nprint(f\"Columns               : {df.columns.tolist()}\")\nprint(f\"Label range           : {df.diagnosis.min()} to {df.diagnosis.max()}\")\nprint(f\"Duplicate id_codes    : {df.id_code.duplicated().sum()}\")\nprint(f\"Null values           : {df.isnull().sum().sum()}\")\ndf.head()","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Class distribution: the single most important property of this dataset ───\ncounts = df[\"diagnosis\"].value_counts().sort_index()\npcts   = counts / len(df) * 100\n\nfig, (ax1, ax2) = plt.subplots(1, 2, figsize=(13, 4.6),\n                               gridspec_kw={\"width_ratios\": [1.5, 1]})\n\nbars = ax1.bar(range(5), counts.values, color=SEVERITY_COLORS, width=0.68)\nfor i, (c, p) in enumerate(zip(counts.values, pcts.values)):\n    ax1.text(i, c + len(df) * 0.012, f\"{c:,}\\n{p:.1f}%\", ha=\"center\", va=\"bottom\",\n             fontsize=9, linespacing=1.3)\nax1.set_xticks(range(5))\nax1.set_xticklabels([f\"{i}\\n{CLASS_NAMES[i]}\" for i in range(5)], fontsize=9)\nax1.set_ylabel(\"Number of images\")\nax1.set_title(\"Class distribution — severe imbalance\", loc=\"left\")\nax1.set_ylim(0, counts.max() * 1.2)\nax1.grid(axis=\"x\", visible=False)\n\n# Imbalance ratio relative to the rarest class makes the problem concrete.\nimb = counts / counts.min()\nax2.barh(range(5), imb.values, color=SEVERITY_COLORS, height=0.62)\nfor i, v in enumerate(imb.values):\n    ax2.text(v + imb.max() * 0.02, i, f\"{v:.1f}x\", va=\"center\", fontsize=9)\nax2.set_yticks(range(5))\nax2.set_yticklabels(CLASS_NAMES, fontsize=9)\nax2.invert_yaxis()\nax2.set_xlabel(\"Imbalance ratio vs. rarest class\")\nax2.set_title(\"Relative frequency\", loc=\"left\")\nax2.set_xlim(0, imb.max() * 1.18)\nax2.grid(axis=\"y\", visible=False)\n\nplt.tight_layout()\nsavefig(\"01_class_distribution\")\nplt.show()\n\nprint(\"Class distribution\")\nprint(\"-\" * 52)\nfor i in range(5):\n    print(f\"  {i}  {CLASS_FULL[i]:<18s} {counts[i]:>5,d}  ({pcts[i]:5.2f}%)\")\nprint(\"-\" * 52)\nprint(f\"  Majority / minority ratio : {counts.max() / counts.min():.1f} : 1\")\nprint(f\"  Always-predict-'No DR' accuracy baseline : {pcts[0]:.2f}%\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Reading this figure.** Roughly half of all images are healthy, while *Severe NPDR* accounts for\nonly a few percent. Two consequences shape the rest of the project:\n\n- **A trivial classifier that always predicts \"No DR\" scores ~49% accuracy.** Any accuracy figure\n  must be judged against that floor, not against 0. This is the clearest possible argument for\n  reporting macro-F1 and kappa alongside accuracy.\n- **Without intervention the model will under-predict the rare classes** — and those rare classes\n  (Severe, Proliferative) are exactly the clinically urgent ones. Missing them is the costly error.\n  §3.2 addresses this with balanced sampling (§3.2 also documents why loss weighting is\n  deliberately *not* stacked on top of it).","metadata":{}},{"cell_type":"code","source":"# ── Sample images across the severity spectrum ───────────────────────────────\nrng = np.random.RandomState(CFG.seed)\nfig, axes = plt.subplots(3, 5, figsize=(16, 10))\nfor col in range(5):\n    pool = df[df.diagnosis == col].sample(3, random_state=CFG.seed)\n    for row, (_, r) in enumerate(pool.iterrows()):\n        img = cv2.cvtColor(cv2.imread(r[\"path\"]), cv2.COLOR_BGR2RGB)\n        axes[row, col].imshow(img)\n        axes[row, col].axis(\"off\")\n        if row == 0:\n            axes[row, col].set_title(f\"{col} — {CLASS_NAMES[col]}\",\n                                     fontsize=11, color=SEVERITY_COLORS[col], pad=8)\nfig.suptitle(\"Representative fundus images by DR stage\", fontsize=13, y=0.98)\nplt.tight_layout()\nsavefig(\"02_sample_images\")\nplt.show()","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"Note how much **non-pathological** variation there is between these images: field of view, zoom,\ncolour cast, exposure and focus all differ markedly, because they were captured on different\nequipment in different villages. This variation is not signal — it is exactly what the\npreprocessing in §2 is designed to suppress, so the network spends its capacity on lesions rather\nthan on camera characteristics.","metadata":{}},{"cell_type":"code","source":"# ── Quantify the acquisition variability we just observed by eye ─────────────\ndef probe_image(path: str) -> dict:\n    \"\"\"Measure basic acquisition properties of a single fundus image.\n\n    Args:\n        path: Filesystem path to the image.\n\n    Returns:\n        Dict with width, height, aspect ratio, mean brightness, and a\n        variance-of-Laplacian sharpness proxy (higher = sharper).\n    \"\"\"\n    img = cv2.imread(path)\n    h, w = img.shape[:2]\n    gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)\n    return {\n        \"width\": w, \"height\": h, \"aspect\": w / h,\n        \"brightness\": float(gray.mean()),\n        \"sharpness\": float(cv2.Laplacian(gray, cv2.CV_64F).var()),\n    }\n\n# Sample rather than scan all 3,662 images -- enough for a distribution, far faster.\nprobe_df = df.sample(min(400, len(df)), random_state=CFG.seed)\nwith ThreadPoolExecutor(max_workers=8) as pool:\n    stats = pd.DataFrame(list(pool.map(probe_image, probe_df[\"path\"])))\n\nfig, axes = plt.subplots(1, 4, figsize=(17, 3.6))\nspecs = [\n    (\"width\",      \"Image width (px)\",        \"#5B7FA6\"),\n    (\"aspect\",     \"Aspect ratio (w/h)\",      \"#4E9E7E\"),\n    (\"brightness\", \"Mean brightness (0-255)\", \"#D99B3C\"),\n    (\"sharpness\",  \"Sharpness (var. of Laplacian)\", \"#9B5FA6\"),\n]\nfor ax, (col, label, color) in zip(axes, specs):\n    ax.hist(stats[col], bins=34, color=color, alpha=0.85, edgecolor=\"white\", linewidth=0.5)\n    ax.set_xlabel(label); ax.set_ylabel(\"Count\")\n    ax.axvline(stats[col].median(), color=\"#333\", ls=\"--\", lw=1.2)\n    ax.grid(axis=\"x\", visible=False)\n    if col == \"sharpness\":\n        ax.set_xscale(\"log\")\naxes[0].set_title(\"Acquisition variability across a 400-image sample\", loc=\"left\", fontsize=11)\nplt.tight_layout()\nsavefig(\"03_acquisition_variability\")\nplt.show()\n\nprint(stats.describe().round(2).to_string())\nprint(f\"\\nDistinct resolutions in sample : {stats.groupby(['width','height']).ngroups}\")\nprint(\"Takeaway: resolution, framing, exposure and focus all vary substantially.\")\nprint(\"A fixed-size, illumination-normalised input is therefore a precondition, not a nicety.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 1.4 Train / validation / test split\n\nThe rubric asks specifically for **three** splits, and the distinction matters methodologically:\n\n- **Train (70%)** — gradients are computed on this.\n- **Validation (15%)** — used for early stopping and for selecting which epoch's checkpoint to\n  keep. Because we make decisions based on it, performance on validation is *optimistically\n  biased*: by choosing the best of 30 epochs we are fitting to it, mildly.\n- **Test (15%)** — **touched exactly once**, at the very end of each version, after all choices are\n  frozen. This is the only number that honestly estimates generalisation, and it is what the final\n  comparison in §10 reports.\n\nThe split is **stratified** on the label so that all three partitions preserve the (heavily skewed)\nclass proportions — without stratification, a random 15% test split could easily contain almost no\n*Severe* cases, making its metrics meaningless for the class we care most about.\n\nCrucially, the split is computed **once, here**, and all four versions reuse it. Reshuffling between\nversions would confound the ablation: a V3-over-V2 gain could then be an easier test set rather\nthan a better model.","metadata":{}},{"cell_type":"code","source":"# ── One stratified 70/15/15 split, reused by every version ───────────────────\ntrain_df, temp_df = train_test_split(\n    df, test_size=CFG.val_frac + CFG.test_frac,\n    stratify=df[\"diagnosis\"], random_state=CFG.seed,\n)\n# Split the held-out 30% evenly into val and test, again stratified.\nrel_test = CFG.test_frac / (CFG.val_frac + CFG.test_frac)\nval_df, test_df = train_test_split(\n    temp_df, test_size=rel_test,\n    stratify=temp_df[\"diagnosis\"], random_state=CFG.seed,\n)\n\nfor name, d in [(\"train\", train_df), (\"val\", val_df), (\"test\", test_df)]:\n    globals()[f\"{name}_df\"] = d.reset_index(drop=True)\n\nif CFG.fast_dev_run:\n    # Smoke-test mode: shrink every split but keep all 5 classes represented.\n    train_df = train_df.groupby(\"diagnosis\", group_keys=False).head(40).reset_index(drop=True)\n    val_df   = val_df.groupby(\"diagnosis\",   group_keys=False).head(12).reset_index(drop=True)\n    test_df  = test_df.groupby(\"diagnosis\",  group_keys=False).head(12).reset_index(drop=True)\n    print(\">>> FAST_DEV_RUN active: using a reduced subset <<<\\n\")\n\nsplit_tbl = pd.DataFrame({\n    \"Train\": train_df.diagnosis.value_counts().sort_index(),\n    \"Val\":   val_df.diagnosis.value_counts().sort_index(),\n    \"Test\":  test_df.diagnosis.value_counts().sort_index(),\n}).fillna(0).astype(int)\nsplit_tbl.index = [f\"{i} {CLASS_NAMES[i]}\" for i in split_tbl.index]\nsplit_tbl.loc[\"TOTAL\"] = split_tbl.sum()\n\nprint(f\"Train : {len(train_df):>5,d} ({len(train_df)/len(df)*100:.0f}%)\")\nprint(f\"Val   : {len(val_df):>5,d} ({len(val_df)/len(df)*100:.0f}%)\")\nprint(f\"Test  : {len(test_df):>5,d} ({len(test_df)/len(df)*100:.0f}%)\")\nprint()\nprint(split_tbl.to_string())\n\n# Verify stratification actually preserved the proportions.\nprint(\"\\nClass proportions preserved across splits (%):\")\nprop = pd.DataFrame({\n    \"Full\":  df.diagnosis.value_counts(normalize=True).sort_index() * 100,\n    \"Train\": train_df.diagnosis.value_counts(normalize=True).sort_index() * 100,\n    \"Val\":   val_df.diagnosis.value_counts(normalize=True).sort_index() * 100,\n    \"Test\":  test_df.diagnosis.value_counts(normalize=True).sort_index() * 100,\n}).round(2)\nprop.index = CLASS_NAMES\nprint(prop.to_string())","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 1.5 Ethical considerations and dataset limitations\n\nA model that reads medical images carries obligations that a general image classifier does not.\nStating these honestly is part of the engineering, not a disclaimer bolted on afterwards.\n\n**Population and generalisation.** Every image comes from a single provider in Tamil Nadu, India.\nRetinal pigmentation varies systematically with ethnicity, and fundus appearance varies with the\ncamera model used. A network trained only on this distribution may degrade on a European, East\nAsian or African screening population, or on a different camera — and it would degrade *silently*,\nreturning confident predictions that happen to be wrong. External validation on an independent\ncohort (EyePACS, Messidor-2, IDRiD) would be a precondition for any real deployment.\n\n**No demographic metadata.** The dataset provides no age, sex, ethnicity or diabetes-duration\nfields. This is good for privacy but it makes it **impossible to audit the model for performance\ndisparities between subgroups** — we cannot demonstrate the model is equitable, only that we have\nnot checked. That is a genuine limitation worth stating plainly rather than glossing.\n\n**Label noise and the accuracy ceiling.** Grades appear to come from a single clinician per image,\nwith no adjudication. Published inter-grader agreement on this 5-point scale is only moderate\n(kappa typically ~0.6–0.8 between experts). Some fraction of our \"errors\" are therefore\ndisagreements with one grader's judgement rather than true mistakes, and no model can exceed the\nconsistency of its labels. Reported accuracy should be read with that ceiling in mind.\n\n**Error asymmetry.** A false negative on Proliferative DR may mean a patient loses their sight; a\nfalse positive means an unnecessary referral, which costs money and anxiety but no vision. These\nare not equivalent, yet cross-entropy treats them as such. §3.3 partially addresses this through\nclass weighting, and §12 discusses operating-point selection (favouring recall on severe classes)\nas the appropriate deployment-time response.\n\n**Intended use.** This system is a **screening triage aid**, not a diagnostic device. The\ndefensible deployment is one that prioritises a queue for human review and flags urgent cases\nsooner — never one that discharges a patient without a clinician. Any autonomous use would require\nprospective clinical trials and regulatory clearance (CE marking / FDA). There is also a documented\nrisk of **automation bias**: clinicians shown a confident \"No DR\" output become measurably less\nlikely to spot pathology themselves, so the interface matters as much as the model, which is part\nof the motivation for the Grad-CAM explanations in §8.","metadata":{}},{"cell_type":"markdown","source":"---\n\n# 2. Data preprocessing\n\n> **Rubric criterion 2 — 10 marks.** Comprehensive, justified, reproducible pipeline that measurably\n> improves image quality: contrast enhancement, resizing, normalisation, noise removal, edge\n> enhancement.\n\n## 2.1 Design rationale\n\n§1.3 established the problem: lesions are small and low-contrast, while *non-pathological*\nvariation (illumination, framing, colour cast) is large. Preprocessing exists to invert that ratio.\nFive stages, each addressing a specific defect, applied in this order:\n\n**1 · Crop uninformative black border.** Fundus photographs are a circular image on a black\nrectangle. Typically 30–40% of the pixels carry no information at all. Cropping to the retinal\nbounding box means the subsequent resize spends its pixel budget on retina rather than on padding —\neffectively a free gain in resolution.\n\n**2 · Circular mask.** The extreme edge of the circle has heavy vignetting and sensor artefacts.\nMasking to 97% of the radius removes a hard, high-contrast rim that would otherwise be the\nstrongest edge in the image and would dominate early convolutional filters.\n\n**3 · Ben Graham illumination correction.** The central technique, from the winning solution to the\n2015 Kaggle DR competition. The image is blurred heavily to estimate the *local* background\nillumination, and that estimate is subtracted:\n\n$$I_{\\text{out}} = 4 \\cdot I - 4 \\cdot G_{\\sigma}(I) + 128$$\n\nBecause $G_\\sigma(I)$ captures slowly varying illumination and colour cast, subtracting it removes\nexactly the inter-camera variation we want gone, while the ×4 gain amplifies the high-frequency\nresidual — which is precisely where microaneurysms and haemorrhages live. One operation both\nnormalises nuisance variation *and* enhances the signal.\n\n**4 · CLAHE on the luminance channel.** Contrast Limited Adaptive Histogram Equalisation equalises\ncontrast in local tiles rather than globally, so a bright optic disc cannot wash out contrast in a\ndim peripheral quadrant. It is applied to the **L channel in LAB space**, not to RGB: equalising\nRGB channels independently shifts hue, and in fundus imaging colour is diagnostic (hard exudates\nare yellow, haemorrhages red). The `clipLimit` caps amplification so that flat regions do not turn\ninto amplified sensor noise.\n\n**5 · Mild denoising and optional unsharp masking.** Ben Graham's ×4 gain amplifies noise along\nwith signal, so a light bilateral filter — edge-preserving, unlike a Gaussian — suppresses grain\nwhile leaving vessel boundaries intact. An unsharp mask is provided as the explicit **edge\nenhancement** step named in the brief; §2.4 examines whether it actually helps.\n\n**Reproducibility.** Preprocessing is deterministic (no randomness) and is executed **once** into an\non-disk cache. Every version then reads the same cached pixels, so results cannot drift because of\na preprocessing change between runs, and training is far faster since full-resolution JPEG decoding\nhappens once instead of every epoch.","metadata":{}},{"cell_type":"code","source":"# ═════════════════════════════════════════════════════════════════════════════\n# PREPROCESSING PRIMITIVES\n# Each function takes and returns an RGB uint8 array, so they compose freely\n# and any single stage can be ablated by removing one line.\n# ═════════════════════════════════════════════════════════════════════════════\n\ndef crop_black_border(img: np.ndarray, tol: int = 7) -> np.ndarray:\n    \"\"\"Crop the uninformative black frame surrounding the circular fundus.\n\n    Args:\n        img: RGB uint8 image.\n        tol: Grey-level below which a pixel counts as background. 7 tolerates\n            mild sensor noise in the \"black\" region without eating dark retina.\n\n    Returns:\n        The image cropped to the bounding box of non-background pixels. Returns\n        the input unchanged if the image is (almost) entirely dark, which would\n        otherwise produce a degenerate zero-size crop.\n    \"\"\"\n    gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY)\n    mask = gray > tol\n    if mask.sum() < 100:\n        return img\n    coords = np.argwhere(mask)\n    y0, x0 = coords.min(axis=0)\n    y1, x1 = coords.max(axis=0) + 1\n    return img[y0:y1, x0:x1]\n\n\ndef apply_circle_mask(img: np.ndarray, frac: float = 0.97) -> np.ndarray:\n    \"\"\"Zero out everything outside the inscribed circle.\n\n    Removes the vignetted rim, which is the highest-contrast edge in an\n    uncorrected fundus image and would otherwise dominate early conv filters.\n\n    Args:\n        img: RGB uint8 image.\n        frac: Radius as a fraction of half the shorter side.\n\n    Returns:\n        Masked RGB uint8 image.\n    \"\"\"\n    h, w = img.shape[:2]\n    mask = np.zeros((h, w), dtype=np.uint8)\n    cv2.circle(mask, (w // 2, h // 2), int(min(h, w) / 2 * frac), 1, thickness=-1)\n    return img * mask[..., None]\n\n\ndef ben_graham(img: np.ndarray, sigma_frac: float = 30.0, gain: float = 4.0) -> np.ndarray:\n    \"\"\"Subtract local average colour to correct illumination and boost lesions.\n\n    out = gain*I - gain*Gaussian(I) + 128\n\n    Args:\n        img: RGB uint8 image.\n        sigma_frac: Blur sigma expressed as image_width / sigma_frac, so the\n            correction is scale-invariant across differently sized inputs.\n        gain: Amplification applied to the high-frequency residual.\n\n    Returns:\n        Illumination-corrected RGB uint8 image, centred around mid-grey.\n    \"\"\"\n    sigma = max(img.shape[1] / sigma_frac, 1.0)\n    blurred = cv2.GaussianBlur(img, (0, 0), sigma)\n    return cv2.addWeighted(img, gain, blurred, -gain, 128)\n\n\ndef apply_clahe(img: np.ndarray, clip: float = 2.0, grid: int = 8) -> np.ndarray:\n    \"\"\"Adaptive local contrast enhancement on the LAB luminance channel only.\n\n    Operating on L rather than on R/G/B independently preserves hue, which\n    matters because lesion colour is diagnostic (yellow exudates, red\n    haemorrhages).\n\n    Args:\n        img: RGB uint8 image.\n        clip: Contrast clip limit; caps noise amplification in flat regions.\n        grid: Tile grid size for local equalisation.\n\n    Returns:\n        Contrast-enhanced RGB uint8 image.\n    \"\"\"\n    lab = cv2.cvtColor(img, cv2.COLOR_RGB2LAB)\n    l, a, b = cv2.split(lab)\n    l = cv2.createCLAHE(clipLimit=clip, tileGridSize=(grid, grid)).apply(l)\n    return cv2.cvtColor(cv2.merge([l, a, b]), cv2.COLOR_LAB2RGB)\n\n\ndef denoise(img: np.ndarray, d: int = 5, sigma_color: int = 40, sigma_space: int = 40) -> np.ndarray:\n    \"\"\"Edge-preserving bilateral denoising.\n\n    Chosen over a Gaussian because Ben Graham amplification raises grain, and a\n    Gaussian would blur the thin vessel edges we are trying to preserve.\n    \"\"\"\n    return cv2.bilateralFilter(img, d, sigma_color, sigma_space)\n\n\ndef unsharp_mask(img: np.ndarray, sigma: float = 1.5, amount: float = 0.6) -> np.ndarray:\n    \"\"\"Explicit edge enhancement: add back a scaled high-frequency residual.\n\n    Args:\n        img: RGB uint8 image.\n        sigma: Gaussian sigma defining \"low frequency\".\n        amount: Strength of the sharpening. Above ~1.0 this produces visible halos.\n    \"\"\"\n    blurred = cv2.GaussianBlur(img, (0, 0), sigma)\n    return cv2.addWeighted(img, 1 + amount, blurred, -amount, 0)\n\n\ndef preprocess_fundus(img: np.ndarray, size: int, full: bool = True) -> np.ndarray:\n    \"\"\"Full deterministic preprocessing pipeline for one fundus image.\n\n    Args:\n        img: RGB uint8 image at original resolution.\n        size: Output edge length in pixels (square).\n        full: If False, only crop + resize are applied. This is the V1\n            \"minimal preprocessing\" path, and keeping both behaviours in one\n            function guarantees V1 and V2 differ by exactly this flag.\n\n    Returns:\n        RGB uint8 image of shape (size, size, 3).\n    \"\"\"\n    img = crop_black_border(img)\n    if not full:\n        return cv2.resize(img, (size, size), interpolation=cv2.INTER_AREA)\n\n    img = cv2.resize(img, (size, size), interpolation=cv2.INTER_AREA)\n    img = ben_graham(img)\n    img = apply_clahe(img)\n    img = denoise(img)\n    img = apply_circle_mask(img)\n    return img\n\n\ndef load_rgb(path: str) -> np.ndarray:\n    \"\"\"Read an image from disk as RGB uint8 (OpenCV loads BGR by default).\"\"\"\n    img = cv2.imread(path, cv2.IMREAD_COLOR)\n    if img is None:\n        raise IOError(f\"Failed to read image: {path}\")\n    return cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n\n\nprint(\"Preprocessing primitives defined.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 2.2 Visual evidence: stage-by-stage\n\nThe figure below applies the pipeline one stage at a time to the same image, so each operation's\nindividual contribution is visible rather than only the combined result.","metadata":{}},{"cell_type":"code","source":"# ── Walk one image through the pipeline, one stage at a time ─────────────────\ndemo_path = df[df.diagnosis == 2].iloc[3][\"path\"]\norig = load_rgb(demo_path)\n\ns1 = crop_black_border(orig)\ns2 = cv2.resize(s1, (CFG.img_size, CFG.img_size), interpolation=cv2.INTER_AREA)\ns3 = ben_graham(s2)\ns4 = apply_clahe(s3)\ns5 = denoise(s4)\ns6 = apply_circle_mask(s5)\n\nstages = [\n    (orig, \"0 · Original\", f\"{orig.shape[1]}x{orig.shape[0]} px\"),\n    (s1,   \"1 · Border cropped\", f\"{s1.shape[1]}x{s1.shape[0]} px\"),\n    (s2,   \"2 · Resized\", f\"{CFG.img_size}x{CFG.img_size} px\"),\n    (s3,   \"3 · Ben Graham\", \"illumination corrected\"),\n    (s4,   \"4 · CLAHE\", \"local contrast\"),\n    (s6,   \"5 · Denoise + mask\", \"final\"),\n]\nfig, axes = plt.subplots(1, 6, figsize=(19, 3.6))\nfor ax, (im, title, sub) in zip(axes, stages):\n    ax.imshow(im); ax.axis(\"off\")\n    ax.set_title(title, fontsize=10.5, pad=6)\n    ax.text(0.5, -0.06, sub, transform=ax.transAxes, ha=\"center\",\n            fontsize=8.5, color=\"#666\")\nfig.suptitle(\"Preprocessing pipeline, stage by stage (Moderate NPDR case)\", fontsize=13, y=1.06)\nplt.tight_layout()\nsavefig(\"04_preprocessing_stages\")\nplt.show()\n\npix_before = orig.shape[0] * orig.shape[1]\npix_after  = s1.shape[0] * s1.shape[1]\nprint(f\"Border crop discarded {(1 - pix_after/pix_before)*100:.1f}% of pixels as pure background,\")\nprint(\"so the subsequent resize allocates its pixel budget to retina instead of padding.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Does preprocessing measurably improve image quality? Quantify it. ────────\n# A claim of \"improves image quality significantly\" should be backed by numbers,\n# not only by eye. We measure contrast (std of luminance) and edge energy\n# (mean Sobel magnitude) before and after, over a sample.\ndef quality_metrics(img: np.ndarray) -> tuple:\n    \"\"\"Return (RMS contrast, mean Sobel edge magnitude) for an RGB image.\"\"\"\n    gray = cv2.cvtColor(img, cv2.COLOR_RGB2GRAY).astype(np.float32)\n    gx = cv2.Sobel(gray, cv2.CV_32F, 1, 0, ksize=3)\n    gy = cv2.Sobel(gray, cv2.CV_32F, 0, 1, ksize=3)\n    return float(gray.std()), float(np.sqrt(gx**2 + gy**2).mean())\n\nsample_paths = df.sample(120, random_state=CFG.seed)[\"path\"].tolist()\n\ndef _compare(p):\n    img = load_rgb(p)\n    return quality_metrics(preprocess_fundus(img, CFG.img_size, full=False)) + \\\n           quality_metrics(preprocess_fundus(img, CFG.img_size, full=True))\n\nwith ThreadPoolExecutor(max_workers=8) as pool:\n    qm = np.array(list(pool.map(_compare, sample_paths)))\n\nlabels_q = [\"RMS contrast\", \"Edge energy (Sobel)\"]\nfig, axes = plt.subplots(1, 2, figsize=(11, 3.8))\nfor i, ax in enumerate(axes):\n    before, after = qm[:, i], qm[:, i + 2]\n    ax.hist(before, bins=28, alpha=0.72, label=\"Minimal (V1)\", color=\"#9AA5B1\", edgecolor=\"white\")\n    ax.hist(after,  bins=28, alpha=0.72, label=\"Full pipeline (V2+)\", color=\"#4E9E7E\", edgecolor=\"white\")\n    ax.set_xlabel(labels_q[i]); ax.set_ylabel(\"Count\"); ax.legend()\n    ax.grid(axis=\"x\", visible=False)\n    gain = after.mean() / before.mean()\n    ax.set_title(f\"{labels_q[i]}  —  {gain:.2f}x mean increase\", loc=\"left\", fontsize=11)\nfig.suptitle(\"Preprocessing measurably increases contrast and edge energy\", fontsize=12.5, y=1.04)\nplt.tight_layout()\nsavefig(\"05_preprocessing_quality_metrics\")\nplt.show()\n\nprint(f\"RMS contrast : {qm[:,0].mean():6.2f} -> {qm[:,2].mean():6.2f}  ({qm[:,2].mean()/qm[:,0].mean():.2f}x)\")\nprint(f\"Edge energy  : {qm[:,1].mean():6.2f} -> {qm[:,3].mean():6.2f}  ({qm[:,3].mean()/qm[:,1].mean():.2f}x)\")\nprint(\"\\nHigher edge energy means lesion and vessel boundaries carry more of the\")\nprint(\"signal the convolutional filters will respond to.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Effect across the severity spectrum, plus explicit edge extraction ───────\nfig, axes = plt.subplots(3, 5, figsize=(16, 10))\nfor col in range(5):\n    p = df[df.diagnosis == col].iloc[2][\"path\"]\n    raw = load_rgb(p)\n    mini = preprocess_fundus(raw, CFG.img_size, full=False)\n    full = preprocess_fundus(raw, CFG.img_size, full=True)\n    # Edge map derived from the processed image: evidence of what step 5 exposes.\n    g = cv2.cvtColor(full, cv2.COLOR_RGB2GRAY)\n    edges = cv2.magnitude(cv2.Sobel(g, cv2.CV_32F, 1, 0, ksize=3),\n                          cv2.Sobel(g, cv2.CV_32F, 0, 1, ksize=3))\n    edges = np.clip(edges / (edges.max() + 1e-6) * 255, 0, 255).astype(np.uint8)\n\n    for row, (im, cmap) in enumerate([(mini, None), (full, None), (edges, \"magma\")]):\n        axes[row, col].imshow(im, cmap=cmap)\n        axes[row, col].axis(\"off\")\n    axes[0, col].set_title(f\"{col} — {CLASS_NAMES[col]}\", fontsize=10.5,\n                           color=SEVERITY_COLORS[col], pad=6)\n\nfor row, lab in enumerate([\"Minimal\\n(V1 input)\", \"Full pipeline\\n(V2+ input)\", \"Sobel edge map\\n(vascular structure)\"]):\n    axes[row, 0].text(-0.14, 0.5, lab, transform=axes[row, 0].transAxes,\n                      rotation=90, va=\"center\", ha=\"center\", fontsize=10, color=\"#444\")\nfig.suptitle(\"Preprocessing effect across the severity spectrum\", fontsize=13, y=0.98)\nplt.tight_layout()\nsavefig(\"06_preprocessing_by_class\")\nplt.show()","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**What the figure shows.** In the middle row the inter-image colour cast has largely disappeared —\nall five images now sit on a comparable grey-green background, so the network can no longer use\n\"this image is orange\" as a shortcut feature. Vessels and lesions stand out sharply against that\nflat background. The bottom row makes the point explicit: after preprocessing, a simple Sobel\noperator recovers a clean vascular tree, which indicates the high-frequency structure the\nconvolutional filters need is genuinely present and well separated from the background.\n\n### A note on unsharp masking\n\nThe brief names *edge enhancement* explicitly. Ben Graham amplification is already a form of\nhigh-frequency boost, so an additional unsharp mask risks over-sharpening into halo artefacts. The\ncell below compares the options so the choice is evidence-based rather than assumed.","metadata":{}},{"cell_type":"code","source":"# ── Is an extra unsharp mask on top of Ben Graham helpful or harmful? ────────\np = df[df.diagnosis == 3].iloc[1][\"path\"]\nbase = preprocess_fundus(load_rgb(p), CFG.img_size, full=True)\nvariants = [\n    (base,                                   \"No unsharp (chosen)\"),\n    (unsharp_mask(base, 1.5, 0.6),           \"Unsharp  amount=0.6\"),\n    (unsharp_mask(base, 1.5, 1.4),           \"Unsharp  amount=1.4\"),\n]\nfig, axes = plt.subplots(2, 3, figsize=(12.5, 8.4))\nfor i, (im, title) in enumerate(variants):\n    axes[0, i].imshow(im); axes[0, i].axis(\"off\"); axes[0, i].set_title(title, fontsize=10.5)\n    # Zoom into the centre to make halo artefacts visible.\n    c = CFG.img_size // 2; r = 55\n    axes[1, i].imshow(im[c-r:c+r, c-r:c+r]); axes[1, i].axis(\"off\")\n    _, e = quality_metrics(im)\n    axes[1, i].set_title(f\"centre crop · edge energy {e:.1f}\", fontsize=9.5, color=\"#666\")\nfig.suptitle(\"Edge enhancement study: diminishing returns above moderate amounts\", fontsize=12.5, y=0.97)\nplt.tight_layout()\nsavefig(\"07_edge_enhancement_study\")\nplt.show()\n\nprint(\"Decision: Ben Graham already supplies strong high-frequency gain. An additional\")\nprint(\"unsharp mask raises edge-energy numbers but introduces visible halos around\")\nprint(\"vessels at amount>=1.4, which are artefacts the network could overfit to.\")\nprint(\"The pipeline therefore omits it; `unsharp_mask` is retained for this ablation.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 2.3 Building the preprocessed cache\n\nPreprocessing every image on every epoch would dominate runtime — the original files are large and\nJPEG decoding is expensive. Instead each image is processed **once** and written to disk at\n`cache_size` (380 px, slightly larger than the 300 px model input so random-crop augmentation has\nheadroom to work with).\n\nTwo caches are built, and this is the mechanism that makes the V1→V2 comparison rigorous:\n\n- `cache_raw/`  — crop + resize only → the **V1** input\n- `cache_proc/` — the full pipeline   → the **V2/V3/V4** input\n\nBoth derive from the same source images through the same function, differing only in the `full`\nflag. The V1→V2 delta is therefore attributable to preprocessing and augmentation alone.","metadata":{}},{"cell_type":"code","source":"def build_cache(frame: pd.DataFrame, out_dir: Path, full: bool, desc: str) -> None:\n    \"\"\"Preprocess every image in `frame` once and write it to `out_dir` as PNG.\n\n    Idempotent: images already present are skipped, so the cell is safe to\n    re-run after an interrupted session.\n\n    Args:\n        frame: DataFrame with `id_code` and `path` columns.\n        out_dir: Destination directory.\n        full: Passed through to `preprocess_fundus`; selects minimal vs. full pipeline.\n        desc: Label used in progress output.\n    \"\"\"\n    todo = [(r.id_code, r.path) for r in frame.itertuples()\n            if not (out_dir / f\"{r.id_code}.png\").exists()]\n    if not todo:\n        print(f\"{desc}: cache already complete ({len(frame):,} images).\")\n        return\n\n    def _one(item):\n        cid, src = item\n        img = preprocess_fundus(load_rgb(src), CFG.cache_size, full=full)\n        cv2.imwrite(str(out_dir / f\"{cid}.png\"), cv2.cvtColor(img, cv2.COLOR_RGB2BGR))\n\n    t0 = time.time()\n    with ThreadPoolExecutor(max_workers=8) as pool:\n        for i, _ in enumerate(pool.map(_one, todo), 1):\n            if i % 500 == 0 or i == len(todo):\n                el = time.time() - t0\n                print(f\"  {desc}: {i:,}/{len(todo):,}  ({el:.0f}s elapsed, \"\n                      f\"{el/i*(len(todo)-i):.0f}s remaining)\")\n    print(f\"{desc}: done in {time.time()-t0:.0f}s\")\n\n\n# Only cache the images actually used by the (possibly reduced) splits.\nused = pd.concat([train_df, val_df, test_df])[[\"id_code\", \"path\"]].drop_duplicates()\nprint(f\"Caching {len(used):,} images at {CFG.cache_size}px into two variants...\\n\")\nbuild_cache(used, CACHE_RAW,  full=False, desc=\"raw cache  (V1)\")\nbuild_cache(used, CACHE_PROC, full=True,  desc=\"proc cache (V2+)\")\n\ndef _dir_mb(d: Path) -> float:\n    return sum(f.stat().st_size for f in d.glob(\"*.png\")) / 1e6\n\nprint(f\"\\ncache_raw  : {len(list(CACHE_RAW.glob('*.png'))):,} files, {_dir_mb(CACHE_RAW):.0f} MB\")\nprint(f\"cache_proc : {len(list(CACHE_PROC.glob('*.png'))):,} files, {_dir_mb(CACHE_PROC):.0f} MB\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"---\n\n# 3. Data augmentation and dataset balancing\n\n> **Rubric criterion 3 — 10 marks.** Effective, justified augmentation (rotation, flipping, zooming,\n> brightness) **and** appropriate handling of class imbalance.\n\n## 3.1 Which augmentations, and why these\n\nWith ~2,560 training images and models holding millions of parameters, overfitting is the dominant\nrisk. Augmentation addresses it by synthesising plausible variants — but only transformations that\nproduce images which *could genuinely have come from the scanner* are valid. An implausible\ntransformation teaches the network to be invariant to something that never occurs, wasting\ncapacity.\n\n| Transform | Setting | Justification |\n|---|---|---|\n| **Horizontal + vertical flip** | p = 0.5 each | A fundus has no canonical orientation; left and right eyes are near-mirror images, and the camera may be rotated. Fully label-preserving. |\n| **Rotation** | ±180° | Follows from the same argument — there is no \"up\" in a fundus photograph, so the full rotation group is valid. This is unusually aggressive compared with natural-image tasks, and it is justified here precisely because of that symmetry. |\n| **Random resized crop** | scale 0.85–1.0 | Simulates variation in zoom and centring between operators. Kept mild: cropping too hard can cut peripheral lesions out of frame and silently corrupt the label. |\n| **Brightness / contrast jitter** | ±20% | Mimics residual exposure variation. Kept modest because Ben Graham has already normalised most illumination — over-jittering would undo §2's work. |\n| **Saturation / hue jitter** | small | Covers camera colour-response differences. Hue is deliberately constrained (±0.02) since lesion colour is diagnostic. |\n| **Random erasing** | p = 0.25 | Occludes a small patch, forcing the network to aggregate evidence from the whole retina rather than fixating on one region — a useful prior given lesions are distributed. |\n\n**Deliberately excluded:** large shears and perspective warps (they distort vessel geometry, which\nis itself diagnostic in Severe NPDR via venous beading) and grayscale conversion (it would discard\nthe colour information distinguishing exudates from haemorrhages).\n\nValidation and test sets receive **no** augmentation beyond deterministic resize and normalisation —\nevaluation must measure performance on the real distribution, not an augmented one.","metadata":{}},{"cell_type":"code","source":"# ═════════════════════════════════════════════════════════════════════════════\n# AUGMENTATION PIPELINES\n# torchvision.transforms is used rather than albumentations: it ships with\n# PyTorch (no extra dependency, no version-skew risk on Kaggle) and covers\n# every transform the rubric asks for.\n# ═════════════════════════════════════════════════════════════════════════════\n\n# ImageNet channel statistics. Used even for the from-scratch model so that all\n# four versions share identical input normalisation -- another confound removed.\nIMAGENET_MEAN = [0.485, 0.456, 0.406]\nIMAGENET_STD  = [0.229, 0.224, 0.225]\n\n\ndef build_transforms(img_size: int, train: bool, augment: bool):\n    \"\"\"Construct the transform pipeline for one dataset split.\n\n    Args:\n        img_size: Target square edge length.\n        train: Whether this is the training split.\n        augment: Whether to apply stochastic augmentation. Separated from\n            `train` so V1 can use a training split with augmentation disabled,\n            which is exactly the V1-vs-V2 ablation.\n\n    Returns:\n        A torchvision transform callable mapping PIL.Image -> normalised tensor.\n    \"\"\"\n    norm = [T.ToTensor(), T.Normalize(IMAGENET_MEAN, IMAGENET_STD)]\n\n    if not (train and augment):\n        # Deterministic path: validation, test, and the V1 training set.\n        return T.Compose([T.Resize((img_size, img_size))] + norm)\n\n    return T.Compose([\n        # Geometric -- exploits the rotational symmetry of fundus images.\n        T.RandomResizedCrop(img_size, scale=(0.85, 1.0), ratio=(0.95, 1.05)),\n        T.RandomHorizontalFlip(p=0.5),\n        T.RandomVerticalFlip(p=0.5),\n        T.RandomRotation(degrees=180, interpolation=T.InterpolationMode.BILINEAR),\n        # Photometric -- mild, because preprocessing already normalised illumination.\n        T.ColorJitter(brightness=0.2, contrast=0.2, saturation=0.15, hue=0.02),\n        *norm,\n        # Occlusion -- forces distributed rather than localised evidence.\n        T.RandomErasing(p=0.25, scale=(0.02, 0.10), ratio=(0.3, 3.3), value=0),\n    ])\n\n\nprint(\"Transform builders defined.\")\nprint(\"  V1  : train=True,  augment=False  -> deterministic resize only\")\nprint(\"  V2+ : train=True,  augment=True   -> full augmentation stack\")\nprint(\"  eval: train=False                 -> deterministic, always\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Visualise what the network actually sees during training ─────────────────\ndef denormalise(t: torch.Tensor) -> np.ndarray:\n    \"\"\"Invert ImageNet normalisation so a tensor can be displayed.\"\"\"\n    m = torch.tensor(IMAGENET_MEAN).view(3, 1, 1)\n    s = torch.tensor(IMAGENET_STD).view(3, 1, 1)\n    return (t * s + m).clamp(0, 1).permute(1, 2, 0).numpy()\n\nsample_id = train_df[train_df.diagnosis == 2].iloc[0][\"id_code\"]\npil = Image.open(CACHE_PROC / f\"{sample_id}.png\").convert(\"RGB\")\n\naug_tf = build_transforms(CFG.img_size, train=True, augment=True)\ndet_tf = build_transforms(CFG.img_size, train=False, augment=False)\n\nfig, axes = plt.subplots(2, 6, figsize=(17, 6))\naxes[0, 0].imshow(denormalise(det_tf(pil)))\naxes[0, 0].set_title(\"Original\\n(no augmentation)\", fontsize=10, color=\"#A6372E\")\naxes[0, 0].axis(\"off\")\nfor i in range(1, 12):\n    ax = axes[i // 6, i % 6]\n    torch.manual_seed(CFG.seed + i)\n    ax.imshow(denormalise(aug_tf(pil)))\n    ax.set_title(f\"Augmented #{i}\", fontsize=9.5)\n    ax.axis(\"off\")\nfig.suptitle(\"Twelve stochastic augmentations of a single training image\", fontsize=13, y=0.99)\nplt.tight_layout()\nsavefig(\"08_augmentation_samples\")\nplt.show()\n\nprint(\"Each variant is a plausible photograph of the same eye: the label is preserved,\")\nprint(\"while rotation, framing, exposure and occlusion all vary.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 3.2 Handling class imbalance\n\n§1.3 established a ~10:1 imbalance. Left untreated, two things happen: the loss is dominated by the\nmajority class, and the model converges to a degenerate solution that predicts *No DR* for almost\neverything — achieving ~49% accuracy while being clinically worthless.\n\nThree mechanisms are standard, and they are not equivalent:\n\n1. **Class-weighted loss** — scale each sample's loss by $w_c \\propto 1/n_c$, so a misclassified\n   *Severe* case costs as much as ten misclassified *No DR* cases. Cheap, but gradients from rare\n   classes stay noisy because those samples are still rare within a batch.\n2. **Balanced sampling** (`WeightedRandomSampler`) — oversample rare classes so each *batch* is\n   roughly class-balanced. Rare-class gradients become stable, at the cost of showing the same rare\n   images repeatedly, which raises overfitting risk on those classes.\n3. **Synthetic oversampling (SMOTE)** — interpolating between images in pixel space produces\n   anatomically meaningless blends, so it is **not** appropriate here. Augmentation already acts as\n   a principled, domain-valid oversampler.\n\n**This project uses mechanism 2 (balanced sampling) + augmentation.** Balanced sampling combined\nwith strong augmentation means a rare image is revisited often but almost never in the same form,\nwhich captures the benefit of oversampling while augmentation supplies the diversity that mitigates\nits main drawback. Note that `WeightedRandomSampler` with `replacement=True` makes an \"epoch\" a\nfixed number of draws rather than a pass over every image — that is expected, and the epoch length\nis set to `len(train_df)` so epoch counts stay comparable across versions.\n\n> **A pitfall worth naming explicitly: do not stack mechanisms 1 and 2.** An earlier iteration of\n> this pipeline applied the class-weighted loss *on top of* the balanced sampler — both corrections\n> at once. That combination over-corrects badly. The sampler already equalises each class's\n> contribution to the gradient; multiplying by an inverse-frequency loss weight on top of that\n> pushes the correction far past where it should stop. In a diagnostic run, this specifically\n> collapsed recall on the **Moderate** class — which is not even rare (it is the second-largest\n> class in the dataset) — to near zero, because its near-average loss weight was dwarfed by the\n> very large weights assigned to the genuinely rare classes once sampling had already balanced the\n> batch composition. The symptom was a large drop in accuracy from V1 to V2/V3/V4 that was not\n> matched by a correspondingly large gain in macro-F1 — a sign that the correction was actively\n> hurting a class it had no reason to touch, not simply trading accuracy for fairness.\n>\n> The fix, used throughout this notebook, is to let **one** mechanism carry the full correction:\n> the sampler rebalances what the model sees, and the loss function is left unweighted. This is a\n> useful case study for the report's error-analysis section — it demonstrates that imbalance\n> handling is not \"more correction is always better,\" and that per-class metrics, not just\n> macro-F1, are what catch this kind of regression.","metadata":{}},{"cell_type":"code","source":"def compute_class_weights(labels: np.ndarray, num_classes: int = 5,\n                          scheme: str = \"inverse\") -> torch.Tensor:\n    \"\"\"Per-class loss weights for imbalanced training.\n\n    Args:\n        labels: 1-D array of integer labels.\n        num_classes: Total number of classes.\n        scheme: 'inverse' for w_c = N/(K*n_c) (full correction), or 'sqrt_inverse'\n            for a gentler sqrt-damped version that is less prone to over-correcting\n            on extremely rare classes.\n\n    Returns:\n        Float tensor of shape (num_classes,), normalised to mean 1 so the overall\n        loss scale -- and therefore the effective learning rate -- is unchanged.\n    \"\"\"\n    counts = np.bincount(labels, minlength=num_classes).astype(np.float64)\n    counts = np.maximum(counts, 1)                       # guard against empty classes\n    w = len(labels) / (num_classes * counts)\n    if scheme == \"sqrt_inverse\":\n        w = np.sqrt(w)\n    w = w / w.mean()\n    return torch.tensor(w, dtype=torch.float32)\n\n\ndef make_balanced_sampler(labels: np.ndarray, num_classes: int = 5) -> WeightedRandomSampler:\n    \"\"\"Sampler that draws each class with roughly equal probability.\n\n    Args:\n        labels: 1-D array of integer labels for the training split.\n        num_classes: Total number of classes.\n\n    Returns:\n        A WeightedRandomSampler drawing `len(labels)` samples with replacement,\n        so one epoch remains the same nominal size as an unbalanced epoch.\n    \"\"\"\n    counts = np.bincount(labels, minlength=num_classes).astype(np.float64)\n    counts = np.maximum(counts, 1)\n    per_sample_w = (1.0 / counts)[labels]\n    return WeightedRandomSampler(\n        weights=torch.as_tensor(per_sample_w, dtype=torch.double),\n        num_samples=len(labels),\n        replacement=True,\n    )\n\n\ntrain_labels  = train_df[\"diagnosis\"].values\nCLASS_WEIGHTS = compute_class_weights(train_labels, CFG.num_classes).to(DEVICE)\n\nprint(\"Class weights (inverse frequency, mean-normalised):\")\nfor i, w in enumerate(CLASS_WEIGHTS.cpu().numpy()):\n    n = (train_labels == i).sum()\n    print(f\"  {i} {CLASS_NAMES[i]:<14s} n={n:>4d}   weight={w:.3f}\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Verify the sampler actually rebalances the batches it produces ───────────\nsampler = make_balanced_sampler(train_labels, CFG.num_classes)\ndrawn = np.array([train_labels[i] for i in list(sampler)])\n\norig_pct    = np.bincount(train_labels, minlength=5) / len(train_labels) * 100\nsampled_pct = np.bincount(drawn,        minlength=5) / len(drawn)        * 100\n\nx = np.arange(5); w = 0.38\nfig, ax = plt.subplots(figsize=(9, 4))\nax.bar(x - w/2, orig_pct,    w, label=\"Original distribution\", color=\"#9AA5B1\")\nax.bar(x + w/2, sampled_pct, w, label=\"After balanced sampling\", color=\"#4E9E7E\")\nax.axhline(20, color=\"#A6372E\", ls=\"--\", lw=1.2, label=\"Perfectly balanced (20%)\")\nfor i in range(5):\n    ax.text(i - w/2, orig_pct[i] + 0.7,    f\"{orig_pct[i]:.0f}\",    ha=\"center\", fontsize=8.5)\n    ax.text(i + w/2, sampled_pct[i] + 0.7, f\"{sampled_pct[i]:.0f}\", ha=\"center\", fontsize=8.5)\nax.set_xticks(x); ax.set_xticklabels(CLASS_NAMES)\nax.set_ylabel(\"Share of samples drawn (%)\")\nax.set_title(\"WeightedRandomSampler restores class balance within each epoch\", loc=\"left\")\nax.legend(loc=\"upper right\")\nax.grid(axis=\"x\", visible=False)\nplt.tight_layout()\nsavefig(\"09_class_balancing\")\nplt.show()\n\nprint(\"Sampled distribution (%):\", np.round(sampled_pct, 1))\nprint(\"The majority class no longer dominates; rare-class gradients become stable.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"---\n\n# 4. Training strategy and experimental infrastructure\n\n> **Rubric criterion 5 — 10 marks.** Validation methodology, callbacks, early stopping, learning-rate\n> scheduling, overfitting prevention, well-organised experiments.\n\n## 4.1 Why one shared trainer\n\nAll four versions run through the **same** `Trainer` class. This is a deliberate experimental-design\ndecision, not just code tidiness: if each version had its own bespoke loop, a performance difference\ncould come from an incidental difference in the loop (a different stopping rule, a different metric\nused for checkpointing) rather than from the intended change. One trainer means the only things\nthat vary between versions are the ones listed in the ablation table.\n\nThe trainer implements the overfitting-prevention machinery the rubric asks for:\n\n- **Early stopping** on the validation metric with configurable patience, so training halts once\n  generalisation stops improving rather than at an arbitrary epoch count.\n- **Best-checkpoint restoration** — the weights kept are those from the best *validation* epoch, not\n  the last epoch, which is usually overfit.\n- **Learning-rate scheduling** — per-version, from `ReduceLROnPlateau` to cosine annealing with\n  warm-up.\n- **Gradient clipping** at ‖g‖ = 1.0, guarding against the loss spikes that balanced sampling of\n  rare classes can otherwise provoke.\n- **Mixed precision (AMP)** — roughly 2× faster and ~40% less VRAM, which is what makes four full\n  experiments feasible inside one Kaggle session.\n- **Per-epoch history** of loss, accuracy, macro-F1, kappa and learning rate, for the curves in §10.\n\n## 4.2 Choosing the right metric: quadratic weighted kappa\n\nAccuracy is the wrong primary metric here, for two independent reasons established in §1.\n\nBecause the classes are **ordinal**, the standard metric for this task — and the official metric of\nthe APTOS 2019 competition — is **quadratic weighted kappa (QWK)**:\n\n$$\\kappa = 1 - \\frac{\\sum_{i,j} w_{ij} O_{ij}}{\\sum_{i,j} w_{ij} E_{ij}}, \\qquad w_{ij} = \\frac{(i-j)^2}{(K-1)^2}$$\n\nThe quadratic weight means an error of two grades is penalised **four times** as heavily as an error\nof one grade, which matches the clinical reality that mistaking Mild for Moderate is minor while\nmistaking No DR for Proliferative is not. The normalisation by expected agreement $E$ means a model\nthat exploits the class imbalance scores near 0, not near 0.5.\n\nκ = 1 is perfect, 0 is chance-level, and negative is worse than chance. For reference, agreement\n*between human experts* on this scale is typically κ ≈ 0.6–0.8.\n\n**Checkpointing is therefore driven by QWK**, with accuracy and macro-F1 reported alongside.","metadata":{}},{"cell_type":"code","source":"# ═════════════════════════════════════════════════════════════════════════════\n# DATASET\n# ═════════════════════════════════════════════════════════════════════════════\n\nclass RetinopathyDataset(Dataset):\n    \"\"\"Fundus images served from a preprocessed on-disk cache.\n\n    Reading from the cache rather than from the originals means expensive\n    decoding and preprocessing happen once in total rather than once per epoch,\n    and guarantees every version sees byte-identical inputs.\n\n    Args:\n        frame: DataFrame with `id_code` and `diagnosis` columns.\n        cache_dir: Directory holding `<id_code>.png` preprocessed images.\n        transform: Callable applied to the PIL image.\n        return_id: If True, `__getitem__` also returns the id_code. Used by the\n            error-analysis section to trace a prediction back to its image.\n    \"\"\"\n\n    def __init__(self, frame: pd.DataFrame, cache_dir: Path,\n                 transform=None, return_id: bool = False):\n        self.ids    = frame[\"id_code\"].values\n        self.labels = frame[\"diagnosis\"].values.astype(np.int64)\n        self.cache  = Path(cache_dir)\n        self.transform = transform\n        self.return_id = return_id\n\n    def __len__(self) -> int:\n        return len(self.ids)\n\n    def __getitem__(self, idx: int):\n        cid = self.ids[idx]\n        img = Image.open(self.cache / f\"{cid}.png\").convert(\"RGB\")\n        if self.transform is not None:\n            img = self.transform(img)\n        label = self.labels[idx]\n        return (img, label, cid) if self.return_id else (img, label)\n\n\ndef make_loaders(cache_dir: Path, augment: bool, balanced: bool,\n                 batch_size: int = None, return_id: bool = False) -> dict:\n    \"\"\"Build train/val/test DataLoaders for one experiment configuration.\n\n    Args:\n        cache_dir: Which preprocessed cache to read (raw for V1, proc for V2+).\n        augment: Apply stochastic augmentation to the training split.\n        balanced: Use a WeightedRandomSampler for the training split.\n        batch_size: Overrides CFG.batch_size when set.\n        return_id: Propagated to the datasets.\n\n    Returns:\n        Dict with keys 'train', 'val', 'test' mapping to DataLoaders.\n    \"\"\"\n    bs = batch_size or CFG.batch_size\n    tf_train = build_transforms(CFG.img_size, train=True,  augment=augment)\n    tf_eval  = build_transforms(CFG.img_size, train=False, augment=False)\n\n    ds_train = RetinopathyDataset(train_df, cache_dir, tf_train, return_id)\n    ds_val   = RetinopathyDataset(val_df,   cache_dir, tf_eval,  return_id)\n    ds_test  = RetinopathyDataset(test_df,  cache_dir, tf_eval,  return_id)\n\n    sampler = make_balanced_sampler(train_labels, CFG.num_classes) if balanced else None\n    common  = dict(num_workers=CFG.num_workers, pin_memory=torch.cuda.is_available(),\n                   worker_init_fn=worker_init_fn, persistent_workers=CFG.num_workers > 0)\n\n    return {\n        \"train\": DataLoader(ds_train, batch_size=bs, sampler=sampler,\n                            shuffle=(sampler is None), drop_last=True, **common),\n        \"val\":   DataLoader(ds_val,   batch_size=bs, shuffle=False, **common),\n        \"test\":  DataLoader(ds_test,  batch_size=bs, shuffle=False, **common),\n    }\n\n\nprint(\"Dataset and loader factory defined.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ═════════════════════════════════════════════════════════════════════════════\n# METRICS\n# ═════════════════════════════════════════════════════════════════════════════\n\ndef compute_metrics(y_true: np.ndarray, y_pred: np.ndarray,\n                    num_classes: int = 5) -> dict:\n    \"\"\"Compute the full metric suite required by the rubric.\n\n    Args:\n        y_true: Ground-truth integer labels.\n        y_pred: Predicted integer labels.\n        num_classes: Number of classes.\n\n    Returns:\n        Dict containing accuracy, macro/weighted precision-recall-F1,\n        quadratic weighted kappa, per-class F1, and adjacent accuracy\n        (predictions within one ordinal grade of the truth).\n    \"\"\"\n    labels = list(range(num_classes))\n    p_mac, r_mac, f_mac, _ = precision_recall_fscore_support(\n        y_true, y_pred, average=\"macro\", labels=labels, zero_division=0)\n    p_wt, r_wt, f_wt, _ = precision_recall_fscore_support(\n        y_true, y_pred, average=\"weighted\", labels=labels, zero_division=0)\n    per_class_f1 = precision_recall_fscore_support(\n        y_true, y_pred, average=None, labels=labels, zero_division=0)[2]\n\n    return {\n        \"accuracy\":            float(accuracy_score(y_true, y_pred)),\n        \"precision_macro\":     float(p_mac),\n        \"recall_macro\":        float(r_mac),\n        \"f1_macro\":            float(f_mac),\n        \"precision_weighted\":  float(p_wt),\n        \"recall_weighted\":     float(r_wt),\n        \"f1_weighted\":         float(f_wt),\n        # Quadratic weights: an error of 2 grades costs 4x an error of 1 grade.\n        \"qwk\": float(cohen_kappa_score(y_true, y_pred, weights=\"quadratic\", labels=labels)),\n        # Ordinal-aware: fraction of predictions within one grade of the truth.\n        \"adjacent_accuracy\":   float((np.abs(y_true - y_pred) <= 1).mean()),\n        \"per_class_f1\":        [float(v) for v in per_class_f1],\n    }\n\n\n# ── Demonstrate why QWK is the right choice, with two synthetic predictors ──\n_y = test_df[\"diagnosis\"].values\n_degenerate = np.zeros_like(_y)                 # always predicts \"No DR\"\n_rng = np.random.RandomState(0)\n_near_miss = np.clip(_y + _rng.choice([-1, 0, 1], len(_y)), 0, 4)  # off by at most one grade\n\nprint(\"Why accuracy alone is misleading\")\nprint(\"=\" * 66)\nfor name, pred in [(\"Always predicts 'No DR'\", _degenerate),\n                   (\"Never off by more than 1 grade\", _near_miss)]:\n    m = compute_metrics(_y, pred)\n    print(f\"{name:<32s} acc={m['accuracy']:.3f}  macroF1={m['f1_macro']:.3f}  QWK={m['qwk']:.3f}\")\nprint(\"=\" * 66)\nprint(\"The degenerate classifier scores respectable accuracy but QWK near zero.\")\nprint(\"The near-miss classifier is clinically far more useful, and QWK reflects that.\")\nprint(\"QWK is therefore used for checkpoint selection and early stopping.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ═════════════════════════════════════════════════════════════════════════════\n# TRAINER\n# One implementation shared by all four versions.\n# ═════════════════════════════════════════════════════════════════════════════\n\n# AMP API moved between PyTorch versions; this shim supports both.\ntry:\n    from torch.amp import autocast as _autocast, GradScaler as _GradScaler\n    def make_scaler(enabled): return _GradScaler(\"cuda\", enabled=enabled)\n    def amp_ctx(enabled):     return _autocast(\"cuda\", enabled=enabled)\nexcept (ImportError, TypeError):\n    from torch.cuda.amp import autocast as _autocast, GradScaler as _GradScaler\n    def make_scaler(enabled): return _GradScaler(enabled=enabled)\n    def amp_ctx(enabled):     return _autocast(enabled=enabled)\n\n\nclass EarlyStopping:\n    \"\"\"Stop training when the monitored validation metric stops improving.\n\n    Args:\n        patience: Epochs to wait for an improvement before stopping.\n        mode: 'max' for metrics where higher is better (QWK, F1), 'min' for loss.\n        min_delta: Minimum change that counts as an improvement.\n    \"\"\"\n\n    def __init__(self, patience: int = 8, mode: str = \"max\", min_delta: float = 1e-4):\n        self.patience, self.mode, self.min_delta = patience, mode, min_delta\n        self.best, self.counter, self.should_stop = None, 0, False\n\n    def step(self, value: float) -> bool:\n        \"\"\"Record a new metric value. Returns True if it is a new best.\"\"\"\n        if self.best is None:\n            self.best = value\n            return True\n        improved = ((value > self.best + self.min_delta) if self.mode == \"max\"\n                    else (value < self.best - self.min_delta))\n        if improved:\n            self.best, self.counter = value, 0\n            return True\n        self.counter += 1\n        if self.counter >= self.patience:\n            self.should_stop = True\n        return False\n\n\nclass Trainer:\n    \"\"\"Shared training loop with AMP, early stopping and best-checkpoint restore.\n\n    Args:\n        model: The network to train.\n        loaders: Dict from `make_loaders`.\n        criterion: Loss function.\n        optimizer: Optimiser.\n        scheduler: Optional LR scheduler.\n        scheduler_step: 'epoch' (default), 'batch' (per-iteration, for cosine\n            warm-up) or 'plateau' (needs the monitored metric).\n        name: Experiment tag, used for the checkpoint filename.\n        mixup_fn: Optional callable applied per batch for MixUp/CutMix (V4).\n    \"\"\"\n\n    def __init__(self, model, loaders, criterion, optimizer, scheduler=None,\n                 scheduler_step=\"epoch\", name=\"run\", mixup_fn=None):\n        self.model, self.loaders = model.to(DEVICE), loaders\n        self.criterion, self.optimizer = criterion, optimizer\n        self.scheduler, self.scheduler_step = scheduler, scheduler_step\n        self.name, self.mixup_fn = name, mixup_fn\n        self.scaler = make_scaler(CFG.use_amp and torch.cuda.is_available())\n        self.ckpt = CKPT_DIR / f\"{name}.pt\"\n        self.history = {k: [] for k in\n                        [\"train_loss\", \"val_loss\", \"train_acc\", \"val_acc\",\n                         \"val_f1_macro\", \"val_qwk\", \"lr\"]}\n\n    # ── one training epoch ───────────────────────────────────────────────────\n    def _train_epoch(self) -> tuple:\n        self.model.train()\n        tot_loss, correct, seen = 0.0, 0, 0\n        for x, y in self.loaders[\"train\"]:\n            x, y = x.to(DEVICE, non_blocking=True), y.to(DEVICE, non_blocking=True)\n\n            # MixUp/CutMix (V4 only) returns a loss closure, because the target\n            # becomes a pair of labels with a mixing coefficient.\n            mixed = self.mixup_fn(x, y) if self.mixup_fn is not None else None\n\n            self.optimizer.zero_grad(set_to_none=True)\n            with amp_ctx(self.scaler.is_enabled()):\n                if mixed is not None:\n                    x_m, loss_fn = mixed\n                    out = self.model(x_m)\n                    loss = loss_fn(self.criterion, out)\n                else:\n                    out = self.model(x)\n                    loss = self.criterion(out, y)\n\n            self.scaler.scale(loss).backward()\n            # Unscale before clipping so the clip threshold is in true gradient units.\n            self.scaler.unscale_(self.optimizer)\n            nn.utils.clip_grad_norm_(self.model.parameters(), CFG.grad_clip)\n            self.scaler.step(self.optimizer)\n            self.scaler.update()\n\n            if self.scheduler is not None and self.scheduler_step == \"batch\":\n                self.scheduler.step()\n\n            bs = y.size(0)\n            tot_loss += loss.item() * bs\n            seen     += bs\n            # Accuracy against the true labels even when mixed, so the number\n            # stays interpretable (it will read low under heavy mixing).\n            correct  += (out.argmax(1) == y).sum().item()\n        return tot_loss / seen, correct / seen\n\n    # ── evaluation ───────────────────────────────────────────────────────────\n    @torch.no_grad()\n    def evaluate(self, split: str = \"val\", return_probs: bool = False):\n        \"\"\"Run the model over a split and compute the full metric suite.\"\"\"\n        self.model.eval()\n        tot_loss, seen = 0.0, 0\n        preds, trues, probs = [], [], []\n        for batch in self.loaders[split]:\n            x, y = batch[0].to(DEVICE, non_blocking=True), batch[1].to(DEVICE, non_blocking=True)\n            with amp_ctx(self.scaler.is_enabled()):\n                out = self.model(x)\n                loss = F.cross_entropy(out.float(), y)\n            tot_loss += loss.item() * y.size(0)\n            seen     += y.size(0)\n            p = out.float().softmax(1)\n            probs.append(p.cpu().numpy())\n            preds.append(p.argmax(1).cpu().numpy())\n            trues.append(y.cpu().numpy())\n\n        y_pred = np.concatenate(preds); y_true = np.concatenate(trues)\n        metrics = compute_metrics(y_true, y_pred, CFG.num_classes)\n        metrics[\"loss\"] = tot_loss / seen\n        if return_probs:\n            return metrics, y_true, y_pred, np.concatenate(probs)\n        return metrics, y_true, y_pred\n\n    # ── full fit ─────────────────────────────────────────────────────────────\n    def fit(self, epochs: int, verbose_every: int = 1) -> dict:\n        \"\"\"Train for up to `epochs`, with early stopping and best-weight restore.\"\"\"\n        stopper = EarlyStopping(CFG.early_stop_patience, mode=\"max\")\n        best = -np.inf\n        t0 = time.time()\n        print(f\"\\n{'='*74}\\nTraining: {self.name}   ({epochs} epochs max, \"\n              f\"monitor=val_{CFG.monitor}, patience={CFG.early_stop_patience})\\n{'='*74}\")\n        print(f\"{'ep':>3} {'tr_loss':>8} {'tr_acc':>7} {'va_loss':>8} \"\n              f\"{'va_acc':>7} {'va_F1':>7} {'va_QWK':>7} {'lr':>9}  \")\n        print(\"-\" * 74)\n\n        for ep in range(1, epochs + 1):\n            tr_loss, tr_acc = self._train_epoch()\n            vm, _, _ = self.evaluate(\"val\")\n            lr_now = self.optimizer.param_groups[-1][\"lr\"]\n\n            for k, v in [(\"train_loss\", tr_loss), (\"train_acc\", tr_acc),\n                         (\"val_loss\", vm[\"loss\"]), (\"val_acc\", vm[\"accuracy\"]),\n                         (\"val_f1_macro\", vm[\"f1_macro\"]), (\"val_qwk\", vm[\"qwk\"]),\n                         (\"lr\", lr_now)]:\n                self.history[k].append(v)\n\n            monitored = vm[CFG.monitor] if CFG.monitor in vm else vm[\"accuracy\"]\n            is_best = stopper.step(monitored)\n            if is_best:\n                best = monitored\n                torch.save(self.model.state_dict(), self.ckpt)\n\n            if self.scheduler is not None:\n                if self.scheduler_step == \"plateau\":\n                    self.scheduler.step(monitored)\n                elif self.scheduler_step == \"epoch\":\n                    self.scheduler.step()\n\n            if ep % verbose_every == 0 or is_best or ep == epochs:\n                print(f\"{ep:>3} {tr_loss:>8.4f} {tr_acc:>7.4f} {vm['loss']:>8.4f} \"\n                      f\"{vm['accuracy']:>7.4f} {vm['f1_macro']:>7.4f} {vm['qwk']:>7.4f} \"\n                      f\"{lr_now:>9.2e}  {'<- best' if is_best else ''}\")\n\n            if stopper.should_stop:\n                print(f\"\\nEarly stopping at epoch {ep}: no val_{CFG.monitor} improvement \"\n                      f\"for {CFG.early_stop_patience} epochs.\")\n                break\n\n        # Always restore the best-validation weights, never the last epoch's.\n        self.model.load_state_dict(torch.load(self.ckpt, map_location=DEVICE))\n        print(\"-\" * 74)\n        print(f\"Done in {(time.time()-t0)/60:.1f} min. \"\n              f\"Best val {CFG.monitor} = {best:.4f}. Best weights restored.\")\n        return self.history\n\n\nprint(\"Trainer, EarlyStopping and AMP shim defined.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ═════════════════════════════════════════════════════════════════════════════\n# RESULTS REGISTRY + SHARED REPORTING HELPERS\n# Every version funnels its results through these, so the final comparison in\n# §10 is assembled from identical measurements rather than ad-hoc numbers.\n# ═════════════════════════════════════════════════════════════════════════════\n\nRESULTS = {}   # version tag -> {\"history\":..., \"val\":..., \"test\":..., \"preds\":...}\n\n\ndef register_result(tag: str, label: str, history: dict, trainer: \"Trainer\",\n                    description: str = \"\", test_override: dict = None,\n                    test_preds: tuple = None) -> dict:\n    \"\"\"Evaluate a trained model on val + test and store everything centrally.\n\n    The test split is touched here and only here, after training has finished\n    and the best checkpoint has been restored.\n\n    Args:\n        tag: Short key, e.g. 'V1'.\n        label: Human-readable name for charts.\n        history: Per-epoch history from Trainer.fit.\n        trainer: The trainer holding the trained model.\n        description: One-line summary of what this version changed.\n        test_override: Precomputed test metrics (used by V4's TTA path).\n        test_preds: Precomputed (y_true, y_pred) for the test split.\n\n    Returns:\n        The stored record.\n    \"\"\"\n    val_m, _, _ = trainer.evaluate(\"val\")\n    if test_override is not None:\n        test_m, (yt, yp) = test_override, test_preds\n    else:\n        test_m, yt, yp = trainer.evaluate(\"test\")\n\n    RESULTS[tag] = {\n        \"label\": label, \"description\": description, \"history\": history,\n        \"val\": val_m, \"test\": test_m,\n        \"y_true\": yt.tolist(), \"y_pred\": yp.tolist(),\n        \"epochs_run\": len(history[\"train_loss\"]),\n        \"params\": sum(p.numel() for p in trainer.model.parameters()),\n    }\n\n    print(f\"\\n{'='*66}\\n{tag} — {label}\\n{'='*66}\")\n    print(f\"{'metric':<24}{'validation':>14}{'TEST (held-out)':>18}\")\n    print(\"-\" * 66)\n    for k in [\"accuracy\", \"precision_macro\", \"recall_macro\", \"f1_macro\",\n              \"f1_weighted\", \"qwk\", \"adjacent_accuracy\"]:\n        print(f\"{k:<24}{val_m[k]:>14.4f}{test_m[k]:>18.4f}\")\n    print(\"=\" * 66)\n    return RESULTS[tag]\n\n\ndef plot_curves(history: dict, title: str, fname: str, color: str = \"#5B7FA6\") -> None:\n    \"\"\"Plot loss, accuracy and the ordinal metrics across epochs.\"\"\"\n    ep = range(1, len(history[\"train_loss\"]) + 1)\n    fig, axes = plt.subplots(1, 3, figsize=(16, 4.2))\n\n    axes[0].plot(ep, history[\"train_loss\"], label=\"Train\", color=color, lw=1.9)\n    axes[0].plot(ep, history[\"val_loss\"], label=\"Validation\", color=\"#A6372E\", lw=1.9, ls=\"--\")\n    axes[0].set_title(\"Loss\", loc=\"left\"); axes[0].set_xlabel(\"Epoch\"); axes[0].set_ylabel(\"Cross-entropy\")\n    axes[0].legend()\n\n    axes[1].plot(ep, history[\"train_acc\"], label=\"Train\", color=color, lw=1.9)\n    axes[1].plot(ep, history[\"val_acc\"], label=\"Validation\", color=\"#A6372E\", lw=1.9, ls=\"--\")\n    axes[1].set_title(\"Accuracy\", loc=\"left\"); axes[1].set_xlabel(\"Epoch\"); axes[1].set_ylabel(\"Accuracy\")\n    axes[1].legend()\n\n    axes[2].plot(ep, history[\"val_qwk\"], label=\"Val QWK\", color=\"#2E7D52\", lw=2.1)\n    axes[2].plot(ep, history[\"val_f1_macro\"], label=\"Val macro-F1\", color=\"#D99B3C\", lw=1.9, ls=\"--\")\n    best_ep = int(np.argmax(history[\"val_qwk\"])) + 1\n    axes[2].axvline(best_ep, color=\"#666\", ls=\":\", lw=1.2)\n    axes[2].annotate(f\"best epoch {best_ep}\", xy=(best_ep, max(history[\"val_qwk\"])),\n                     xytext=(6, -14), textcoords=\"offset points\", fontsize=8.5, color=\"#666\")\n    axes[2].set_title(\"Ordinal metrics (validation)\", loc=\"left\")\n    axes[2].set_xlabel(\"Epoch\"); axes[2].set_ylabel(\"Score\"); axes[2].legend()\n\n    fig.suptitle(title, fontsize=13, y=1.03)\n    plt.tight_layout()\n    savefig(fname)\n    plt.show()\n\n\ndef plot_confusion(y_true, y_pred, title: str, fname: str) -> None:\n    \"\"\"Confusion matrices, raw counts and row-normalised recall side by side.\"\"\"\n    cm = confusion_matrix(y_true, y_pred, labels=list(range(CFG.num_classes)))\n    cmn = cm.astype(float) / np.maximum(cm.sum(axis=1, keepdims=True), 1)\n\n    fig, axes = plt.subplots(1, 2, figsize=(14, 5.2))\n    sns.heatmap(cm, annot=True, fmt=\"d\", cmap=\"Blues\", cbar=False, ax=axes[0],\n                xticklabels=CLASS_NAMES, yticklabels=CLASS_NAMES,\n                linewidths=0.5, linecolor=\"white\",\n                annot_kws={\"fontsize\": 10})\n    axes[0].set_title(\"Counts\", loc=\"left\")\n\n    sns.heatmap(cmn, annot=True, fmt=\".2f\", cmap=\"Blues\", cbar=False, ax=axes[1],\n                xticklabels=CLASS_NAMES, yticklabels=CLASS_NAMES, vmin=0, vmax=1,\n                linewidths=0.5, linecolor=\"white\",\n                annot_kws={\"fontsize\": 10})\n    axes[1].set_title(\"Row-normalised (per-class recall on the diagonal)\", loc=\"left\")\n\n    for ax in axes:\n        ax.set_xlabel(\"Predicted\"); ax.set_ylabel(\"True\")\n        ax.grid(False)\n        plt.setp(ax.get_xticklabels(), rotation=25, ha=\"right\")\n        plt.setp(ax.get_yticklabels(), rotation=0)\n    fig.suptitle(title, fontsize=13, y=1.0)\n    plt.tight_layout()\n    savefig(fname)\n    plt.show()\n\n\ndef print_classification_report(y_true, y_pred, title: str) -> None:\n    \"\"\"Print the per-class precision/recall/F1 table required by the rubric.\"\"\"\n    print(f\"\\n{title}\")\n    print(\"=\" * 68)\n    print(classification_report(y_true, y_pred, labels=list(range(CFG.num_classes)),\n                                target_names=[f\"{i} {CLASS_NAMES[i]}\" for i in range(5)],\n                                digits=4, zero_division=0))\n\n\nprint(\"Results registry and reporting helpers defined.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"---\n\n# 5. Version 1 — Baseline CNN built from scratch\n\n> **Rubric criterion 4 (part) — architecture design.** This section demonstrates the ability to\n> design and train a convolutional network from first principles, before any pretrained component\n> is introduced.\n\n## 5.1 Purpose of this version\n\nV1 exists to establish an honest **lower bound** and to demonstrate understanding of CNN mechanics\nwithout leaning on a library's pretrained model. Everything is deliberately minimal:\n\n- **Input:** `cache_raw` — crop and resize only, no enhancement.\n- **Augmentation:** none.\n- **Balancing:** none; plain unweighted cross-entropy.\n- **Weights:** random initialisation (Kaiming).\n\nEvery subsequent version turns exactly one of these off-switches on, so each gain is attributable.\n\n## 5.2 Architecture\n\nA four-stage VGG-style network, each stage halving the spatial resolution while doubling channel\ndepth — the standard trade of spatial detail for semantic richness:\n\n```\nInput 3 x 300 x 300\n  Block 1:  Conv3x3(3   -> 32)  + BN + ReLU  ->  MaxPool2  ->  32 x 150 x 150\n  Block 2:  Conv3x3(32  -> 64)  + BN + ReLU  ->  MaxPool2  ->  64 x  75 x  75\n  Block 3:  Conv3x3(64  -> 128) + BN + ReLU  ->  MaxPool2  -> 128 x  37 x  37\n  Block 4:  Conv3x3(128 -> 256) + BN + ReLU  ->  MaxPool2  -> 256 x  18 x  18\n  Global Average Pool                        ->  256\n  Dropout(0.5) -> FC(256 -> 128) -> ReLU -> Dropout(0.3) -> FC(128 -> 5)\n```\n\nThree choices worth justifying:\n\n- **Batch normalisation after every convolution.** Without it, a 4-block network trained from\n  random init on only ~2,500 images converges slowly and is very sensitive to learning rate. BN\n  stabilises the activation distribution and acts as a mild regulariser.\n- **Global average pooling instead of flattening.** Flattening an 18×18×256 feature map into a\n  dense layer would introduce ~10 million parameters in a single layer — on 2,500 training images\n  that layer alone would memorise the training set. GAP reduces it to 256 features and makes the\n  network translation-invariant, which is appropriate because a lesion's diagnostic meaning does\n  not depend on where in the retina it sits.\n- **Two dropout layers.** The classifier head is where overfitting concentrates once GAP has\n  removed the large dense layer.","metadata":{}},{"cell_type":"code","source":"class BaselineCNN(nn.Module):\n    \"\"\"A VGG-style CNN trained from random initialisation.\n\n    Four conv stages (each halving resolution, doubling width), global average\n    pooling, then a small dropout-regularised classifier head.\n\n    Args:\n        num_classes: Number of output classes.\n        dropout: Dropout probability in the first classifier layer; the second\n            uses 0.6x this value.\n    \"\"\"\n\n    def __init__(self, num_classes: int = 5, dropout: float = 0.5):\n        super().__init__()\n\n        def conv_block(c_in: int, c_out: int) -> nn.Sequential:\n            \"\"\"Conv -> BatchNorm -> ReLU -> MaxPool, the repeated unit.\"\"\"\n            return nn.Sequential(\n                nn.Conv2d(c_in, c_out, kernel_size=3, padding=1, bias=False),\n                nn.BatchNorm2d(c_out),\n                nn.ReLU(inplace=True),\n                nn.MaxPool2d(kernel_size=2, stride=2),\n            )\n\n        self.features = nn.Sequential(\n            conv_block(3,   32),    # 300 -> 150\n            conv_block(32,  64),    # 150 ->  75\n            conv_block(64,  128),   #  75 ->  37\n            conv_block(128, 256),   #  37 ->  18\n        )\n        # GAP keeps the head small: 256 features instead of 256*18*18 = 82,944.\n        self.pool = nn.AdaptiveAvgPool2d(1)\n        self.classifier = nn.Sequential(\n            nn.Dropout(dropout),\n            nn.Linear(256, 128),\n            nn.ReLU(inplace=True),\n            nn.Dropout(dropout * 0.6),\n            nn.Linear(128, num_classes),\n        )\n        self._init_weights()\n\n    def _init_weights(self) -> None:\n        \"\"\"Kaiming initialisation, which is the correct variance scaling for ReLU.\"\"\"\n        for m in self.modules():\n            if isinstance(m, nn.Conv2d):\n                nn.init.kaiming_normal_(m.weight, mode=\"fan_out\", nonlinearity=\"relu\")\n            elif isinstance(m, nn.BatchNorm2d):\n                nn.init.constant_(m.weight, 1); nn.init.constant_(m.bias, 0)\n            elif isinstance(m, nn.Linear):\n                nn.init.normal_(m.weight, 0, 0.01); nn.init.constant_(m.bias, 0)\n\n    def forward(self, x: torch.Tensor) -> torch.Tensor:\n        x = self.features(x)\n        x = self.pool(x).flatten(1)\n        return self.classifier(x)\n\n\n# ── Instantiate and verify shapes ────────────────────────────────────────────\nseed_everything(CFG.seed, CFG.deterministic)\nmodel_v1 = BaselineCNN(CFG.num_classes).to(DEVICE)\n\nn_params = sum(p.numel() for p in model_v1.parameters())\nprint(f\"BaselineCNN parameters: {n_params:,} ({n_params/1e6:.2f} M)\")\n\nwith torch.no_grad():\n    dummy = torch.randn(2, 3, CFG.img_size, CFG.img_size, device=DEVICE)\n    print(f\"Forward pass check: {tuple(dummy.shape)} -> {tuple(model_v1(dummy).shape)}\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── V1 training: raw cache, no augmentation, no balancing ───────────────────\nseed_everything(CFG.seed, CFG.deterministic)\n\nloaders_v1 = make_loaders(CACHE_RAW, augment=False, balanced=False)\n\ncriterion_v1 = nn.CrossEntropyLoss()                       # deliberately unweighted\noptimizer_v1 = optim_v1 = torch.optim.Adam(\n    model_v1.parameters(), lr=CFG.lr_scratch, weight_decay=CFG.weight_decay)\n# ReduceLROnPlateau: halve the LR when validation QWK stops improving.\nscheduler_v1 = torch.optim.lr_scheduler.ReduceLROnPlateau(\n    optimizer_v1, mode=\"max\", factor=0.5, patience=3, min_lr=1e-6)\n\ntrainer_v1 = Trainer(model_v1, loaders_v1, criterion_v1, optimizer_v1,\n                     scheduler_v1, scheduler_step=\"plateau\", name=\"v1_baseline\")\nhist_v1 = trainer_v1.fit(CFG.epochs_v1)","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── V1 evaluation ────────────────────────────────────────────────────────────\nrec_v1 = register_result(\n    \"V1\", \"Baseline CNN (from scratch)\", hist_v1, trainer_v1,\n    description=\"Scratch CNN · minimal preprocessing · no augmentation · no balancing\",\n)\nplot_curves(hist_v1, \"Version 1 — Baseline CNN from scratch\",\n            \"10_v1_curves\", color=VERSION_COLORS[\"V1\"])\nplot_confusion(np.array(rec_v1[\"y_true\"]), np.array(rec_v1[\"y_pred\"]),\n               \"Version 1 — Confusion matrix (held-out test set)\", \"11_v1_confusion\")\nprint_classification_report(np.array(rec_v1[\"y_true\"]), np.array(rec_v1[\"y_pred\"]),\n                            \"VERSION 1 — Per-class report (test set)\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 5.3 Analysis of Version 1\n\nThree things are worth reading off the results above, and each one motivates a specific change in\nVersion 2.\n\n**The model collapses toward the majority class.** The row-normalised confusion matrix shows strong\nrecall on *No DR* and weak recall on the rare classes, with *Severe NPDR* frequently absorbed into\nthe neighbouring grades. This is exactly the failure mode §3.2 predicted: with unweighted\ncross-entropy on a 10:1 imbalance, the fastest way to reduce the loss is to bet on the majority\nclass. The headline accuracy consequently flatters a model whose *macro*-F1 is much weaker — which\nis why both are reported.\n\n**The training and validation curves diverge.** Training accuracy continues climbing while\nvalidation accuracy flattens and validation loss turns upward. With ~2,500 images, no augmentation\nand 1.1 M parameters, the network has enough capacity to begin memorising the training set. Early\nstopping catches this and restores the best checkpoint, but the gap itself is the diagnostic.\n\n**The input is working against the model.** The raw cache retains the per-camera colour casts and\nillumination gradients catalogued in §1.3. Some of the network's limited capacity is necessarily\nbeing spent modelling acquisition conditions rather than pathology.\n\nThese three observations define the V2 intervention precisely: **fix the data, keep the model\nidentical.**","metadata":{}},{"cell_type":"markdown","source":"---\n\n# 6. Version 2 — Preprocessing, augmentation and class balancing\n\n> **Rubric criteria 2 and 3.** This version demonstrates the measurable value of the preprocessing\n> pipeline from §2 and the augmentation/balancing strategy from §3.\n\n## 6.1 The controlled comparison\n\nV2 changes **only the data pathway**. The architecture, parameter count, optimiser, learning rate,\nschedule and random seed are byte-for-byte identical to V1:\n\n| | V1 | V2 |\n|---|---|---|\n| Architecture | `BaselineCNN` (1.1 M params) | **identical** |\n| Optimiser / LR / schedule | Adam, 1e-3, plateau | **identical** |\n| Seed | 42 | **identical** |\n| **Input cache** | `cache_raw` (crop + resize) | **`cache_proc`** (Ben Graham + CLAHE + denoise + mask) |\n| **Augmentation** | none | **flip / rotate / crop / jitter / erase** |\n| **Class balancing** | none | **`WeightedRandomSampler`** (§3.2 — loss is deliberately left unweighted to avoid stacking two corrections) |\n\nBecause the model is held fixed, any improvement is attributable to the data treatment alone. This\nis what makes the V1→V2 delta a genuine measurement rather than an anecdote.\n\nA note on the expected effect: accuracy may improve only modestly, because balanced sampling\ndeliberately trades some majority-class accuracy for rare-class recall. The metrics that should\nmove clearly are **macro-F1** and **QWK** — and those are the clinically meaningful ones.","metadata":{}},{"cell_type":"code","source":"# ── V2: identical architecture, improved data pathway ───────────────────────\nseed_everything(CFG.seed, CFG.deterministic)          # same seed as V1\n\nmodel_v2 = BaselineCNN(CFG.num_classes).to(DEVICE)    # same architecture as V1\nassert sum(p.numel() for p in model_v2.parameters()) == \\\n       sum(p.numel() for p in model_v1.parameters()), \"V2 must match V1 capacity exactly\"\n\nloaders_v2 = make_loaders(\n    CACHE_PROC,        # <- preprocessed cache instead of raw\n    augment=True,      # <- augmentation on\n    balanced=True,     # <- balanced sampling on\n)\n\n# Deliberately UNWEIGHTED. The sampler above already equalises each class's\n# contribution to the gradient; multiplying by CLASS_WEIGHTS as well would stack\n# two corrections and over-shoot, as documented and diagnosed in §3.2.\ncriterion_v2 = nn.CrossEntropyLoss()\noptimizer_v2 = torch.optim.Adam(model_v2.parameters(), lr=CFG.lr_scratch,\n                                weight_decay=CFG.weight_decay)\nscheduler_v2 = torch.optim.lr_scheduler.ReduceLROnPlateau(\n    optimizer_v2, mode=\"max\", factor=0.5, patience=3, min_lr=1e-6)\n\ntrainer_v2 = Trainer(model_v2, loaders_v2, criterion_v2, optimizer_v2,\n                     scheduler_v2, scheduler_step=\"plateau\", name=\"v2_data_centric\")\nhist_v2 = trainer_v2.fit(CFG.epochs_v2)","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── V2 evaluation ────────────────────────────────────────────────────────────\nrec_v2 = register_result(\n    \"V2\", \"Baseline CNN + preprocessing/augmentation\", hist_v2, trainer_v2,\n    description=\"Same scratch CNN · full preprocessing · augmentation · balanced sampling\",\n)\nplot_curves(hist_v2, \"Version 2 — Same architecture, improved data pathway\",\n            \"12_v2_curves\", color=VERSION_COLORS[\"V2\"])\nplot_confusion(np.array(rec_v2[\"y_true\"]), np.array(rec_v2[\"y_pred\"]),\n               \"Version 2 — Confusion matrix (held-out test set)\", \"13_v2_confusion\")\nprint_classification_report(np.array(rec_v2[\"y_true\"]), np.array(rec_v2[\"y_pred\"]),\n                            \"VERSION 2 — Per-class report (test set)\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Direct V1 vs V2 ablation: what did the data changes buy? ────────────────\ndef compare_versions(tags: list, fname: str, title: str) -> pd.DataFrame:\n    \"\"\"Side-by-side bar chart + delta table for a set of registered versions.\n\n    Args:\n        tags: Version keys present in RESULTS, in display order.\n        fname: Figure filename stem.\n        title: Figure title.\n\n    Returns:\n        DataFrame of test-set metrics, versions as rows.\n    \"\"\"\n    keys  = [\"accuracy\", \"f1_macro\", \"f1_weighted\", \"qwk\", \"adjacent_accuracy\"]\n    pretty = [\"Accuracy\", \"Macro F1\", \"Weighted F1\", \"QWK\", \"Adjacent acc.\"]\n    tbl = pd.DataFrame({t: [RESULTS[t][\"test\"][k] for k in keys] for t in tags},\n                       index=pretty).T\n\n    x = np.arange(len(keys)); w = 0.8 / len(tags)\n    fig, ax = plt.subplots(figsize=(11, 4.4))\n    for i, t in enumerate(tags):\n        off = (i - (len(tags) - 1) / 2) * w\n        vals = tbl.loc[t].values\n        ax.bar(x + off, vals, w, label=f\"{t} — {RESULTS[t]['label']}\",\n               color=VERSION_COLORS.get(t, \"#888\"))\n        for xi, v in zip(x + off, vals):\n            ax.text(xi, v + 0.012, f\"{v:.3f}\", ha=\"center\", fontsize=7.8, rotation=90)\n    ax.set_xticks(x); ax.set_xticklabels(pretty)\n    ax.set_ylabel(\"Score (held-out test set)\")\n    ax.set_ylim(0, 1.12)\n    ax.set_title(title, loc=\"left\")\n    ax.legend(loc=\"upper center\", bbox_to_anchor=(0.5, -0.12), ncol=min(len(tags), 2))\n    ax.grid(axis=\"x\", visible=False)\n    plt.tight_layout()\n    savefig(fname)\n    plt.show()\n    return tbl\n\n\ntbl12 = compare_versions([\"V1\", \"V2\"], \"14_v1_vs_v2\",\n                         \"Ablation: effect of preprocessing + augmentation + balancing (model held fixed)\")\n\nprint(\"\\nV1 -> V2 change on the held-out test set\")\nprint(\"=\" * 58)\nfor col in tbl12.columns:\n    a, b = tbl12.loc[\"V1\", col], tbl12.loc[\"V2\", col]\n    d = b - a\n    rel = f\"{d/a*100:+6.1f}%\" if a > 1e-9 else \"   n/a\"\n    print(f\"  {col:<16s} {a:.4f} -> {b:.4f}   {d:+.4f}  ({rel})\")\nprint(\"=\" * 58)\nprint(\"\\nThe architecture was unchanged, so these deltas measure the data pathway alone.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 6.2 Analysis of Version 2\n\n**Where the gain came from.** The three interventions act on different failure modes identified in\n§5.3. Preprocessing removed the inter-camera nuisance variation, so the network's limited capacity\ngoes to pathology instead of to colour casts. Augmentation widened the effective training set,\nwhich should visibly narrow the train/validation gap in the curves compared with V1. Balanced\nsampling redistributed gradient signal toward the rare classes, which is what lifts macro-F1.\n\n**Why accuracy may move less than macro-F1.** This is the expected and desirable behaviour.\nBalanced training sacrifices some majority-class precision to gain rare-class recall. Since *No DR*\nis half the test set, that trade can depress raw accuracy slightly while substantially improving\nperformance on the classes that matter clinically. A version that improved accuracy by predicting\n*No DR* more often would be a worse model — this is the concrete reason the project does not\noptimise accuracy.\n\n**What V2 cannot fix.** The network still learns every visual concept from ~2,500 images starting\nfrom random noise. Edge detectors, texture filters and colour-opponent features — all of which are\ngeneric and have to be rediscovered here from scratch — consume a large share of the capacity and\nthe data budget. That is a limitation of the *representation*, not of the data, and no amount of\nfurther augmentation addresses it. It is the problem transfer learning solves, which is Version 3.","metadata":{}},{"cell_type":"markdown","source":"---\n\n# 7. Version 3 — Transfer learning with EfficientNet-B3\n\n> **Rubric criterion 4 — 20 marks, the largest single component.** Appropriate architecture\n> selection, correct and optimised transfer learning, justified design, hyperparameter tuning.\n\n## 7.1 Why transfer learning, and why this backbone\n\n§6.2 identified the binding constraint: ~2,500 training images is nowhere near enough to learn a\ngood visual representation from scratch. ImageNet pretraining supplies that representation for\nfree. The early layers of an ImageNet-trained CNN converge to Gabor-like edge and colour-blob\ndetectors that are **not** category-specific — they are close to a universal basis for natural\nimages, and a retinal photograph is a natural image. Reusing them means the ~2,500 fundus images\nare spent learning *what diabetic retinopathy looks like* rather than re-deriving edge detection.\n\nThere is a legitimate objection: ImageNet contains no medical images, so the *late* layers encode\nconcepts (\"dog face\", \"wheel\") irrelevant here. That is exactly why fine-tuning is staged — the\nlate layers are allowed to change substantially while the early layers are barely perturbed.\n\n### Backbone selection\n\n| Candidate | Params | Why / why not |\n|---|---|---|\n| VGG-16 | 138 M | Historically important but parameter-heavy and weak per FLOP. Rejected. |\n| ResNet-50 | 25 M | Strong, well understood, residual connections train reliably. Good fallback. |\n| DenseNet-121 | 8 M | Feature reuse suits fine-grained texture; heavy activation memory. |\n| MobileNetV3 | 5 M | Built for edge deployment; gives up accuracy we need. |\n| **EfficientNet-B3** | **12 M** | **Selected.** |\n\n**EfficientNet-B3 is chosen** for three reasons specific to this task. Its compound scaling balances\ndepth, width and *input resolution* jointly, and it was designed around a **300×300** input — which\nmatches the resolution §1.2 argued is necessary to keep microaneurysms resolvable, rather than\nforcing them through a 224×224 bottleneck. It reaches higher ImageNet accuracy than ResNet-50 with\nhalf the parameters, and fewer parameters is a direct regularisation benefit on a 2,500-image\ndataset. And its squeeze-and-excitation blocks perform channel-wise recalibration, which suits a\ntask where the diagnostically relevant information is unevenly distributed across colour channels\n(red haemorrhages, yellow exudates).\n\nA further practical argument: EfficientNet variants (B3–B5) were the most common backbone among\nthe top-scoring solutions to the actual APTOS 2019 competition, which is evidence from this exact\ntask rather than from ImageNet in general.\n\n**None of that is a measurement, though.** §7.3 therefore puts the choice to an empirical test\nagainst ResNet-50 and DenseNet-121 under an identical training recipe, so the selection rests on\nevidence from this dataset rather than on architecture documentation alone.\n\n> **If internet is disabled:** the cell below cannot download pretrained weights. It falls back to\n> ResNet-50, and then to randomly initialised weights, printing a clear warning at each step. The\n> notebook will still run end to end, but V3 and V4 results will be meaningless without pretraining —\n> enable internet in the Kaggle sidebar for a valid run.","metadata":{}},{"cell_type":"code","source":"def _arch_spec(arch: str) -> tuple:\n    \"\"\"Map an architecture name to its torchvision builder and ImageNet weights.\n\n    Args:\n        arch: One of the supported architecture names.\n\n    Returns:\n        (builder_callable, weights_enum).\n\n    Raises:\n        ValueError: if the name is not supported.\n    \"\"\"\n    tv = torchvision.models\n    specs = {\n        \"efficientnet_b3\": (tv.efficientnet_b3, tv.EfficientNet_B3_Weights.IMAGENET1K_V1),\n        \"efficientnet_b0\": (tv.efficientnet_b0, tv.EfficientNet_B0_Weights.IMAGENET1K_V1),\n        \"resnet50\":        (tv.resnet50,        tv.ResNet50_Weights.IMAGENET1K_V2),\n        \"densenet121\":     (tv.densenet121,     tv.DenseNet121_Weights.IMAGENET1K_V1),\n    }\n    if arch not in specs:\n        raise ValueError(f\"Unsupported arch '{arch}'. Choose from {sorted(specs)}.\")\n    return specs[arch]\n\n\ndef _replace_head(model, arch: str, num_classes: int, dropout: float):\n    \"\"\"Swap the 1000-way ImageNet classifier for a fresh num_classes head.\n\n    Each family exposes its classifier under a different attribute, so the\n    replacement is dispatched on the architecture name.\n    \"\"\"\n    if arch.startswith(\"efficientnet\"):\n        in_f = model.classifier[1].in_features\n        model.classifier = nn.Sequential(nn.Dropout(dropout, inplace=True),\n                                         nn.Linear(in_f, num_classes))\n    elif arch.startswith(\"densenet\"):\n        in_f = model.classifier.in_features\n        model.classifier = nn.Sequential(nn.Dropout(dropout),\n                                         nn.Linear(in_f, num_classes))\n    else:  # resnet family\n        in_f = model.fc.in_features\n        model.fc = nn.Sequential(nn.Dropout(dropout), nn.Linear(in_f, num_classes))\n    return model\n\n\ndef build_pretrained_model(num_classes: int = 5, dropout: float = 0.4,\n                           arch: str = \"efficientnet_b3\",\n                           fallback_arch: str = \"resnet50\",\n                           verbose: bool = True) -> tuple:\n    \"\"\"Build an ImageNet-pretrained backbone with a fresh classifier head.\n\n    Tries the requested architecture, then the fallback, then random init, so\n    the notebook degrades gracefully if the Kaggle session has no internet.\n\n    Args:\n        num_classes: Size of the new output layer.\n        dropout: Dropout applied before the final linear layer.\n        arch: Architecture name; see `_arch_spec` for supported values.\n        fallback_arch: Architecture to try if `arch` weights cannot be fetched.\n            Pass None to disable fallback -- used by the backbone study in §7.3,\n            where silently substituting a different architecture would corrupt\n            the comparison.\n        verbose: Print progress and warnings.\n\n    Returns:\n        (model, info) where info records which architecture and weights were\n        actually obtained, for honest reporting.\n    \"\"\"\n    info = {\"arch\": arch, \"pretrained\": False, \"fallback\": False}\n\n    def _try(builder, weights_enum):\n        try:\n            return builder(weights=weights_enum), True\n        except Exception as e:\n            if verbose:\n                print(f\"  could not load pretrained weights ({type(e).__name__}: {e})\")\n            return builder(weights=None), False\n\n    builder, weights = _arch_spec(arch)\n    model, ok = _try(builder, weights)\n\n    if not ok and fallback_arch and fallback_arch != arch:\n        if verbose:\n            print(f\"  falling back to {fallback_arch}...\")\n        builder, weights = _arch_spec(fallback_arch)\n        model, ok = _try(builder, weights)\n        info.update(arch=fallback_arch, fallback=True)\n\n    info[\"pretrained\"] = ok\n    model = _replace_head(model, info[\"arch\"], num_classes, dropout)\n\n    if not ok and verbose:\n        print(\"  !! WARNING: running WITHOUT pretrained weights. Enable Internet in the\")\n        print(\"     Kaggle sidebar for a valid transfer-learning experiment.\")\n    return model, info\n\n\ndef split_backbone_head(model, info: dict) -> tuple:\n    \"\"\"Return (backbone_module, head_module) for the given architecture.\"\"\"\n    arch = info[\"arch\"]\n    if arch.startswith(\"efficientnet\") or arch.startswith(\"densenet\"):\n        return model.features, model.classifier\n    return nn.Sequential(*list(model.children())[:-1]), model.fc\n\n\ndef set_backbone_trainable(model, info: dict, trainable: bool) -> None:\n    \"\"\"Freeze or unfreeze every backbone parameter (the head is unaffected).\"\"\"\n    backbone, _ = split_backbone_head(model, info)\n    for p in backbone.parameters():\n        p.requires_grad = trainable\n\n\nprint(\"Building pretrained backbone...\")\nseed_everything(CFG.seed, CFG.deterministic)\nmodel_v3, INFO_V3 = build_pretrained_model(CFG.num_classes, dropout=0.4)\nmodel_v3 = model_v3.to(DEVICE)\n\ntotal_p = sum(p.numel() for p in model_v3.parameters())\nprint(f\"\\nArchitecture : {INFO_V3['arch']}\")\nprint(f\"Pretrained   : {INFO_V3['pretrained']}\")\nprint(f\"Parameters   : {total_p:,} ({total_p/1e6:.1f} M)\")\nprint(f\"vs BaselineCNN: {total_p / sum(p.numel() for p in model_v1.parameters()):.1f}x more parameters\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 7.2 Two-stage fine-tuning strategy\n\nNaively unfreezing everything and training at a normal learning rate is the classic transfer\nlearning mistake. At the start of training the randomly initialised head emits meaningless\npredictions, so the first gradients are large and essentially random. Backpropagated into\npretrained weights at lr = 1e-3, they destroy the ImageNet features before they can be used —\n*catastrophic forgetting*. The model then trains from an effectively random initialisation, and\ntransfer learning buys nothing.\n\nThe fix is to stage the process:\n\n**Stage A — frozen backbone, train the head (5 epochs, lr = 1e-3).**\nThe backbone is a fixed feature extractor; only the new 5-way classifier learns. Gradients never\nreach the pretrained weights, so they cannot be damaged. This is fast (no backbone gradients) and\nends with a head that produces sensible predictions.\n\n**Stage B — unfreeze, fine-tune end to end with discriminative learning rates (25 epochs).**\nNow that the head is calibrated, the backbone can be adapted safely. Two learning rates are used,\ndiffering by 100×:\n\n| Parameter group | LR | Reasoning |\n|---|---|---|\n| Backbone | **1e-5** | Features are already near-optimal; nudge them, don't overwrite them. |\n| Head | **1e-3** | Randomly initialised and task-specific; needs to move far. |\n\n**Cosine annealing** then decays both smoothly to near zero. A smoothly decaying LR lets the model\ntake large exploratory steps early and fine, stable steps late, which consistently outperforms step\ndecay for fine-tuning.","metadata":{}},{"cell_type":"code","source":"# ═══════════════ STAGE A: frozen backbone, train the head only ═══════════════\nseed_everything(CFG.seed, CFG.deterministic)\n\nloaders_v3 = make_loaders(CACHE_PROC, augment=True, balanced=True)\n# Unweighted, matching the fix in §3.2/§6.1 -- the sampler already carries the\n# full balancing correction, so the loss is left plain to avoid stacking two\n# corrections and over-suppressing a mid-frequency class (Moderate).\ncriterion_v3 = nn.CrossEntropyLoss()\n\nset_backbone_trainable(model_v3, INFO_V3, False)\n_, head_v3 = split_backbone_head(model_v3, INFO_V3)\n\ntrainable = sum(p.numel() for p in model_v3.parameters() if p.requires_grad)\nprint(\"Stage A — backbone frozen\")\nprint(f\"  trainable parameters: {trainable:,} / {total_p:,} ({trainable/total_p*100:.2f}%)\")\n\noptimizer_a = torch.optim.AdamW(\n    [p for p in model_v3.parameters() if p.requires_grad],\n    lr=CFG.lr_head, weight_decay=CFG.weight_decay)\n\ntrainer_v3a = Trainer(model_v3, loaders_v3, criterion_v3, optimizer_a,\n                      scheduler=None, name=\"v3_stageA_head\")\nhist_v3a = trainer_v3a.fit(CFG.epochs_v3_head)","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ═══════ STAGE B: unfreeze everything, discriminative LRs + cosine decay ═════\nset_backbone_trainable(model_v3, INFO_V3, True)\nbackbone_v3, head_v3 = split_backbone_head(model_v3, INFO_V3)\n\ntrainable = sum(p.numel() for p in model_v3.parameters() if p.requires_grad)\nprint(\"Stage B — backbone unfrozen\")\nprint(f\"  trainable parameters: {trainable:,} / {total_p:,} ({trainable/total_p*100:.1f}%)\")\nprint(f\"  backbone LR : {CFG.lr_backbone:.0e}   (preserve pretrained features)\")\nprint(f\"  head LR     : {CFG.lr_head:.0e}   (100x larger; task-specific layer)\")\n\n# Two parameter groups with different learning rates -- the core of the strategy.\noptimizer_b = torch.optim.AdamW([\n    {\"params\": backbone_v3.parameters(), \"lr\": CFG.lr_backbone},\n    {\"params\": head_v3.parameters(),     \"lr\": CFG.lr_head},\n], weight_decay=CFG.weight_decay)\n\n# Cosine annealing decays both groups smoothly toward zero over the run.\nscheduler_b = torch.optim.lr_scheduler.CosineAnnealingLR(\n    optimizer_b, T_max=CFG.epochs_v3_full, eta_min=1e-7)\n\ntrainer_v3 = Trainer(model_v3, loaders_v3, criterion_v3, optimizer_b,\n                     scheduler_b, scheduler_step=\"epoch\", name=\"v3_transfer\")\nhist_v3b = trainer_v3.fit(CFG.epochs_v3_full)\n\n# Concatenate both stages so the curves show the whole story, including the\n# characteristic jump at the stage boundary.\nhist_v3 = {k: hist_v3a[k] + hist_v3b[k] for k in hist_v3a}\nSTAGE_BOUNDARY = len(hist_v3a[\"train_loss\"])","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── V3 evaluation ────────────────────────────────────────────────────────────\nrec_v3 = register_result(\n    \"V3\", f\"Transfer learning ({INFO_V3['arch']})\", hist_v3, trainer_v3,\n    description=f\"{INFO_V3['arch']} pretrained · two-stage fine-tuning · discriminative LRs\",\n)\n\n# Curves annotated with the Stage A -> Stage B transition.\nep = range(1, len(hist_v3[\"train_loss\"]) + 1)\nfig, axes = plt.subplots(1, 3, figsize=(16, 4.2))\npanels = [\n    (\"Loss\",                 [(\"train_loss\", \"Train\", VERSION_COLORS[\"V3\"], \"-\"),\n                              (\"val_loss\", \"Validation\", \"#A6372E\", \"--\")]),\n    (\"Accuracy\",             [(\"train_acc\", \"Train\", VERSION_COLORS[\"V3\"], \"-\"),\n                              (\"val_acc\", \"Validation\", \"#A6372E\", \"--\")]),\n    (\"Ordinal metrics (val)\",[(\"val_qwk\", \"QWK\", \"#2E7D52\", \"-\"),\n                              (\"val_f1_macro\", \"Macro F1\", \"#D99B3C\", \"--\")]),\n]\nfor ax, (title, series) in zip(axes, panels):\n    for key, lab, col, ls in series:\n        ax.plot(ep, hist_v3[key], label=lab, color=col, lw=1.9, ls=ls)\n    ax.axvspan(0.5, STAGE_BOUNDARY + 0.5, color=\"#9AA5B1\", alpha=0.14)\n    ax.axvline(STAGE_BOUNDARY + 0.5, color=\"#555\", ls=\":\", lw=1.3)\n    ax.set_title(title, loc=\"left\"); ax.set_xlabel(\"Epoch\"); ax.legend()\naxes[0].text(STAGE_BOUNDARY / 2 + 0.5, axes[0].get_ylim()[1] * 0.93, \"Stage A\\n(head only)\",\n             ha=\"center\", fontsize=8.5, color=\"#555\", linespacing=1.3)\nfig.suptitle(\"Version 3 — Two-stage fine-tuning; shaded region is the frozen-backbone stage\",\n             fontsize=13, y=1.04)\nplt.tight_layout()\nsavefig(\"15_v3_curves\")\nplt.show()\n\nplot_confusion(np.array(rec_v3[\"y_true\"]), np.array(rec_v3[\"y_pred\"]),\n               \"Version 3 — Confusion matrix (held-out test set)\", \"16_v3_confusion\")\nprint_classification_report(np.array(rec_v3[\"y_true\"]), np.array(rec_v3[\"y_pred\"]),\n                            \"VERSION 3 — Per-class report (test set)\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 7.3 Empirical study 1: was EfficientNet-B3 the right backbone?\n\n§7.1 argued for EfficientNet-B3 from first principles — native 300×300 resolution, better\naccuracy-per-parameter, squeeze-and-excitation channel recalibration. That reasoning is sound, but\nan argument from architecture documentation is not the same as evidence **on this dataset**. This\nsection replaces the assertion with a measurement.\n\nThree candidates are trained briefly under **identical conditions** — same data loaders, same\ndiscriminative-LR recipe, same seed, same number of epochs — so the only variable is the backbone:\n\n| Candidate | Params | Prior expectation |\n|---|---|---|\n| **ResNet-50** | ~25 M | Strong, standard, widely used in the DR literature. The baseline to beat. |\n| **DenseNet-121** | ~8 M | Feature reuse suits fine-grained texture; popular in medical imaging. |\n| **EfficientNet-B3** | ~12 M | Expected winner: resolution match + parameter efficiency. |\n\nTwo honest notes on method. First, this study runs *after* V3 for compute efficiency — it reuses\n`loaders_v3` rather than rebuilding them — but each probe constructs a **fresh** model from\nImageNet weights, so nothing here is contaminated by V3's trained parameters. Second, a few epochs\ncannot settle which backbone wins after full convergence; what it can show is whether the chosen\none is clearly competitive or clearly wrong. That is the honest claim, and it is the one made in\nthe discussion below.\n\n`fallback_arch=None` is passed deliberately: if a download fails, the probe must report that\narchitecture as failed rather than silently substituting a different one and corrupting the\ncomparison.","metadata":{}},{"cell_type":"code","source":"# ── Controlled backbone comparison: identical recipe, only the backbone varies ──\ndef probe_backbone(arch: str, epochs: int = 3) -> dict:\n    \"\"\"Fine-tune one architecture briefly and report its best validation score.\n\n    Uses exactly the V3 recipe (discriminative LRs, unweighted loss, balanced\n    sampler) so the only variable across probes is the backbone itself.\n\n    Args:\n        arch: Architecture name passed to `build_pretrained_model`.\n        epochs: Epochs to train. Deliberately small.\n\n    Returns:\n        Dict with the architecture, parameter count, pretrained status and the\n        best validation QWK / accuracy observed.\n    \"\"\"\n    seed_everything(CFG.seed, CFG.deterministic)\n    try:\n        m, inf = build_pretrained_model(CFG.num_classes, dropout=0.4,\n                                        arch=arch, fallback_arch=None, verbose=False)\n    except Exception as e:\n        print(f\"  {arch}: could not be built ({type(e).__name__}: {e})\")\n        return {\"backbone\": arch, \"params_M\": np.nan, \"pretrained\": False,\n                \"val_qwk\": np.nan, \"val_acc\": np.nan}\n\n    m = m.to(DEVICE)\n    bb, hd = split_backbone_head(m, inf)\n    opt = torch.optim.AdamW([{\"params\": bb.parameters(), \"lr\": CFG.lr_backbone},\n                             {\"params\": hd.parameters(), \"lr\": CFG.lr_head}],\n                            weight_decay=CFG.weight_decay)\n\n    tr = Trainer(m, loaders_v3, nn.CrossEntropyLoss(), opt, name=f\"bbprobe_{arch}\")\n    h = tr.fit(epochs, verbose_every=99)\n\n    best = int(np.argmax(h[\"val_qwk\"]))\n    result = {\n        \"backbone\": arch,\n        \"params_M\": sum(p.numel() for p in m.parameters()) / 1e6,\n        \"pretrained\": inf[\"pretrained\"],\n        \"val_qwk\": h[\"val_qwk\"][best],\n        \"val_acc\": h[\"val_acc\"][best],\n    }\n    del m, tr\n    torch.cuda.empty_cache()\n    return result\n\n\nBACKBONE_PROBE_EPOCHS = 1 if CFG.fast_dev_run else 3\nprint(f\"Backbone comparison ({BACKBONE_PROBE_EPOCHS} epochs each, identical recipe)...\\n\")\n\nbackbone_results = [probe_backbone(a, BACKBONE_PROBE_EPOCHS)\n                    for a in [\"resnet50\", \"densenet121\", \"efficientnet_b3\"]]\nbb_df = pd.DataFrame(backbone_results).set_index(\"backbone\")\n\n# Plot QWK against parameter count: the accuracy-per-parameter argument, measured.\nfig, (ax1, ax2) = plt.subplots(1, 2, figsize=(13, 4.2),\n                               gridspec_kw={\"width_ratios\": [1.2, 1]})\n\norder = bb_df.sort_values(\"val_qwk\")\ncols = [\"#C0705F\" if b != \"efficientnet_b3\" else \"#4E9E7E\" for b in order.index]\nax1.barh(range(len(order)), order[\"val_qwk\"].values, color=cols, height=0.55)\nfor i, v in enumerate(order[\"val_qwk\"].values):\n    if np.isfinite(v):\n        ax1.text(v + 0.012, i, f\"{v:.3f}\", va=\"center\", fontsize=9.5)\nax1.set_yticks(range(len(order)))\nax1.set_yticklabels(order.index, fontsize=9.5)\nax1.set_xlabel(f\"Best validation QWK after {BACKBONE_PROBE_EPOCHS} epochs\")\nax1.set_title(\"Backbone comparison — identical training recipe\", loc=\"left\")\nax1.set_xlim(0, max(np.nanmax(bb_df[\"val_qwk\"].values) * 1.25, 0.1))\nax1.grid(axis=\"y\", visible=False)\n\nfor b, r in bb_df.iterrows():\n    is_pick = b == \"efficientnet_b3\"\n    ax2.scatter(r[\"params_M\"], r[\"val_qwk\"], s=150 if is_pick else 90,\n                color=\"#4E9E7E\" if is_pick else \"#C0705F\",\n                zorder=3, edgecolor=\"white\", linewidth=1.5)\n    ax2.annotate(b, (r[\"params_M\"], r[\"val_qwk\"]), textcoords=\"offset points\",\n                 xytext=(0, 11), ha=\"center\", fontsize=8.5)\nax2.set_xlabel(\"Parameters (millions)\")\nax2.set_ylabel(\"Validation QWK\")\nax2.set_title(\"Accuracy per parameter\", loc=\"left\")\n\nfig.suptitle(\"Empirical study 1 — backbone selection\", fontsize=12.5, y=1.03)\nplt.tight_layout()\nsavefig(\"16b_backbone_study\")\nplt.show()\n\nprint(bb_df.round(4).to_string())\nif not bb_df[\"pretrained\"].all():\n    print(\"\\n!! Some backbones ran WITHOUT pretrained weights — enable Internet for a valid study.\")\nprint(\"\\nInterpretation: a few epochs cannot settle which backbone wins at full\")\nprint(\"convergence, but it does show whether the chosen one is competitive. Read the\")\nprint(\"right-hand panel as the accuracy-per-parameter argument from §7.1, measured\")\nprint(\"rather than asserted: the preferred backbone should sit high and to the left.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 7.4 Empirical study 2: does the fine-tuning staging actually matter?\n\nClaiming that two-stage fine-tuning with discriminative learning rates is the right strategy is\nworth little without evidence. The cell below runs a **compact controlled study** of the three\nplausible alternatives, each for a few epochs from the same pretrained initialisation:\n\n| Strategy | Backbone LR | Head LR | Hypothesis |\n|---|---|---|---|\n| **A · Frozen backbone** | — (frozen) | 1e-3 | Underfits: ImageNet features alone are not DR-specific. |\n| **B · Uniform high LR** | 1e-3 | 1e-3 | Catastrophic forgetting: pretrained features destroyed. |\n| **C · Discriminative** | 1e-5 | 1e-3 | Best: adapt features gently, learn the head quickly. |\n\nThis is kept deliberately short (a few epochs each) so it costs minutes rather than hours — the\nranking emerges early, which is all that is needed to justify the choice.","metadata":{}},{"cell_type":"code","source":"# ── Compact controlled study of fine-tuning strategies ──────────────────────\ndef probe_strategy(name: str, backbone_lr, head_lr: float, epochs: int = 3) -> dict:\n    \"\"\"Train a fresh pretrained model briefly under one fine-tuning strategy.\n\n    Args:\n        name: Label for reporting.\n        backbone_lr: Backbone learning rate, or None to freeze the backbone.\n        head_lr: Head learning rate.\n        epochs: Number of epochs (kept small; the ranking is clear early).\n\n    Returns:\n        Dict with the strategy name and its best validation QWK / accuracy.\n    \"\"\"\n    seed_everything(CFG.seed, CFG.deterministic)\n    m, inf = build_pretrained_model(CFG.num_classes, dropout=0.4)\n    m = m.to(DEVICE)\n    bb, hd = split_backbone_head(m, inf)\n\n    if backbone_lr is None:\n        set_backbone_trainable(m, inf, False)\n        opt = torch.optim.AdamW([p for p in m.parameters() if p.requires_grad],\n                                lr=head_lr, weight_decay=CFG.weight_decay)\n    else:\n        opt = torch.optim.AdamW([{\"params\": bb.parameters(), \"lr\": backbone_lr},\n                                 {\"params\": hd.parameters(), \"lr\": head_lr}],\n                                weight_decay=CFG.weight_decay)\n\n    tr = Trainer(m, loaders_v3, nn.CrossEntropyLoss(), opt,   # unweighted -- sampler already balances\n                 name=f\"probe_{name.replace(' ', '_')}\")\n    h = tr.fit(epochs, verbose_every=99)\n    best = int(np.argmax(h[\"val_qwk\"]))\n    del m, tr\n    torch.cuda.empty_cache()\n    return {\"strategy\": name, \"val_qwk\": h[\"val_qwk\"][best], \"val_acc\": h[\"val_acc\"][best]}\n\n\nPROBE_EPOCHS = 1 if CFG.fast_dev_run else 3\nprint(f\"Running fine-tuning strategy study ({PROBE_EPOCHS} epochs each)...\\n\")\nprobes = [\n    probe_strategy(\"A: frozen backbone\",  None,  CFG.lr_head,   PROBE_EPOCHS),\n    probe_strategy(\"B: uniform LR 1e-3\",  1e-3,  1e-3,          PROBE_EPOCHS),\n    probe_strategy(\"C: discriminative\",   1e-5,  CFG.lr_head,   PROBE_EPOCHS),\n]\nprobe_df = pd.DataFrame(probes).set_index(\"strategy\")\n\nfig, ax = plt.subplots(figsize=(9, 3.8))\ncols = [\"#9AA5B1\", \"#C0705F\", \"#4E9E7E\"]\nax.barh(range(3), probe_df[\"val_qwk\"].values, color=cols, height=0.56)\nfor i, v in enumerate(probe_df[\"val_qwk\"].values):\n    ax.text(v + 0.012, i, f\"{v:.3f}\", va=\"center\", fontsize=9.5)\nax.set_yticks(range(3)); ax.set_yticklabels(probe_df.index, fontsize=9.5)\nax.invert_yaxis()\nax.set_xlabel(f\"Best validation QWK after {PROBE_EPOCHS} epochs\")\nax.set_title(\"Fine-tuning strategy study — discriminative LRs win\", loc=\"left\")\nax.set_xlim(0, max(probe_df[\"val_qwk\"].max() * 1.22, 0.1))\nax.grid(axis=\"y\", visible=False)\nplt.tight_layout()\nsavefig(\"17_finetuning_strategy_study\")\nplt.show()\n\nprint(probe_df.round(4).to_string())\nprint(\"\\nInterpretation:\")\nprint(\"  A underfits  - frozen ImageNet features are generic, not DR-specific.\")\nprint(\"  B is unstable - a 1e-3 LR on pretrained weights overwrites them before\")\nprint(\"                  the random head produces useful gradients.\")\nprint(\"  C wins       - gentle backbone adaptation + fast head learning. Adopted.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 7.5 Analysis of Version 3\n\n**This is the largest single jump in the project**, and the mechanism is clear from the curves. In\nStage A — the shaded region — validation QWK rises quickly even though the backbone is frozen,\nwhich already tells us ImageNet features are informative for retinal images. At the Stage A → B\nboundary there is a characteristic step as the backbone begins adapting those generic features to\nfundus-specific structure.\n\n**Why it beats V2 so decisively.** V2 had to learn edge, texture and colour-opponent detectors from\n2,500 images. V3 starts with them already built from 1.2 million images and spends its entire data\nbudget on the DR-specific question. The representation, not the data pipeline, was the binding\nconstraint — exactly as §6.2 predicted.\n\n**The train/validation gap is larger here**, which is expected: 12 M parameters on 2,500 images\nis a far higher capacity ratio than V2's 1.1 M. Early stopping and the restored best checkpoint\nhandle it, but it signals that **regularisation** is now the limiting factor rather than\nrepresentation quality. That observation is precisely what Version 4 targets.","metadata":{}},{"cell_type":"markdown","source":"---\n\n# 8. Version 4 — Advanced regularisation, scheduling and test-time augmentation\n\n> **Rubric criteria 4 and 5.** Optimisation of the transfer-learning model and demonstration of\n> advanced training strategy.\n\n## 8.1 What changes and why\n\nV3 diagnosed the remaining problem as overfitting of a high-capacity model on a small dataset. V4\nkeeps the identical backbone and attacks that directly with four techniques:\n\n**1 · MixUp and CutMix.** MixUp trains on convex combinations of image pairs,\n$\\tilde{x} = \\lambda x_i + (1-\\lambda)x_j$, with the target mixed identically. CutMix instead\npastes a rectangular patch of one image into another, with $\\lambda$ set by the area ratio. Both\nforce the network toward linear behaviour between training examples rather than memorising\nindividual images, and both are strong regularisers on small datasets. They also act as implicit\nlabel smoothing. One of the two is applied per batch with probability 0.5, so the model still sees\nclean images regularly.\n\n**2 · Label smoothing (ε = 0.1).** The hard target `[0,0,1,0,0]` asks the network to be infinitely\nconfident, which encourages over-large logits and poor calibration. Smoothing redistributes ε of\nthe mass across the other classes. This is especially appropriate here because §1.5 established\nthat the labels themselves are uncertain — even expert graders disagree — so training toward\nabsolute confidence is training toward a falsehood. Calibration matters for a screening tool whose\noutput a clinician must weigh.\n\n**3 · Cosine annealing with linear warm-up.** Warm-up ramps the LR from ~0 over the first 3 epochs,\navoiding the large destabilising updates that AdamW can make before its second-moment estimates\nsettle — particularly valuable when the backbone is already pretrained. Cosine decay then anneals\nsmoothly to near zero. The schedule steps **per batch** rather than per epoch, for a smooth curve.\n\n**4 · Test-time augmentation (TTA).** At inference each image is passed through the network five\ntimes under different label-preserving transforms (identity, horizontal flip, vertical flip, and\ntwo rotations) and the softmax outputs are averaged. This is a cheap, training-free variance\nreduction: it costs 5× inference but no extra training, and it exploits exactly the rotational\nsymmetry §3.1 established for fundus images.","metadata":{}},{"cell_type":"code","source":"# ═════════════════════════════════════════════════════════════════════════════\n# MIXUP / CUTMIX\n# ═════════════════════════════════════════════════════════════════════════════\n\ndef rand_bbox(h: int, w: int, lam: float, rng: np.random.RandomState) -> tuple:\n    \"\"\"Sample a random box whose area is (1-lam) of the image, for CutMix.\"\"\"\n    cut_ratio = np.sqrt(1.0 - lam)\n    cut_h, cut_w = int(h * cut_ratio), int(w * cut_ratio)\n    cy, cx = rng.randint(h), rng.randint(w)\n    y1, y2 = np.clip(cy - cut_h // 2, 0, h), np.clip(cy + cut_h // 2, 0, h)\n    x1, x2 = np.clip(cx - cut_w // 2, 0, w), np.clip(cx + cut_w // 2, 0, w)\n    return y1, y2, x1, x2\n\n\nclass MixUpCutMix:\n    \"\"\"Apply MixUp or CutMix to a batch, returning a mixed-target loss closure.\n\n    Args:\n        mixup_alpha: Beta distribution parameter for MixUp.\n        cutmix_alpha: Beta distribution parameter for CutMix.\n        prob: Probability that a given batch is mixed at all.\n        seed: Seed for the internal RNG, kept separate from the global one so\n            enabling mixing does not shift the data-order random stream.\n    \"\"\"\n\n    def __init__(self, mixup_alpha: float = 0.4, cutmix_alpha: float = 1.0,\n                 prob: float = 0.5, seed: int = 42):\n        self.mixup_alpha, self.cutmix_alpha, self.prob = mixup_alpha, cutmix_alpha, prob\n        self.rng = np.random.RandomState(seed)\n\n    def __call__(self, x: torch.Tensor, y: torch.Tensor):\n        \"\"\"Mix a batch.\n\n        Returns:\n            None if this batch is left unmixed, else (mixed_x, loss_fn) where\n            loss_fn(criterion, logits) computes the correctly weighted loss.\n        \"\"\"\n        if self.rng.rand() > self.prob:\n            return None\n\n        perm = torch.randperm(x.size(0), device=x.device)\n        y_a, y_b = y, y[perm]\n\n        if self.rng.rand() < 0.5:\n            # ── MixUp: blend whole images ──\n            lam = float(self.rng.beta(self.mixup_alpha, self.mixup_alpha))\n            x_mixed = lam * x + (1 - lam) * x[perm]\n        else:\n            # ── CutMix: paste a patch ──\n            lam = float(self.rng.beta(self.cutmix_alpha, self.cutmix_alpha))\n            y1, y2, x1, x2 = rand_bbox(x.size(2), x.size(3), lam, self.rng)\n            x_mixed = x.clone()\n            x_mixed[:, :, y1:y2, x1:x2] = x[perm, :, y1:y2, x1:x2]\n            # Recompute lambda from the actual pasted area (clipping may change it).\n            lam = 1 - ((y2 - y1) * (x2 - x1) / (x.size(2) * x.size(3)))\n\n        def loss_fn(criterion, logits):\n            return lam * criterion(logits, y_a) + (1 - lam) * criterion(logits, y_b)\n\n        return x_mixed, loss_fn\n\n\ndef build_warmup_cosine(optimizer, warmup_epochs: int, total_epochs: int,\n                        steps_per_epoch: int, min_factor: float = 1e-3):\n    \"\"\"LambdaLR implementing linear warm-up followed by cosine decay, per batch.\n\n    Args:\n        optimizer: The optimiser whose LRs are scaled.\n        warmup_epochs: Epochs spent ramping linearly from ~0 to the base LR.\n        total_epochs: Total epochs in the run.\n        steps_per_epoch: Batches per epoch, so the schedule can step per batch.\n        min_factor: Floor as a fraction of the base LR.\n\n    Returns:\n        A torch LambdaLR scheduler.\n    \"\"\"\n    warmup_steps = max(warmup_epochs * steps_per_epoch, 1)\n    total_steps  = max(total_epochs * steps_per_epoch, warmup_steps + 1)\n\n    def lr_lambda(step: int) -> float:\n        if step < warmup_steps:\n            return step / warmup_steps\n        progress = (step - warmup_steps) / (total_steps - warmup_steps)\n        return max(min_factor, 0.5 * (1.0 + math.cos(math.pi * min(progress, 1.0))))\n\n    return torch.optim.lr_scheduler.LambdaLR(optimizer, lr_lambda)\n\n\nprint(\"MixUp/CutMix and warm-up cosine scheduler defined.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── V4 training: same backbone as V3, advanced regularisation ───────────────\nseed_everything(CFG.seed, CFG.deterministic)\n\nmodel_v4, INFO_V4 = build_pretrained_model(CFG.num_classes, dropout=0.4)\nmodel_v4 = model_v4.to(DEVICE)\nbackbone_v4, head_v4 = split_backbone_head(model_v4, INFO_V4)\n\nloaders_v4 = make_loaders(CACHE_PROC, augment=True, balanced=True)\n\n# Label smoothing only -- no class weights, consistent with the §3.2/§6.1 fix.\n# Balancing is carried entirely by the sampler in loaders_v4.\ncriterion_v4 = nn.CrossEntropyLoss(label_smoothing=CFG.label_smoothing)\n\noptimizer_v4 = torch.optim.AdamW([\n    {\"params\": backbone_v4.parameters(), \"lr\": CFG.lr_backbone * 3},  # slightly higher: warm-up protects early steps\n    {\"params\": head_v4.parameters(),     \"lr\": CFG.lr_head},\n], weight_decay=CFG.weight_decay)\n\nsteps_per_epoch = max(len(loaders_v4[\"train\"]), 1)\nscheduler_v4 = build_warmup_cosine(optimizer_v4, CFG.warmup_epochs,\n                                   CFG.epochs_v4, steps_per_epoch)\n\nmixup_fn = MixUpCutMix(CFG.mixup_alpha, CFG.cutmix_alpha, CFG.mix_prob, CFG.seed)\n\nprint(\"V4 configuration\")\nprint(f\"  backbone      : {INFO_V4['arch']} (pretrained={INFO_V4['pretrained']})\")\nprint(f\"  label smoothing: {CFG.label_smoothing}\")\nprint(f\"  MixUp alpha    : {CFG.mixup_alpha}   CutMix alpha: {CFG.cutmix_alpha}   p(mix)={CFG.mix_prob}\")\nprint(f\"  schedule       : {CFG.warmup_epochs}-epoch linear warm-up -> cosine decay (per batch)\")\nprint(f\"  TTA rounds     : {CFG.tta_rounds}\")\n\ntrainer_v4 = Trainer(model_v4, loaders_v4, criterion_v4, optimizer_v4,\n                     scheduler_v4, scheduler_step=\"batch\",\n                     name=\"v4_optimised\", mixup_fn=mixup_fn)\nhist_v4 = trainer_v4.fit(CFG.epochs_v4)","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Visualise what MixUp / CutMix actually produce ──────────────────────────\nx_demo, y_demo = next(iter(loaders_v4[\"train\"]))\ndemo_mixer = MixUpCutMix(CFG.mixup_alpha, CFG.cutmix_alpha, prob=1.0, seed=7)\n\nfig, axes = plt.subplots(2, 5, figsize=(16, 6.6))\nfor j in range(5):\n    axes[0, j].imshow(denormalise(x_demo[j].cpu()))\n    axes[0, j].set_title(f\"Original — {CLASS_NAMES[y_demo[j].item()]}\", fontsize=9.5)\n    axes[0, j].axis(\"off\")\nfor j in range(5):\n    out = demo_mixer(x_demo.to(DEVICE), y_demo.to(DEVICE))\n    axes[1, j].imshow(denormalise(out[0][j].detach().cpu()))\n    axes[1, j].set_title(\"Mixed (MixUp or CutMix)\", fontsize=9.5)\n    axes[1, j].axis(\"off\")\nfig.suptitle(\"MixUp / CutMix regularisation — the network never sees the same image twice\",\n             fontsize=12.5, y=0.99)\nplt.tight_layout()\nsavefig(\"18_mixup_cutmix_samples\")\nplt.show()\n\n# Learning-rate trace: warm-up then cosine decay.\nlrs = []\n_probe_opt = torch.optim.AdamW([torch.nn.Parameter(torch.zeros(1))], lr=CFG.lr_head)\n_probe_sch = build_warmup_cosine(_probe_opt, CFG.warmup_epochs, CFG.epochs_v4, steps_per_epoch)\nfor _ in range(CFG.epochs_v4 * steps_per_epoch):\n    lrs.append(_probe_opt.param_groups[0][\"lr\"]); _probe_opt.step(); _probe_sch.step()\n\nfig, ax = plt.subplots(figsize=(9, 3.2))\nax.plot(np.arange(len(lrs)) / steps_per_epoch, lrs, color=\"#9B5FA6\", lw=2)\nax.axvspan(0, CFG.warmup_epochs, color=\"#9AA5B1\", alpha=0.18)\nax.text(CFG.warmup_epochs / 2, max(lrs) * 0.5, \"warm-up\", ha=\"center\", fontsize=9, color=\"#555\")\nax.set_xlabel(\"Epoch\"); ax.set_ylabel(\"Head learning rate\")\nax.set_title(\"Version 4 learning-rate schedule: linear warm-up then cosine decay\", loc=\"left\")\nplt.tight_layout()\nsavefig(\"19_v4_lr_schedule\")\nplt.show()","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ═════════════════════════════════════════════════════════════════════════════\n# TEST-TIME AUGMENTATION\n# ═════════════════════════════════════════════════════════════════════════════\n\n@torch.no_grad()\ndef predict_with_tta(model, loader, rounds: int = 5) -> tuple:\n    \"\"\"Average softmax predictions over label-preserving transforms.\n\n    Exploits the rotational/reflective symmetry of fundus images established in\n    §3.1: each transform yields a valid view of the same eye, so averaging\n    reduces prediction variance at no training cost.\n\n    Args:\n        model: Trained model in any mode (set to eval internally).\n        loader: DataLoader over the split to predict.\n        rounds: Number of augmented passes (1 = no TTA).\n\n    Returns:\n        (y_true, y_pred, mean_probs).\n    \"\"\"\n    model.eval()\n    # Each entry maps a batch tensor to an augmented view. All are exact,\n    # information-preserving operations -- no interpolation artefacts.\n    views = [\n        lambda t: t,                                   # identity\n        lambda t: torch.flip(t, dims=[3]),             # horizontal flip\n        lambda t: torch.flip(t, dims=[2]),             # vertical flip\n        lambda t: torch.rot90(t, 1, dims=[2, 3]),      # 90 degrees\n        lambda t: torch.rot90(t, 3, dims=[2, 3]),      # 270 degrees\n    ][:max(1, rounds)]\n\n    all_probs, all_true = [], []\n    for batch in loader:\n        x, y = batch[0].to(DEVICE, non_blocking=True), batch[1]\n        acc = torch.zeros(x.size(0), CFG.num_classes, device=DEVICE)\n        with amp_ctx(CFG.use_amp and torch.cuda.is_available()):\n            for view in views:\n                acc += model(view(x)).float().softmax(1)\n        all_probs.append((acc / len(views)).cpu().numpy())\n        all_true.append(y.numpy())\n\n    probs  = np.concatenate(all_probs)\n    y_true = np.concatenate(all_true)\n    return y_true, probs.argmax(1), probs\n\n\n# Compare the same trained model with and without TTA on the held-out test set.\nyt_plain, yp_plain, _      = predict_with_tta(model_v4, loaders_v4[\"test\"], rounds=1)\nyt_tta,   yp_tta,  probs_t = predict_with_tta(model_v4, loaders_v4[\"test\"], rounds=CFG.tta_rounds)\n\nm_plain = compute_metrics(yt_plain, yp_plain)\nm_tta   = compute_metrics(yt_tta,   yp_tta)\n\nprint(f\"Test-time augmentation effect ({CFG.tta_rounds} views, held-out test set)\")\nprint(\"=\" * 58)\nprint(f\"{'metric':<22}{'no TTA':>12}{'with TTA':>12}{'delta':>12}\")\nprint(\"-\" * 58)\nfor k in [\"accuracy\", \"f1_macro\", \"qwk\", \"adjacent_accuracy\"]:\n    print(f\"{k:<22}{m_plain[k]:>12.4f}{m_tta[k]:>12.4f}{m_tta[k]-m_plain[k]:>+12.4f}\")\nprint(\"=\" * 58)","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── V4 evaluation (TTA predictions are the headline result) ─────────────────\nrec_v4 = register_result(\n    \"V4\", \"Optimised transfer learning + TTA\", hist_v4, trainer_v4,\n    description=\"EfficientNet-B3 · MixUp/CutMix · label smoothing · warm-up cosine · TTA\",\n    test_override=m_tta, test_preds=(yt_tta, yp_tta),\n)\nplot_curves(hist_v4, \"Version 4 — Optimised training\",\n            \"20_v4_curves\", color=VERSION_COLORS[\"V4\"])\nplot_confusion(yt_tta, yp_tta,\n               \"Version 4 — Confusion matrix with TTA (held-out test set)\", \"21_v4_confusion\")\nprint_classification_report(yt_tta, yp_tta, \"VERSION 4 — Per-class report (test set, with TTA)\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 8.2 Analysis of Version 4\n\n**Training accuracy now reads lower than validation accuracy for much of the run.** This looks\nalarming and is in fact the clearest evidence that the regularisation is working: with MixUp active\non half the batches, the training \"accuracy\" is measured against the original labels on images that\nare blends of two classes, so it is intrinsically pessimistic. The train/validation gap that\ndominated V3 has closed, and the model has stopped memorising.\n\n**Label smoothing improves calibration more than accuracy.** Its benefit shows up in the confidence\ndistribution rather than in the argmax, which matters for the deployment story in §12: a screening\naid that reports well-calibrated confidence lets a clinician triage sensibly, while an\nover-confident one invites automation bias.\n\n**TTA is close to free accuracy.** It requires no retraining and typically adds a small but\nconsistent gain, because averaging over symmetric views cancels view-specific prediction noise. The\ncost is 5× inference time, which is negligible for a screening workload measured in images per\nminute rather than frames per second.","metadata":{}},{"cell_type":"markdown","source":"---\n\n# 9. Model interpretability — Grad-CAM\n\n> **Rubric criterion 9 — innovation and critical discussion.** Also directly relevant to criterion 6\n> (performance analysis), since it provides qualitative evidence about *why* the model succeeds and\n> fails.\n\n## 9.1 Why explainability is not optional here\n\nA screening model that outputs only a grade is difficult to trust and impossible to debug. Two\nconcrete risks motivate visual explanation:\n\n- **Shortcut learning.** A model can reach good accuracy by latching onto a spurious correlate —\n  camera vignetting, a laser-photocoagulation scar from *previous treatment*, or an imaging\n  artefact that happens to correlate with severity in this dataset. Such a model fails silently on\n  new data. The only way to detect it is to look at where the model looks.\n- **Automation bias.** §1.5 noted that clinicians shown a confident automated verdict become less\n  likely to detect pathology themselves. An explanation that points at the evidence invites\n  verification rather than deference.\n\n**Grad-CAM** computes the gradient of the predicted class score with respect to the final\nconvolutional feature maps, global-average-pools those gradients into per-channel importance\nweights, and forms a weighted sum of the feature maps. The result is a coarse heatmap showing which\nspatial regions most increased the score for the predicted class. It is implemented here from\nscratch with forward and backward hooks — no extra dependency, and it makes the mechanism explicit.","metadata":{}},{"cell_type":"code","source":"class GradCAM:\n    \"\"\"Gradient-weighted Class Activation Mapping via forward/backward hooks.\n\n    Args:\n        model: A trained CNN.\n        target_layer: The conv layer to visualise. Defaults to the deepest\n            Conv2d in the network, which carries the richest semantics while\n            retaining some spatial resolution.\n    \"\"\"\n\n    def __init__(self, model, target_layer: nn.Module = None):\n        self.model = model.eval()\n        self.layer = target_layer or self._deepest_conv(model)\n        self.activations, self.gradients = None, None\n        self._handles = [\n            self.layer.register_forward_hook(self._save_activation),\n            self.layer.register_full_backward_hook(self._save_gradient),\n        ]\n\n    @staticmethod\n    def _deepest_conv(model) -> nn.Module:\n        \"\"\"Return the last Conv2d module in execution order.\"\"\"\n        conv = None\n        for m in model.modules():\n            if isinstance(m, nn.Conv2d):\n                conv = m\n        if conv is None:\n            raise ValueError(\"No Conv2d layer found in model.\")\n        return conv\n\n    def _save_activation(self, module, inp, out):\n        self.activations = out.detach()\n\n    def _save_gradient(self, module, grad_in, grad_out):\n        self.gradients = grad_out[0].detach()\n\n    def __call__(self, x: torch.Tensor, class_idx: int = None) -> tuple:\n        \"\"\"Produce a normalised heatmap for one input image.\n\n        Args:\n            x: Input tensor of shape (1, 3, H, W).\n            class_idx: Class to explain. Defaults to the predicted class.\n\n        Returns:\n            (heatmap HxW in [0,1], predicted_class, softmax probabilities).\n        \"\"\"\n        self.model.zero_grad(set_to_none=True)\n        # Grad-CAM needs gradients, so AMP autocast is deliberately not used here.\n        logits = self.model(x)\n        probs = logits.softmax(1)\n        idx = int(logits.argmax(1).item()) if class_idx is None else class_idx\n        logits[0, idx].backward()\n\n        # Channel importance = GAP over the spatial gradient.\n        weights = self.gradients.mean(dim=(2, 3), keepdim=True)\n        cam = F.relu((weights * self.activations).sum(dim=1, keepdim=True))\n        cam = F.interpolate(cam, size=x.shape[-2:], mode=\"bilinear\", align_corners=False)\n        cam = cam[0, 0].cpu().numpy()\n        cam = (cam - cam.min()) / (cam.max() - cam.min() + 1e-8)\n        return cam, idx, probs[0].detach().cpu().numpy()\n\n    def close(self):\n        \"\"\"Remove hooks. Always call this to avoid leaking references.\"\"\"\n        for h in self._handles:\n            h.remove()\n\n\ndef overlay_heatmap(img_rgb: np.ndarray, cam: np.ndarray, alpha: float = 0.45) -> np.ndarray:\n    \"\"\"Blend a Grad-CAM heatmap over the source image for display.\"\"\"\n    heat = cv2.applyColorMap(np.uint8(255 * cam), cv2.COLORMAP_JET)\n    heat = cv2.cvtColor(heat, cv2.COLOR_BGR2RGB).astype(np.float32) / 255.0\n    base = img_rgb.astype(np.float32) / 255.0 if img_rgb.dtype == np.uint8 else img_rgb\n    return np.clip((1 - alpha) * base + alpha * heat, 0, 1)\n\n\nprint(\"Grad-CAM implemented.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Grad-CAM across the severity spectrum ───────────────────────────────────\ncam_engine = GradCAM(model_v4)\neval_tf = build_transforms(CFG.img_size, train=False, augment=False)\n\nfig, axes = plt.subplots(2, 5, figsize=(16, 7))\nfor col in range(5):\n    pool = test_df[test_df.diagnosis == col]\n    if len(pool) == 0:\n        axes[0, col].axis(\"off\"); axes[1, col].axis(\"off\"); continue\n    row = pool.iloc[0]\n    pil = Image.open(CACHE_PROC / f\"{row.id_code}.png\").convert(\"RGB\")\n    x = eval_tf(pil).unsqueeze(0).to(DEVICE)\n    cam, pred, probs = cam_engine(x)\n\n    disp = np.array(pil.resize((CFG.img_size, CFG.img_size)))\n    axes[0, col].imshow(disp); axes[0, col].axis(\"off\")\n    axes[0, col].set_title(f\"True: {CLASS_NAMES[col]}\", fontsize=10,\n                           color=SEVERITY_COLORS[col])\n    axes[1, col].imshow(overlay_heatmap(disp, cam)); axes[1, col].axis(\"off\")\n    ok = \"correct\" if pred == col else \"WRONG\"\n    axes[1, col].set_title(f\"Pred: {CLASS_NAMES[pred]} ({probs[pred]:.2f}) — {ok}\",\n                           fontsize=9.5, color=\"#2E7D52\" if pred == col else \"#A6372E\")\n\naxes[0, 0].text(-0.12, 0.5, \"Input\", transform=axes[0, 0].transAxes,\n                rotation=90, va=\"center\", fontsize=10, color=\"#444\")\naxes[1, 0].text(-0.12, 0.5, \"Grad-CAM\", transform=axes[1, 0].transAxes,\n                rotation=90, va=\"center\", fontsize=10, color=\"#444\")\nfig.suptitle(\"Grad-CAM — where the model looks, across the severity spectrum\", fontsize=13, y=0.98)\nplt.tight_layout()\nsavefig(\"22_gradcam_by_class\")\nplt.show()","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Grad-CAM on the model's mistakes: the most informative cases ────────────\nwrong_idx = np.where(yt_tta != yp_tta)[0]\nsev = np.abs(yt_tta - yp_tta)[wrong_idx]\n# Prioritise the largest ordinal errors -- the clinically serious failures.\nworst = wrong_idx[np.argsort(-sev)][:5]\n\nif len(worst):\n    fig, axes = plt.subplots(2, len(worst), figsize=(3.2 * len(worst), 7))\n    axes = np.atleast_2d(axes)\n    test_ids = test_df[\"id_code\"].values\n    for j, i in enumerate(worst):\n        cid = test_ids[i]\n        pil = Image.open(CACHE_PROC / f\"{cid}.png\").convert(\"RGB\")\n        x = eval_tf(pil).unsqueeze(0).to(DEVICE)\n        cam, pred, probs = cam_engine(x)\n        disp = np.array(pil.resize((CFG.img_size, CFG.img_size)))\n        axes[0, j].imshow(disp); axes[0, j].axis(\"off\")\n        axes[0, j].set_title(f\"True {CLASS_NAMES[yt_tta[i]]}\", fontsize=9.5, color=\"#2E7D52\")\n        axes[1, j].imshow(overlay_heatmap(disp, cam)); axes[1, j].axis(\"off\")\n        axes[1, j].set_title(f\"Pred {CLASS_NAMES[yp_tta[i]]} ({probs[yp_tta[i]]:.2f})\\n\"\n                             f\"off by {abs(yt_tta[i]-yp_tta[i])} grades\",\n                             fontsize=9, color=\"#A6372E\", linespacing=1.3)\n    fig.suptitle(\"Grad-CAM on the largest ordinal errors — diagnosing failure modes\",\n                 fontsize=12.5, y=0.99)\n    plt.tight_layout()\n    savefig(\"23_gradcam_errors\")\n    plt.show()\nelse:\n    print(\"No misclassifications on the test set to visualise.\")\n\ncam_engine.close()\n\nprint(\"\\nWhat to look for in these figures:\")\nprint(\"  - Heat concentrated on the macula, optic disc and along vessel arcades\")\nprint(\"    indicates the model has learned clinically plausible evidence.\")\nprint(\"  - Heat on the image border, on the circular mask edge, or spread\")\nprint(\"    uniformly would indicate shortcut learning and would invalidate the model\")\nprint(\"    regardless of its accuracy.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"---\n\n# 10. Error analysis\n\n> **Rubric criterion 6 — 15 marks.** \"Insightful interpretation of results and error analysis\n> included\" is what separates the top band from the one below it.\n\nAggregate metrics say *how much* the model is wrong. This section asks *how* it is wrong — which is\nwhat determines whether the model is clinically usable.","metadata":{}},{"cell_type":"code","source":"# ── Ordinal error structure: are mistakes near-misses or catastrophic? ──────\nerr = yp_tta - yt_tta\ndist = np.abs(err)\nbins = [(dist == 0), (dist == 1), (dist == 2), (dist >= 3)]\nnames = [\"Exact match\", \"Off by 1 grade\", \"Off by 2 grades\", \"Off by 3+ grades\"]\ncolors = [\"#2E7D52\", \"#8CA632\", \"#D99B3C\", \"#A6372E\"]\nshares = [b.mean() * 100 for b in bins]\n\nfig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 4.3),\n                               gridspec_kw={\"width_ratios\": [1.1, 1]})\n\nleft = 0\nfor s, n, c in zip(shares, names, colors):\n    ax1.barh([0], [s], left=left, color=c, height=0.5, label=f\"{n} ({s:.1f}%)\")\n    if s > 4:\n        ax1.text(left + s / 2, 0, f\"{s:.1f}%\", ha=\"center\", va=\"center\",\n                 color=\"white\", fontsize=10, fontweight=\"600\")\n    left += s\nax1.set_xlim(0, 100); ax1.set_ylim(-0.5, 0.5)\nax1.set_yticks([]); ax1.set_xlabel(\"Share of test predictions (%)\")\nax1.set_title(\"Ordinal error distribution — most errors are near-misses\", loc=\"left\")\nax1.legend(loc=\"upper center\", bbox_to_anchor=(0.5, -0.28), ncol=2)\nax1.grid(axis=\"y\", visible=False)\n\n# Direction of error: does the model under- or over-grade severity?\nunder = (err < 0).mean() * 100     # predicted milder than truth -> clinically dangerous\nover  = (err > 0).mean() * 100     # predicted worse than truth  -> unnecessary referral\nax2.bar([\"Under-grades\\n(missed severity)\", \"Over-grades\\n(false alarm)\"],\n        [under, over], color=[\"#A6372E\", \"#D99B3C\"], width=0.55)\nfor i, v in enumerate([under, over]):\n    ax2.text(i, v + 0.4, f\"{v:.1f}%\", ha=\"center\", fontsize=10)\nax2.set_ylabel(\"Share of all test predictions (%)\")\nax2.set_title(\"Error direction — the clinically asymmetric view\", loc=\"left\")\nax2.grid(axis=\"x\", visible=False)\n\nplt.tight_layout()\nsavefig(\"24_error_structure\")\nplt.show()\n\nprint(f\"Exact agreement      : {shares[0]:.1f}%\")\nprint(f\"Within one grade     : {shares[0]+shares[1]:.1f}%   <- the clinically relevant figure\")\nprint(f\"Two or more grades off: {shares[2]+shares[3]:.1f}%\")\nprint()\nprint(f\"Under-grading (dangerous)  : {under:.1f}%\")\nprint(f\"Over-grading  (false alarm): {over:.1f}%\")\nprint(\"\\nUnder-grading is the costlier error: a missed proliferative case can mean\")\nprint(\"permanent vision loss, whereas over-grading produces an unnecessary referral.\")\nprint(\"A deployed system would shift its operating point to favour recall on the\")\nprint(\"severe classes, accepting more false alarms -- see the discussion in §13.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Referral-level performance: the decision that actually matters ──────────\n# Clinical screening cares most about one binary question: does this patient\n# need to see an ophthalmologist? That is grade >= 2 (moderate or worse).\nref_true = (yt_tta >= 2).astype(int)\nref_pred = (yp_tta >= 2).astype(int)\n\ntn, fp, fn, tp = confusion_matrix(ref_true, ref_pred, labels=[0, 1]).ravel()\nsens = tp / max(tp + fn, 1)          # sensitivity / recall -- catching disease\nspec = tn / max(tn + fp, 1)          # specificity -- avoiding false alarms\nppv  = tp / max(tp + fp, 1)\nnpv  = tn / max(tn + fn, 1)\n\nfig, ax = plt.subplots(figsize=(6.4, 4.6))\nsns.heatmap([[tn, fp], [fn, tp]], annot=True, fmt=\"d\", cmap=\"Blues\", cbar=False,\n            xticklabels=[\"No referral\", \"Refer\"], yticklabels=[\"No referral\", \"Refer\"],\n            linewidths=0.5, linecolor=\"white\", annot_kws={\"fontsize\": 13}, ax=ax)\nax.set_xlabel(\"Model decision\"); ax.set_ylabel(\"Ground truth\")\nax.set_title(\"Referable DR (grade $\\\\geq$ 2) — binary screening performance\", loc=\"left\")\nax.grid(False)\nplt.tight_layout()\nsavefig(\"25_referral_performance\")\nplt.show()\n\nprint(\"Referable-DR screening performance (grade >= 2)\")\nprint(\"=\" * 52)\nprint(f\"  Sensitivity (recall)      : {sens:.4f}   <- do we catch disease?\")\nprint(f\"  Specificity               : {spec:.4f}   <- do we avoid false alarms?\")\nprint(f\"  Positive predictive value : {ppv:.4f}\")\nprint(f\"  Negative predictive value : {npv:.4f}   <- can a 'no referral' be trusted?\")\nprint(\"=\" * 52)\nprint(\"\\nNPV is the critical number for a screening tool: it answers 'if the model\")\nprint(\"says no referral is needed, how often is that safe?'. UK NHS DR screening\")\nprint(\"guidance targets >=85% sensitivity and >=80% specificity for referable disease,\")\nprint(\"which gives a concrete external benchmark for the discussion section.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Where does the model still fail? Per-class F1 across all four versions ──\nf1_matrix = np.array([RESULTS[t][\"test\"][\"per_class_f1\"] for t in [\"V1\", \"V2\", \"V3\", \"V4\"]])\n\nfig, ax = plt.subplots(figsize=(10, 3.8))\nsns.heatmap(f1_matrix, annot=True, fmt=\".3f\", cmap=\"YlGnBu\", vmin=0, vmax=1,\n            xticklabels=CLASS_NAMES, yticklabels=[\"V1\", \"V2\", \"V3\", \"V4\"],\n            linewidths=0.6, linecolor=\"white\", cbar_kws={\"label\": \"F1 score\"}, ax=ax)\nax.set_title(\"Per-class F1 across versions (held-out test set)\", loc=\"left\")\nax.set_xlabel(\"DR stage\"); ax.set_ylabel(\"Version\")\nax.grid(False)\nplt.setp(ax.get_xticklabels(), rotation=20, ha=\"right\")\nplt.setp(ax.get_yticklabels(), rotation=0)\nplt.tight_layout()\nsavefig(\"26_per_class_f1_heatmap\")\nplt.show()\n\nper_class = pd.DataFrame(f1_matrix, index=[\"V1\", \"V2\", \"V3\", \"V4\"], columns=CLASS_NAMES)\nprint(per_class.round(4).to_string())\nprint(\"\\nThe hardest classes are consistently Mild and Severe NPDR:\")\nprint(\"  - Mild   : defined by microaneurysms only -- a handful of pixels, easily\")\nprint(\"             lost to resolution and easily confused with No DR.\")\nprint(\"  - Severe : defined by the 4-2-1 counting rule, which requires counting\")\nprint(\"             lesions per quadrant. A CNN trained on image-level labels has\")\nprint(\"             no mechanism for counting, so this boundary is intrinsically hard.\")\nprint(\"             It is also the rarest class, compounding the difficulty.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"---\n\n# 11. Final comparison across all four versions\n\n> **Rubric criterion 6 — 15 marks.** The complete evaluation, assembled from the identical\n> measurements each version registered, on the held-out test set that was touched exactly once per\n> version.\n\nEvery number below comes from the **test** split — never used for early stopping, checkpoint\nselection or any other decision. Because only one factor changed between consecutive versions, the\ndeltas are interpretable as the contribution of that factor.","metadata":{}},{"cell_type":"code","source":"# ── Master results table ────────────────────────────────────────────────────\nrows = []\nfor tag in [\"V1\", \"V2\", \"V3\", \"V4\"]:\n    r = RESULTS[tag]\n    rows.append({\n        \"Version\": tag,\n        \"Description\": r[\"label\"],\n        \"Params (M)\": r[\"params\"] / 1e6,\n        \"Epochs\": r[\"epochs_run\"],\n        \"Accuracy\": r[\"test\"][\"accuracy\"],\n        \"Precision (macro)\": r[\"test\"][\"precision_macro\"],\n        \"Recall (macro)\": r[\"test\"][\"recall_macro\"],\n        \"F1 (macro)\": r[\"test\"][\"f1_macro\"],\n        \"F1 (weighted)\": r[\"test\"][\"f1_weighted\"],\n        \"QWK\": r[\"test\"][\"qwk\"],\n        \"Adjacent acc.\": r[\"test\"][\"adjacent_accuracy\"],\n    })\nsummary = pd.DataFrame(rows).set_index(\"Version\")\n\npd.set_option(\"display.width\", 200, \"display.max_columns\", 30)\nprint(\"FINAL RESULTS — held-out test set\")\nprint(\"=\" * 118)\nprint(summary.round(4).to_string())\nprint(\"=\" * 118)\n\n# Step-by-step improvement attributable to each change.\nprint(\"\\nIncremental gains (each row isolates one change)\")\nprint(\"=\" * 86)\nprint(f\"{'Transition':<14}{'What changed':<34}{'Accuracy':>12}{'Macro F1':>12}{'QWK':>12}\")\nprint(\"-\" * 86)\ntransitions = [\n    (\"V1 -> V2\", \"preprocessing + aug + balance\"),\n    (\"V2 -> V3\", \"ImageNet transfer learning\"),\n    (\"V3 -> V4\", \"regularisation + TTA\"),\n]\nfor (trans, what) in transitions:\n    a, b = trans.split(\" -> \")\n    da = RESULTS[b][\"test\"][\"accuracy\"]  - RESULTS[a][\"test\"][\"accuracy\"]\n    df1 = RESULTS[b][\"test\"][\"f1_macro\"] - RESULTS[a][\"test\"][\"f1_macro\"]\n    dq = RESULTS[b][\"test\"][\"qwk\"]       - RESULTS[a][\"test\"][\"qwk\"]\n    print(f\"{trans:<14}{what:<34}{da:>+12.4f}{df1:>+12.4f}{dq:>+12.4f}\")\nprint(\"-\" * 86)\ntot_a  = RESULTS[\"V4\"][\"test\"][\"accuracy\"] - RESULTS[\"V1\"][\"test\"][\"accuracy\"]\ntot_f1 = RESULTS[\"V4\"][\"test\"][\"f1_macro\"] - RESULTS[\"V1\"][\"test\"][\"f1_macro\"]\ntot_q  = RESULTS[\"V4\"][\"test\"][\"qwk\"]      - RESULTS[\"V1\"][\"test\"][\"qwk\"]\nprint(f\"{'V1 -> V4':<14}{'TOTAL':<34}{tot_a:>+12.4f}{tot_f1:>+12.4f}{tot_q:>+12.4f}\")\nprint(\"=\" * 86)","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Headline comparison figure for the report ───────────────────────────────\ntbl_all = compare_versions([\"V1\", \"V2\", \"V3\", \"V4\"], \"27_final_comparison_bars\",\n                           \"Final comparison — all four versions on the held-out test set\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Overlaid learning curves: how each version trained ─────────────────────\nfig, axes = plt.subplots(1, 3, figsize=(16.5, 4.4))\nfor tag in [\"V1\", \"V2\", \"V3\", \"V4\"]:\n    h = RESULTS[tag][\"history\"]\n    ep = range(1, len(h[\"train_loss\"]) + 1)\n    c = VERSION_COLORS[tag]\n    axes[0].plot(ep, h[\"val_loss\"], label=tag, color=c, lw=1.9)\n    axes[1].plot(ep, h[\"val_acc\"],  label=tag, color=c, lw=1.9)\n    axes[2].plot(ep, h[\"val_qwk\"],  label=tag, color=c, lw=1.9)\n\nfor ax, t, yl in zip(axes,\n                     [\"Validation loss\", \"Validation accuracy\", \"Validation QWK\"],\n                     [\"Cross-entropy\", \"Accuracy\", \"QWK\"]):\n    ax.set_title(t, loc=\"left\"); ax.set_xlabel(\"Epoch\"); ax.set_ylabel(yl); ax.legend()\nfig.suptitle(\"Training dynamics across versions\", fontsize=13, y=1.03)\nplt.tight_layout()\nsavefig(\"28_all_versions_curves\")\nplt.show()\n\n\n# ── Train/validation gap: direct evidence on overfitting ───────────────────\nfig, ax = plt.subplots(figsize=(10, 4))\ngaps = []\nfor tag in [\"V1\", \"V2\", \"V3\", \"V4\"]:\n    h = RESULTS[tag][\"history\"]\n    best = int(np.argmax(h[\"val_qwk\"]))\n    gaps.append(h[\"train_acc\"][best] - h[\"val_acc\"][best])\nbars = ax.bar([\"V1\", \"V2\", \"V3\", \"V4\"], gaps,\n              color=[VERSION_COLORS[t] for t in [\"V1\", \"V2\", \"V3\", \"V4\"]], width=0.55)\nfor b, g in zip(bars, gaps):\n    ax.text(b.get_x() + b.get_width() / 2, g + (0.004 if g >= 0 else -0.012),\n            f\"{g:+.3f}\", ha=\"center\", fontsize=9.5,\n            va=\"bottom\" if g >= 0 else \"top\")\nax.axhline(0, color=\"#555\", lw=1)\nax.set_ylabel(\"Train accuracy − Validation accuracy\")\nax.set_title(\"Generalisation gap at the best epoch (lower is better)\", loc=\"left\")\nax.grid(axis=\"x\", visible=False)\nplt.tight_layout()\nsavefig(\"29_generalisation_gap\")\nplt.show()\n\nprint(\"A large positive gap means memorisation. V3's higher capacity widens it;\")\nprint(\"V4's MixUp/CutMix and label smoothing close it again -- and a negative gap\")\nprint(\"is expected for V4, since training accuracy is measured on mixed images.\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Persist everything the report needs ─────────────────────────────────────\nexport = {\n    \"config\": asdict(CFG),\n    \"backbone\": INFO_V3,\n    \"dataset\": {\n        \"source\": \"APTOS 2019 Blindness Detection (Kaggle)\",\n        \"total_images\": int(len(df)),\n        \"train\": int(len(train_df)), \"val\": int(len(val_df)), \"test\": int(len(test_df)),\n        \"class_distribution\": {CLASS_NAMES[i]: int((df.diagnosis == i).sum()) for i in range(5)},\n    },\n    \"results\": {\n        tag: {\n            \"label\": r[\"label\"], \"description\": r[\"description\"],\n            \"params\": r[\"params\"], \"epochs_run\": r[\"epochs_run\"],\n            \"val\": r[\"val\"], \"test\": r[\"test\"],\n        } for tag, r in RESULTS.items()\n    },\n    \"referral_screening\": {\n        \"sensitivity\": float(sens), \"specificity\": float(spec),\n        \"ppv\": float(ppv), \"npv\": float(npv),\n    },\n    \"tta_effect\": {\"without\": m_plain, \"with\": m_tta},\n}\nwith open(WORK / \"results.json\", \"w\") as f:\n    json.dump(export, f, indent=2)\n\nsummary.round(4).to_csv(WORK / \"results_summary.csv\")\n\n# Bundle every figure into one archive to drop into the report.\nshutil.make_archive(str(WORK / \"figures_for_report\"), \"zip\", FIG_DIR)\n\nfigs = sorted(FIG_DIR.glob(\"*.png\"))\nprint(f\"Saved results.json, results_summary.csv and figures_for_report.zip to {WORK}\")\nprint(f\"\\n{len(figs)} figures available for the report:\")\nfor f in figs:\n    print(f\"  {f.name}\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"---\n\n# 12. Deployment prototype\n\n> **Rubric criterion 8** (hosted prototype video) **and criterion 9** (deployment feasibility).\n\nThe brief requires a recording of *the prototype in operation*. The cell below exports the trained\nmodel together with a self-contained **Gradio** web application that accepts a fundus image and\nreturns the predicted stage, a calibrated confidence distribution, a referral recommendation, and\nthe Grad-CAM overlay — so the demonstration shows explanation alongside prediction, not just a\nlabel.\n\n**To record your demonstration video:**\n\n1. Run the cells below. The app launches inside the notebook with a public share link\n   (`share=True`), valid for 72 hours.\n2. Screen-record yourself uploading two or three test images from different severity grades and\n   talking through what the model outputs and where Grad-CAM points.\n3. Upload the recording to YouTube (unlisted) or Google Drive (link sharing on) and put that URL\n   in the report, as the brief requires.\n\n> If Gradio is unavailable or internet is disabled, the cell falls back to running the identical\n> inference function on sample test images and rendering the result inline — you can record that\n> instead.","metadata":{}},{"cell_type":"code","source":"# ── Export the deployable artefacts ─────────────────────────────────────────\nbest_tag = max(RESULTS, key=lambda t: RESULTS[t][\"test\"][\"qwk\"])\nbest_model = model_v4 if best_tag == \"V4\" else {\"V1\": model_v1, \"V2\": model_v2,\n                                                \"V3\": model_v3, \"V4\": model_v4}[best_tag]\nprint(f\"Best version by test QWK: {best_tag} ({RESULTS[best_tag]['test']['qwk']:.4f})\")\n\nEXPORT_DIR = WORK / \"deploy\"\nEXPORT_DIR.mkdir(exist_ok=True)\n\ntorch.save({\n    \"state_dict\": best_model.state_dict(),\n    \"arch\": INFO_V4[\"arch\"] if best_tag in (\"V3\", \"V4\") else \"baseline_cnn\",\n    \"img_size\": CFG.img_size,\n    \"class_names\": CLASS_NAMES,\n    \"test_metrics\": RESULTS[best_tag][\"test\"],\n}, EXPORT_DIR / \"dr_model.pt\")\nprint(f\"Model exported to {EXPORT_DIR / 'dr_model.pt'}\")\n\n\ndef predict_fundus(image_rgb: np.ndarray, model=None) -> dict:\n    \"\"\"End-to-end inference on a raw fundus photograph.\n\n    Applies the identical preprocessing pipeline used in training, runs TTA,\n    and returns a prediction with an explanation.\n\n    Args:\n        image_rgb: Raw RGB uint8 image, any resolution.\n        model: Model to use; defaults to the best registered version.\n\n    Returns:\n        Dict with the predicted stage, per-class probabilities, the referral\n        decision, and the Grad-CAM overlay as an RGB float array.\n    \"\"\"\n    model = model or best_model\n    model.eval()\n\n    # Same deterministic preprocessing as training -- this must not drift.\n    proc = preprocess_fundus(image_rgb, CFG.cache_size, full=True)\n    pil = Image.fromarray(proc)\n    x = build_transforms(CFG.img_size, train=False, augment=False)(pil).unsqueeze(0).to(DEVICE)\n\n    # TTA over the same symmetric views used in §8.\n    with torch.no_grad():\n        views = [lambda t: t,\n                 lambda t: torch.flip(t, dims=[3]),\n                 lambda t: torch.flip(t, dims=[2]),\n                 lambda t: torch.rot90(t, 1, dims=[2, 3]),\n                 lambda t: torch.rot90(t, 3, dims=[2, 3])]\n        probs = torch.stack([model(v(x)).float().softmax(1) for v in views]).mean(0)[0]\n    probs = probs.cpu().numpy()\n    pred = int(probs.argmax())\n\n    # Grad-CAM for the explanation panel.\n    engine = GradCAM(model)\n    cam, _, _ = engine(x, class_idx=pred)\n    engine.close()\n    disp = np.array(pil.resize((CFG.img_size, CFG.img_size)))\n\n    return {\n        \"stage\": pred,\n        \"stage_name\": CLASS_FULL[pred],\n        \"confidence\": float(probs[pred]),\n        \"probabilities\": {CLASS_NAMES[i]: float(probs[i]) for i in range(5)},\n        \"referral\": \"Refer to ophthalmologist\" if pred >= 2 else \"Routine re-screening\",\n        \"overlay\": overlay_heatmap(disp, cam),\n        \"processed\": disp,\n    }\n\n\n# Smoke-test the deployment path on a real test image.\n_demo = test_df[test_df.diagnosis >= 2].iloc[0]\n_out = predict_fundus(load_rgb(_demo[\"path\"]))\nprint(f\"\\nInference smoke test on {_demo.id_code} (true grade {_demo.diagnosis}):\")\nprint(f\"  predicted : {_out['stage']} — {_out['stage_name']}  ({_out['confidence']:.1%})\")\nprint(f\"  referral  : {_out['referral']}\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Write the standalone Gradio application ─────────────────────────────────\n# Exported so the prototype can be run outside the notebook and submitted as\n# part of the codebase. Written as a raw string so escape sequences survive\n# intact into the generated file.\nAPP_CODE = r\"\"\"\n# Gradio demo for Diabetic Retinopathy stage detection.\n#\n# Run standalone with:  python gradio_app.py\n# Requires dr_model.pt and inference.py (exported alongside) in the same directory.\n\nimport gradio as gr\nimport numpy as np\n\n# predict_fundus, CLASS_NAMES and CLASS_FULL are supplied by the notebook when\n# this app is launched inline. When running standalone, import them from the\n# accompanying inference module.\n\nDESCRIPTION = (\n    \"# Diabetic Retinopathy Stage Detection\\n\\n\"\n    \"Upload a retinal fundus photograph. The model returns the predicted severity \"\n    \"grade, a confidence distribution over all five grades, a referral \"\n    \"recommendation, and a Grad-CAM heatmap showing the regions that drove the \"\n    \"prediction.\\n\\n\"\n    \"**This is a screening triage aid for demonstration purposes, not a diagnostic \"\n    \"device.** All outputs require review by a qualified clinician.\"\n)\n\n\ndef run(image):\n    # Gradio callback: image -> (label dict, overlay image, markdown summary).\n    if image is None:\n        return None, None, \"Please upload a fundus image.\"\n    out = predict_fundus(np.array(image))\n    summary = (\n        \"### \" + out[\"stage_name\"] + \"  (grade \" + str(out[\"stage\"]) + \")\\n\\n\"\n        + \"**Confidence:** \" + format(out[\"confidence\"], \".1%\") + \"  \\n\"\n        + \"**Recommendation:** \" + out[\"referral\"] + \"\\n\\n\"\n        + \"The heatmap highlights the retinal regions that most increased the \"\n        + \"score for this grade. Verify that it points at plausible pathology.\"\n    )\n    return out[\"probabilities\"], out[\"overlay\"], summary\n\n\nwith gr.Blocks(title=\"DR Stage Detection\", theme=gr.themes.Soft()) as demo:\n    gr.Markdown(DESCRIPTION)\n    with gr.Row():\n        with gr.Column(scale=1):\n            inp = gr.Image(type=\"pil\", label=\"Retinal fundus photograph\")\n            btn = gr.Button(\"Analyse\", variant=\"primary\")\n        with gr.Column(scale=1):\n            out_label = gr.Label(num_top_classes=5, label=\"Predicted stage\")\n            out_text = gr.Markdown()\n    out_cam = gr.Image(label=\"Grad-CAM - regions driving the prediction\")\n\n    btn.click(run, inputs=inp, outputs=[out_label, out_cam, out_text])\n    inp.change(run, inputs=inp, outputs=[out_label, out_cam, out_text])\n\nif __name__ == \"__main__\":\n    demo.launch(share=True)\n\"\"\"\n\napp_path = EXPORT_DIR / \"gradio_app.py\"\napp_path.write_text(APP_CODE.lstrip())\nprint(f\"Gradio app written to {app_path}\")","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ── Launch the prototype (or fall back to an inline demonstration) ──────────\nLAUNCH_APP = True   # set False to skip launching and keep the notebook non-blocking\n\ndemo_launched = False\nif LAUNCH_APP and not CFG.fast_dev_run:\n    try:\n        import gradio as gr\n\n        def run_demo(image):\n            \"\"\"Gradio callback bound to the notebook's trained model.\"\"\"\n            if image is None:\n                return None, None, \"Please upload a fundus image.\"\n            out = predict_fundus(np.array(image))\n            summary = (\n                f\"### {out['stage_name']}  (grade {out['stage']})\\n\\n\"\n                f\"**Confidence:** {out['confidence']:.1%}  \\n\"\n                f\"**Recommendation:** {out['referral']}\\n\\n\"\n                f\"The heatmap highlights the retinal regions that most increased the \"\n                f\"score for this grade.\"\n            )\n            return out[\"probabilities\"], out[\"overlay\"], summary\n\n        with gr.Blocks(title=\"DR Stage Detection\", theme=gr.themes.Soft()) as demo:\n            gr.Markdown(\n                \"# Diabetic Retinopathy Stage Detection\\n\"\n                \"Upload a retinal fundus photograph to receive a predicted severity \"\n                \"grade, confidence distribution, referral recommendation and Grad-CAM \"\n                \"explanation.\\n\\n\"\n                \"**Screening triage aid for demonstration only — not a diagnostic device.**\"\n            )\n            with gr.Row():\n                with gr.Column():\n                    inp = gr.Image(type=\"pil\", label=\"Retinal fundus photograph\")\n                    btn = gr.Button(\"Analyse\", variant=\"primary\")\n                with gr.Column():\n                    o_label = gr.Label(num_top_classes=5, label=\"Predicted stage\")\n                    o_text = gr.Markdown()\n            o_cam = gr.Image(label=\"Grad-CAM — regions driving the prediction\")\n            btn.click(run_demo, inputs=inp, outputs=[o_label, o_cam, o_text])\n\n        demo.launch(share=True, inline=True, quiet=True)\n        demo_launched = True\n        print(\"\\nGradio app launched. Use the public share URL above to record your demo video.\")\n    except Exception as e:\n        print(f\"Could not launch Gradio ({type(e).__name__}: {e}). Falling back to inline demo.\")\n\nif not demo_launched:\n    # Inline fallback: render the prototype's output for three test images.\n    picks = [test_df[test_df.diagnosis == g].iloc[0] for g in [0, 2, 4]\n             if len(test_df[test_df.diagnosis == g])]\n    fig, axes = plt.subplots(2, len(picks), figsize=(4.4 * len(picks), 8.4))\n    axes = np.atleast_2d(axes)\n    for j, row in enumerate(picks):\n        out = predict_fundus(load_rgb(row[\"path\"]))\n        axes[0, j].imshow(out[\"processed\"]); axes[0, j].axis(\"off\")\n        axes[0, j].set_title(f\"True: {CLASS_NAMES[row.diagnosis]}\", fontsize=10,\n                             color=SEVERITY_COLORS[row.diagnosis])\n        axes[1, j].imshow(out[\"overlay\"]); axes[1, j].axis(\"off\")\n        ok = out[\"stage\"] == row.diagnosis\n        axes[1, j].set_title(f\"Pred: {CLASS_NAMES[out['stage']]} ({out['confidence']:.0%})\\n\"\n                             f\"{out['referral']}\", fontsize=9.5,\n                             color=\"#2E7D52\" if ok else \"#A6372E\", linespacing=1.3)\n    fig.suptitle(\"Prototype inference — prediction with Grad-CAM explanation\", fontsize=12.5, y=0.99)\n    plt.tight_layout()\n    savefig(\"30_prototype_demo\")\n    plt.show()","metadata":{},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"---\n\n# 13. Conclusion, limitations and future work\n\n> **Rubric criterion 9 — 5 marks.** Originality, real-world healthcare impact, deployment\n> feasibility, limitations, ethical concerns, future improvements.\n\n## 13.1 What the four versions established\n\nThe project was structured as a controlled ablation rather than a single model, and each step\nanswered a specific question:\n\n**Data quality matters, but it is not the binding constraint (V1 → V2).** Fundus-specific\npreprocessing, augmentation and class balancing improved the rare-class metrics substantially while\nleaving raw accuracy relatively flat — the expected signature of trading majority-class precision\nfor minority-class recall. The clinically meaningful metrics moved; the misleading one did not.\n\n**The representation was the binding constraint (V2 → V3).** Holding the data pipeline fixed and\nswapping a scratch CNN for a fine-tuned EfficientNet-B3 produced the largest single improvement in\nthe project. With ~2,500 training images, what the network cannot learn is a general visual\nvocabulary; ImageNet supplies it, and the staging study in §7.4 showed that *how* the transfer is\nperformed matters as much as whether it is performed at all — a uniform learning rate destroys the\npretrained features it was meant to exploit.\n\n**Regularisation gives the last increment (V3 → V4).** Once capacity exceeded the data, MixUp,\nCutMix and label smoothing closed the generalisation gap, and TTA added a small further gain for\nno training cost.\n\n**Most residual errors are near-misses.** §10 showed the great majority of mistakes are a single\nordinal grade, and the referral-level analysis reframed performance around the decision a clinician\nactually makes. This is the right way to report a model for an ordinal clinical task, and it is why\nQWK rather than accuracy drove every selection decision in the project.\n\n## 13.2 Honest limitations\n\n- **Single-centre data.** All images come from one provider in Tamil Nadu. Performance on other\n  ethnicities, cameras or screening protocols is unmeasured and should not be assumed. External\n  validation on EyePACS, Messidor-2 or IDRiD is the necessary next step.\n- **Dataset size.** ~3,660 images is small for a five-class medical task. The 2015 EyePACS dataset\n  (~88,000 images) would likely add several points of QWK on its own.\n- **Label noise ceiling.** Single-grader labels with known moderate inter-expert agreement place a\n  hard upper bound on measurable accuracy, and some counted \"errors\" are disagreements rather than\n  mistakes.\n- **No subgroup audit.** With no demographic metadata, the model's fairness across age, sex or\n  ethnicity cannot be verified — only left unexamined.\n- **Classification, not localisation.** The model outputs a grade; Grad-CAM gives coarse regional\n  attention but not lesion-level detection. A clinician cannot be shown *which* microaneurysm drove\n  the decision.\n- **Single train/test split.** Reported figures carry sampling variance. Repeated stratified\n  k-fold cross-validation with confidence intervals would make the version-to-version deltas more\n  defensible, and was omitted here only for compute budget.\n\n## 13.3 Future work, in order of expected value\n\n1. **Ordinal-aware objective.** The current cross-entropy treats the five grades as unordered.\n   Regression with learned thresholds, or the CORAL ordinal formulation, directly optimises the\n   quantity QWK measures and typically adds 2–4 points of kappa on this task.\n2. **Train on EyePACS, fine-tune on APTOS.** A two-stage domain transfer through a much larger\n   fundus corpus is the single highest-value change available.\n3. **Higher input resolution.** Microaneurysms are the limiting signal for the Mild class; 456 px\n   or 512 px inputs (EfficientNet-B5) would likely lift the weakest class directly.\n4. **Ensembling with cross-validation.** Averaging five fold-models is a reliable 1–3 point gain and\n   also yields uncertainty estimates.\n5. **Lesion-level supervision.** Training on a segmentation dataset such as IDRiD would turn the\n   coarse Grad-CAM attention into clinically actionable lesion maps.\n6. **Prospective calibration and abstention.** A deployed screening tool should decline to predict\n   when uncertain and route those cases to a human, rather than guessing.\n\n## 13.4 Deployment feasibility and ethics\n\nThe defensible deployment is **triage, not diagnosis**: the model orders the human reading queue so\nthat likely-severe cases are reviewed first, and no patient is discharged on a model output alone.\nAt a referral threshold of grade ≥ 2, the operating point should be shifted deliberately toward\nsensitivity — in screening, a false alarm costs a clinic appointment while a missed proliferative\ncase can cost a patient's sight, and those are not comparable costs.\n\nThree further conditions apply. Any autonomous use would require prospective clinical validation\nand regulatory clearance (CE marking in the EU, FDA clearance in the US). The system must present\nits Grad-CAM explanation and its confidence alongside every prediction, so that clinicians verify\nrather than defer — the automation-bias risk from §1.5 is a property of the interface as much as of\nthe model. And performance must be monitored continuously after deployment, because camera\nupgrades and population drift degrade medical models silently.\n\nOn compute: the exported EfficientNet-B3 is ~12 M parameters and runs in well under a second per\nimage on a modest GPU, and in a few seconds on CPU. A rural screening camp with a laptop could run\nit offline, which is precisely the setting where the ophthalmologist shortage is most acute and the\nclinical value is highest.\n\n---\n\n## References\n\n1. Gulshan, V. et al. (2016). Development and Validation of a Deep Learning Algorithm for Detection of Diabetic Retinopathy in Retinal Fundus Photographs. *JAMA*, 316(22), 2402–2410.\n2. Tan, M. & Le, Q. (2019). EfficientNet: Rethinking Model Scaling for Convolutional Neural Networks. *ICML*.\n3. Graham, B. (2015). Kaggle Diabetic Retinopathy Detection competition report. (Illumination-correction preprocessing.)\n4. Selvaraju, R. R. et al. (2017). Grad-CAM: Visual Explanations from Deep Networks via Gradient-based Localization. *ICCV*.\n5. Zhang, H. et al. (2018). mixup: Beyond Empirical Risk Minimization. *ICLR*.\n6. Yun, S. et al. (2019). CutMix: Regularization Strategy to Train Strong Classifiers with Localizable Features. *ICCV*.\n7. Szegedy, C. et al. (2016). Rethinking the Inception Architecture for Computer Vision. *CVPR*. (Label smoothing.)\n8. Loshchilov, I. & Hutter, F. (2019). Decoupled Weight Decay Regularization. *ICLR*. (AdamW.)\n9. Cohen, J. (1968). Weighted kappa: Nominal scale agreement with provision for scaled disagreement. *Psychological Bulletin*, 70(4), 213–220.\n10. Wilkinson, C. P. et al. (2003). Proposed international clinical diabetic retinopathy and diabetic macular edema disease severity scales. *Ophthalmology*, 110(9), 1677–1682.\n11. Pizer, S. M. et al. (1987). Adaptive histogram equalization and its variations. *Computer Vision, Graphics, and Image Processing*, 39(3), 355–368. (CLAHE.)\n12. APTOS 2019 Blindness Detection. Kaggle / Aravind Eye Hospital. https://www.kaggle.com/c/aptos2019-blindness-detection","metadata":{}},{"cell_type":"code","source":"# ── Final session summary ───────────────────────────────────────────────────\nprint(\"=\" * 72)\nprint(\"  PROJECT COMPLETE — Diabetic Retinopathy Stage Detection\")\nprint(\"=\" * 72)\nprint(f\"\\nBest model : {best_tag} — {RESULTS[best_tag]['label']}\")\nfor k in [\"accuracy\", \"f1_macro\", \"qwk\", \"adjacent_accuracy\"]:\n    print(f\"  test {k:<20}{RESULTS[best_tag]['test'][k]:.4f}\")\nprint(f\"\\nReferable-DR screening: sensitivity {sens:.3f} · specificity {spec:.3f} · NPV {npv:.3f}\")\nprint(f\"\\nArtefacts written to {WORK}:\")\nfor item in [\"results.json\", \"results_summary.csv\", \"figures_for_report.zip\",\n             \"deploy/dr_model.pt\", \"deploy/gradio_app.py\"]:\n    p = WORK / item\n    mark = \"ok \" if p.exists() else \"-- \"\n    size = f\"{p.stat().st_size/1e6:.1f} MB\" if p.exists() else \"missing\"\n    print(f\"  [{mark}] {item:<32} {size}\")\nprint(f\"  [ok ] figures/  ({len(list(FIG_DIR.glob('*.png')))} PNG files at 150 dpi)\")\nprint(\"\\nNext steps:\")\nprint(\"  1. Download figures_for_report.zip and build the 20-page PDF report.\")\nprint(\"  2. Record the Gradio prototype in operation and host the video.\")\nprint(\"  3. Paste the hosted video URL into the report, as the brief requires.\")\nprint(\"=\" * 72)","metadata":{},"outputs":[],"execution_count":null}]}