Deep Learning for Medical Imaging

Contents

• Introduction and overview
• 0. Review: deep learning and statistics for medical images
• 1. Introduction to medical imaging
• 2. Classification
• 3. Detection
• 4. Segmentation
• 5. Registration
• 6. Denoising and reconstruction
• 7. Generative models: autoencoders and GANs
• 8. Validation, interpretation and reproducibility
• 9. Research directions and open problems
• Exercises
• Projects and assessment
• References

Introduction and overview

Audience: a student with a good M1 in mathematics, computer science or physics who has trained a convolutional network once, and wants to know what changes when the images are three-dimensional scans of patients and the predictions may reach a clinician.

Empirical first: a chest X-ray classifier reaches an excellent AUC in its own hospital and loses much of it in another one, because it had partly learned which hospital, scanner or portable device produced the image. A brain MRI classifier for Alzheimer’s disease reports near-perfect accuracy because slices of the same subject were put in both the training and the test set.

The course covers the main tasks of medical image analysis: classification, detection, segmentation, registration, denoising and reconstruction, and generation. For each task it presents the state-of-the-art deep learning methods next to the classical image-processing and machine-learning approaches they replaced or extended, and it ends with how to validate, interpret and reproduce such models. It continues Deep Learning and Medical Image Analysis Based on Generative, Geometric and Biophysical Models, and runs alongside Representation Learning for Computer Vision and Medical Imaging and Functional Brain Imaging and Brain–Computer Interfaces.

What you will be able to do

• Explain how CT, MRI, PET and ultrasound images are formed, and read a scan correctly: voxel spacing, world coordinates, intensity units and normalisation.
• Train a 2D or 3D convolutional classifier on medical images, handle class imbalance and prevalence, and evaluate it with sensitivity, specificity, ROC curves and AUC.
• Detect lesions and landmarks with anchor-based, focal-loss and heatmap detectors, and evaluate them with free-response ROC analysis.
• Segment organs and lesions with U-Net-type networks, choose between cross-entropy and Dice losses, and report Dice and Hausdorff distances.
• Formulate image registration as an optimisation problem, build diffeomorphic transformations by scaling and squaring, and train an unsupervised registration network.
• Reconstruct undersampled MRI and denoise low-dose images with learned denoisers, data-consistency layers and self-supervised training.
• Derive the training objectives of autoencoders, variational autoencoders and GANs, and use them for anomaly detection and image-to-image translation.
• Design a leakage-free validation protocol, attach confidence intervals to a performance estimate, and test the difference between two models.
• Compute and criticise saliency and Grad-CAM explanations, and write a reproducible, well-documented experiment.

Official objectives and syllabus (2026–2027)

The course description on master-mva.com starts from the observation that medical imaging gives unparalleled means to study the structure and function of the human body in vivo, but that interpretation is difficult: it must take into account three-dimensional, time-varying information coming from several types of images. Artificial intelligence holds great promise for assisting interpretation, and deep learning has led to impressive advances in medical image processing in recent years.

The course covers both theoretical and practical aspects of deep learning for medical imaging. For the main tasks of medical image analysis (classification, segmentation, registration, generative models…) it presents state-of-the-art deep learning techniques, alongside more traditional image-processing and machine-learning approaches, with examples from several applications (brain, cardiac…).

By the end of the course, the description expects students to know state-of-the-art deep learning techniques for medical imaging, to have a deeper understanding of deep learning methods that also applies to other types of data, and to know how to build and validate deep learning models for medical images. The course is organised as 9 sessions of 3 h, each either a 3 h lecture or a 1 h 30 lecture followed by a 1 h 30 lab, at CentraleSupélec. The topics of master-mva.com and of the Université Paris-Saclay programme are the same eight items; each is matched below.

Extensions beyond the syllabus: the prior-shift correction of probabilities (ch. 2), self-supervised denoising (ch. 6), cycle-consistent translation and its hallucinations (ch. 7), calibration and statistical comparison of models (ch. 8), and a self-study chapter on foundation models, federated learning and uncertainty (ch. 9).

Where this course fits in the programme

Official prerequisites (master-mva.com): a background in statistics, computer programming (Python) and some basic knowledge of deep learning. Chapter 0 reviews exactly what is used.

• Refreshers — Probability and Statistics (Bayes’ rule, likelihood, confidence intervals), Linear Algebra, Analysis and Optimization (SVD, gradients, Jacobians), Scientific Python, NumPy and PyTorch (every code block).
• First-semester courses that prepare it — Deep Learning, Introduction to Statistical Learning, Object Recognition and Computer Vision (detectors, segmentation networks), Medical Image Analysis Based on Generative, Geometric and Biophysical Models, Image Denoising: the Human–Machine Competition, Introduction to Probabilistic Graphical Models and Deep Generative Models (variational inference), Convex Optimization.
• Related second-semester courses — Representation Learning for Computer Vision and Medical Imaging, Deformable Models and Geodesic Methods for Image Analysis, Geometry and Shape Spaces, Inverse Problems and Imaging, Deep Learning for Image Restoration and Synthesis, Generative Models for Images, Functional Brain Imaging and Brain–Computer Interfaces, Biostatistics, Responsible Machine Learning, Multimodal Explainable AI.
• Throughout the year and afterwards — the Hospital Immersion pairs students with a physician, and the Research Internship is often a medical-imaging project in a hospital, laboratory or company.

Course outline

• Review: deep learning and statistics for medical images — losses as likelihoods, convolutional layers in 2D and 3D, receptive fields, generalisation and the test set.
• Introduction to medical imaging — how each modality forms its image, voxels and world coordinates, intensity units, preprocessing, the tasks and the public data.
• Classification — 2D and 3D classifiers, transfer learning, prevalence and class imbalance, ROC and AUC, multiple-instance learning, and the radiomics baseline.
• Detection — bounding boxes and IoU, region proposals and one-stage detectors, focal loss, non-maximum suppression, landmark heatmaps, FROC analysis.
• Segmentation — classical segmentation, U-Net and its 3D variants, Dice and cross-entropy losses, overlap and boundary metrics, nnU-Net.
• Registration — transformation models, similarity measures, diffeomorphisms, classical optimisation and learned registration with VoxelMorph.
• Denoising and reconstruction — CT and MRI forward models, filtered backprojection and compressed sensing, learned and self-supervised denoising, unrolled reconstruction networks.
• Generative models: autoencoders and GANs — PCA as the optimal linear autoencoder, the ELBO, anomaly detection, the GAN game and image-to-image translation.
• Validation, interpretation and reproducibility — leakage, confidence intervals, model comparison, calibration, saliency methods and their sanity checks, reproducible research.
• Research directions and open problems — foundation models, self-supervision, federated learning, uncertainty, diffusion models, fairness (self-study).

Session plan

27 h in class = 9 sessions of 3 h. The volume is the official one of the Université Paris-Saclay programme table (27 h of lectures), and master-mva.com announces the same organisation: 9 sessions of 3 h, each a lecture or a lecture followed by a lab. The assignment of chapters to sessions below is a proposal; it follows the order of the official topics.

Level 1 = foundations (reviews, definitions, first worked examples), 2 = core methods with their proofs, 3 = advanced material close to current research; the level never goes down from one session to the next.

Self-study, outside the 27 h (about 5 h in total): ch. 9 Research directions and open problems (≈ 3 h, with Ex. 15) · the code blocks of ch. 6–8, run and modified at home (≈ 2 h).

How to study this module

• Read ch. 0 and ch. 1 before the first session: the conventions on voxels, units and losses are used everywhere afterwards without warning.
• Do the exercises of the Practice column right after each session; the ★ ones check the session, the ★★ ones prepare the next.
• Run every code block, change one parameter (noise level, lesion contrast, number of iterations, split strategy) and predict the effect before looking.
• Read ch. 8 twice: once early, with ch. 2, to set up your project’s protocol; once at the end, when you interpret your results.
• Prepare the article presentation as you go: the papers cited in each chapter are good candidates, and the projects section lists more.
• Environment: Python 3 with NumPy, SciPy, scikit-learn and PyTorch on CPU; every code block uses synthetic data and runs in well under a minute on a laptop. For real data, the libraries nibabel, TorchIO and MONAI read and transform medical volumes.

0. Review: deep learning and statistics for medical images

Every chapter trains a network by minimising a loss, on images that are often volumes, and then measures its performance on a finite test set. This chapter reviews the three facts behind these steps: what a loss estimates, what a convolution costs in 3D, and how precise a test-set estimate is.

0.1 Losses are negative log-likelihoods

A classifier with parameters \theta outputs a probability vector q_\theta(x) over K classes. Training minimises the empirical risk over pairs (x_i, y_i).

Two consequences recur in this course. A well-trained classifier outputs probabilities, which can be corrected when the class frequencies change (ch. 2.3) and must be checked for calibration (ch. 8.3). The squared loss is likewise minimised by the conditional mean, the basis of denoising (ch. 6).

0.2 Convolutions in two and three dimensions

Medical images are often volumes: a brain MRI has about 256^3 voxels. A 3D convolutional layer slides a kernel along three axes instead of two.

The layer is equivariant to translations: shifting the input shifts the output. Stacking L layers of size 3 and stride 1 gives a receptive field of 1 + 2L voxels per axis; pooling and strides make it grow faster.

import torch
import torch.nn as nn
torch.set_num_threads(1)

conv2d = nn.Conv2d(32, 64, kernel_size=3, padding=1)    # 32 -> 64 channels, 3 x 3
conv3d = nn.Conv3d(32, 64, kernel_size=3, padding=1)    # same, 3 x 3 x 3
count = lambda m: sum(p.numel() for p in m.parameters())
print("parameters  2D:", count(conv2d), "  3D:", count(conv3d))

# the 2D layer on a 128 x 128 slice (batch of 1): shape is kept by padding=1
y = conv2d(torch.zeros(1, 32, 128, 128))
print("2D output", tuple(y.shape), f"{y.numel() * 4 / 2**20:.0f} MiB in float32")
# the same layer in 3D on a 128^3 volume: 128 times more voxels
n3 = 64 * 128 ** 3
print("3D output (1, 64, 128, 128, 128)", f"{n3 * 4 / 2**20:.0f} MiB in float32")

# receptive field of L stacked 3x3 convolutions with stride 1: 1 + 2L
for L in (5, 10, 20):
    print(f"{L} layers -> receptive field {1 + 2 * L} voxels per axis")

The 3D layer has three times more parameters (55 360 against 18 496), but its activations take 512 MiB against 4 MiB for one slice. Memory, not parameters, is why 3D networks are trained on patches and small batches (ch. 4.3).

0.3 Generalisation and the test set

Networks are trained with stochastic gradient methods (SGD with momentum, Adam), regularised by weight decay, data augmentation and early stopping. Hyperparameters are chosen on a validation set. The test set is used once, at the end, to estimate the risk of the final model.

With n = 100 test patients, a 95 % guarantee needs 2 e^{-200 t^2} = 0.05, that is t = \sqrt{\log(40)/200} \approx 0.136 .

An accuracy measured on 100 patients is therefore only known to within about 14 points by this bound, and about 8 points by the normal approximation of ch. 8.2 when it is near 80 %. Medical test sets are often that small.

1. Introduction to medical imaging

A medical image is a measurement, not a photograph: each modality measures a different physical quantity, on a 3D grid with its own spacing, units and noise. This chapter explains how the main modalities form their images, how a scan is stored and preprocessed, and which tasks the rest of the course addresses.

1.1 What makes medical images different

Four features separate medical images from the natural images of computer vision.

• Dimension: CT and MRI are volumes, sometimes with a time axis (cine cardiac MRI) or several contrasts per patient; whole-slide microscopy images have billions of pixels.
• Physics: the intensity has a physical meaning in CT (attenuation) and PET (tracer concentration), and an arbitrary scale in MRI.
• Data: data sets have hundreds or thousands of patients, not millions; labels come from experts, who disagree, and pathologies are rare.
• Stakes: predictions concern patients, so validation, interpretation and reproducibility (ch. 8) are part of the method.

1.2 How images are formed

X-ray radiography and CT. An X-ray beam crossing tissue of linear attenuation coefficient \mu loses intensity at a rate proportional to itself.

A radiograph shows \log(I_0 / I), the line integral of \mu along each ray, projected on a plane. Computed tomography (CT) measures these line integrals along many directions and inverts them to recover \mu on a 3D grid (ch. 6.2).

Magnetic resonance imaging (MRI). Hydrogen nuclei, aligned by a strong magnetic field, are excited by radio-frequency pulses and relax with tissue-dependent time constants T_1 and T_2. Magnetic field gradients encode position in the frequency and phase of the signal, so that the scanner samples the Fourier transform of the image: s(k) = \int \rho(x)\, e^{-2 i \pi \langle k, x \rangle}\, dx .

The space of k is called k-space. The weighting \rho depends on the sequence: T_1-weighted images show anatomy, T_2-weighted and FLAIR images show oedema and lesions, diffusion images show the motion of water. MRI intensities have no absolute unit and drift with a smooth multiplicative bias field.

Other modalities. Positron emission tomography (PET) counts photon pairs emitted by a radioactive tracer: its data are Poisson counts, and its images are noisy and of low resolution. Ultrasound images the echoes of sound waves, in real time, with speckle noise. Digital pathology scans stained tissue slides into RGB whole-slide images; ophthalmology and dermatology use photographs.

1.3 Volumes, voxels and world coordinates

A scan is stored as an array together with the geometry that places it in the patient. Formats such as DICOM (the clinical standard, one file per slice with a header) and NIfTI (common in research) both store this geometry.

The code below builds a sphere of radius 30 mm on such a grid, resamples it to 2 mm isotropic voxels, measures its volume and applies a CT soft-tissue window.

import numpy as np
from scipy.ndimage import zoom

# voxel -> world (mm) affine of a volume with spacing 0.8 x 0.8 x 5 mm (axial CT)
spacing = np.array([0.8, 0.8, 5.0])
A = np.diag(np.append(spacing, 1.0))
A[:3, 3] = [-200.0, -200.0, -100.0]                 # world position of voxel (0, 0, 0)
ijk = np.array([250, 250, 20, 1.0])
print("voxel (250, 250, 20) is at", A @ ijk, "mm")

# a synthetic volume: a sphere of radius 30 mm, 40 HU (soft tissue) in air (-1000 HU)
shape = (500, 500, 40)
grid = np.indices(shape).reshape(3, -1).T * spacing + A[:3, 3]
hu = np.where(np.linalg.norm(grid, axis=1) < 30, 40.0, -1000.0).reshape(shape)

# resample to isotropic 2 mm voxels (linear interpolation)
iso = zoom(hu, spacing / 2.0, order=1)
print("shape", hu.shape, "->", iso.shape)
vol_mm3 = (iso > -480).sum() * 2.0 ** 3
print(f"sphere volume: {vol_mm3 / 1e3:.1f} cm3 (exact {4 / 3 * np.pi * 30 ** 3 / 1e3:.1f} cm3)")

# a soft-tissue window [level 40, width 400] maps HU to [0, 1]
level, width = 40.0, 400.0
win = np.clip((iso - (level - width / 2)) / width, 0, 1)
print("window values: air", win.min(), " tissue", win.max())

The measured volume is 113.6 cm³ against 113.1 cm³ exactly: the error comes from the 5 mm slices, not from the method.

1.4 Preprocessing

Most pipelines apply the same steps before any learning.

1. Reorientation and resampling to a common orientation and spacing, with linear interpolation for images and nearest-neighbour interpolation for label maps.
2. Intensity normalisation: a fixed window in HU for CT; for MRI, a bias-field correction (the N4 algorithm of Tustison et al., 2010), then a z-score inside the body or brain mask, or a clipping between percentiles.
3. Cropping and masking, for instance skull stripping in brain MRI.
4. Spatial normalisation to a template by registration (ch. 5), when the analysis compares anatomies across subjects.

Data augmentation then generates plausible variants during training: flips (if the anatomy allows), small rotations and scalings, elastic deformations, intensity and contrast changes, simulated noise and bias fields. Libraries such as TorchIO and MONAI implement these operations for 3D volumes.

1.5 Tasks, data and challenges

The course follows the main tasks of medical image analysis.

Public data are numerous; the official references point to curated lists and to grand-challenge.org, which hosts many benchmark challenges. Examples used in this course: the Medical Segmentation Decathlon, BraTS (brain tumours), ACDC (cardiac MRI), LUNA16 (lung nodules), fastMRI (raw MRI data), OASIS and ADNI (ageing and Alzheimer’s disease).

Medical data are personal health data. They are de-identified before release, which includes removing names from DICOM headers and sometimes the face from head scans, and their use is governed by data-use agreements.

2. Classification

Classification assigns one label to an image or to a patient: diseased or healthy, a disease subtype, a risk group. This chapter adapts convolutional classifiers to medical data, then studies the two statistical issues that dominate in practice: the class frequencies and the choice of a metric.

2.1 Architectures for 2D and 3D images

For 2D images (radiographs, fundus photographs, dermoscopy, pathology patches), the standard choice is a 2D residual network or one of its successors, often pretrained on ImageNet. For volumes there are three options.

• 3D networks take the whole volume, or a large crop, as input. They see the 3D context but are limited by memory (ch. 0.2) and need more data.
• 2D networks on slices classify each slice, then aggregate the slice scores by averaging or by a small attention module. They reuse 2D pretrained weights.
• 2.5D networks feed a few neighbouring slices, or three orthogonal slices, as input channels.

Images over time. Cine cardiac MRI and follow-up scans add a time axis. The frames can be stacked as input channels, convolved over space and time together, or encoded one by one and then combined by a recurrent or attention module across time steps. Several contrasts of the same patient (T1, T2, FLAIR) are stacked as channels once they are registered (ch. 5).

Transfer learning. Initialising from ImageNet weights then fine-tuning is the default for 2D images. Raghu et al. (2019) found, on large retinal and chest X-ray tasks, that pretraining brought little accuracy over well-chosen small networks trained from scratch, but sped up convergence. Its benefit is largest with small data sets.

2.2 Sensitivity, specificity and the ROC curve

A binary classifier outputs a score S and decides “positive” when S > t. The clinical vocabulary names the error rates.

With ties, the AUC counts each tied pair as one half: this is the Mann–Whitney statistic divided by the number of pairs. The AUC does not depend on the threshold or on the prevalence, which makes it a good summary of ranking quality, but not of the decisions taken at one threshold.

2.3 Prevalence and class imbalance

Diseases are rare in screening populations and over-represented in research data sets. The predictive value depends on the prevalence \pi of the disease, by Bayes’ rule.

A classifier trained with cross-entropy estimates the posterior of its training distribution (Proposition 0.2). When the class frequencies differ at deployment while the images of each class do not change (a prior shift), the probabilities can be corrected.

A network trained on balanced data (\pi_0 = 1/2) that outputs q = 0.9 means, at \pi = 1\%, a probability of \frac{0.01 \times 0.9}{0.01 \times 0.9 + 0.99 \times 0.1} \approx 0.083, the number of Example 2.3. Weighting the loss of class k by w_k is equivalent to training under a prevalence proportional to w_k \pi_k, so the same correction undoes class weights. Weighted losses and oversampling of the rare class change the effective prior. The focal loss of ch. 3.3 is not a prior change, but it also distorts the outputs, which must then be recalibrated (ch. 8.3) before being read as probabilities.

2.4 Weak labels and multiple-instance learning

A whole-slide pathology image has about 10^{10} pixels and a single diagnosis. It is cut into N tiles, and only the label of the whole bag of tiles is known.

Ilse, Tomczak and Welling (2018) pool with attention weights computed by a small network with parameters w and V: z = \sum_{n=1}^N a_n h_n, \qquad a_n = \frac{\exp\big(w^\top \tanh(V h_n)\big)}{\sum_{m} \exp\big(w^\top \tanh(V h_m)\big)} .

The pooled vector z is invariant to permutations of the tiles, and the weights a_n show which tiles drove the decision. Campanella et al. (2019) trained weakly supervised multiple-instance models (max-pooling over tiles, then a recurrent aggregator) on more than 40 000 slides, with diagnoses as the only labels.

2.5 The radiomics baseline and a worked example

Before deep learning, and still today in many clinical studies, classification used radiomics (Lambin et al., 2012): hand-crafted features of a segmented region (volume, shape, intensity statistics, co-occurrence texture), followed by a penalised logistic regression or a random forest. On small data sets this baseline is hard to beat, and every deep model should be compared with it.

The code trains a small CNN to detect a faint blob (a “lesion”) in a textured background, at a prevalence of 30 %, with a weighted loss. It compares the CNN with a hand-crafted feature, the maximum of the smoothed image, and checks Theorem 2.2.

import numpy as np
import torch
import torch.nn as nn
from scipy.ndimage import gaussian_filter
from sklearn.metrics import roc_auc_score

rng = np.random.default_rng(0)
torch.manual_seed(0)
torch.set_num_threads(1)

def make(n, prevalence=0.3, size=24):
    """Smooth random 'tissue' texture; positives carry a faint blob."""
    x = gaussian_filter(rng.standard_normal((n, size, size)), sigma=(0, 1.5, 1.5)) * 3
    y = (rng.random(n) < prevalence).astype(np.float32)
    ii, jj = np.mgrid[:size, :size]
    for k in np.flatnonzero(y):
        ci, cj = rng.integers(5, size - 5, 2)
        x[k] += 1.5 * np.exp(-((ii - ci) ** 2 + (jj - cj) ** 2) / (2 * 2.0 ** 2))
    x += 0.3 * rng.standard_normal(x.shape)
    return torch.tensor(x[:, None], dtype=torch.float32), torch.tensor(y)

xtr, ytr = make(1200)
xte, yte = make(600)

net = nn.Sequential(
    nn.Conv2d(1, 8, 3, padding=1), nn.ReLU(), nn.MaxPool2d(2),
    nn.Conv2d(8, 16, 3, padding=1), nn.ReLU(), nn.MaxPool2d(2),
    nn.Conv2d(16, 16, 3, padding=1), nn.ReLU(),
    nn.AdaptiveMaxPool2d(1), nn.Flatten(), nn.Linear(16, 1))
pos_weight = (1 - ytr.mean()) / ytr.mean()           # re-weights the rare class
loss_fn = nn.BCEWithLogitsLoss(pos_weight=pos_weight)
opt = torch.optim.Adam(net.parameters(), lr=3e-3)
for epoch in range(6):
    perm = torch.randperm(len(xtr))
    for b in range(0, len(xtr), 64):
        idx = perm[b:b + 64]
        opt.zero_grad()
        loss = loss_fn(net(xtr[idx]).squeeze(1), ytr[idx])
        loss.backward()
        opt.step()

with torch.no_grad():
    s = net(xte).squeeze(1).numpy()
y = yte.numpy()
# AUC as the Mann-Whitney probability P(score of a positive > score of a negative)
pos, neg = s[y == 1], s[y == 0]
mw = (pos[:, None] > neg[None, :]).mean() + 0.5 * (pos[:, None] == neg[None, :]).mean()
print(f"CNN test AUC: {roc_auc_score(y, s):.3f} (Mann-Whitney: {mw:.3f})")
# baseline 'radiomic' feature: maximum of the image after a smoothing filter
feat = np.array([gaussian_filter(im[0], 2).max() for im in xte.numpy()])
print(f"hand-crafted feature AUC: {roc_auc_score(y, feat):.3f}")

The CNN reaches an AUC of about 0.76 and the hand-crafted feature about 0.60; the Mann–Whitney count agrees with the library to three decimals. The blob is barely visible, and the CNN has learned a better matched filter than the smoothing kernel chosen by hand.

3. Detection

Detection answers “where are the lesions, and how many?”: lung nodules in CT, polyps in colonoscopy, microbleeds in brain MRI, cells in microscopy, anatomical landmarks. The output is a set of boxes or points with confidence scores, and the evaluation must count both missed lesions and false alarms.

3.1 Boxes and overlap

A 2D box is b = (x_0, y_0, x_1, y_1); a 3D box adds a third axis. A predicted box is matched to a ground-truth box when they overlap enough.

IoU is very sensitive to localisation for small objects: shifting a 4 \times 4 box by one pixel along one axis gives an IoU of 12 / 20 = 0.6.

3.2 Region proposals and one-stage detectors

Modern detectors place many reference boxes, called anchors, at every position of a feature map, at several scales and aspect ratios. For each anchor, the network predicts an objectness score and a correction of the box.

• Two-stage detectors (Faster R-CNN, Ren et al., 2015): a region proposal network scores the anchors; the best few hundred proposals are cropped from the feature map and passed to a second head that classifies them and refines the boxes.
• One-stage detectors (YOLO, Redmon et al., 2016; RetinaNet, Lin et al., 2017) classify and regress all anchors (all grid cells, for the first YOLO) in one pass. They are faster, but face an extreme imbalance: tens of thousands of background anchors per lesion.
• Medical variants add 3D anchors, and Retina U-Net (Jaeger et al., 2020) adds a segmentation decoder whose voxel-level supervision helps the detector when lesions are few.

3.3 The focal loss

The background anchors are mostly easy: the network quickly gives them a probability close to 1 of being background. Their sum still dominates the cross-entropy and its gradient.

With \gamma = 2, an easy background anchor at p_t = 0.99 is down-weighted by 10^{-4}, a hard one at p_t = 0.5 by only 0.25. A thousand easy anchors then weigh less than one hard anchor.

3.4 Non-maximum suppression

A detector fires on many overlapping anchors around each lesion. Non-maximum suppression keeps one box per object.

3.5 Candidates, false positive reduction and landmarks

The classical detection pipeline has two stages. A sensitive candidate generator proposes locations, for instance with Hessian-based blob filters or intensity thresholds. A classifier then rejects the false candidates. Setio et al. (2016) kept this structure for lung nodules in CT and replaced the second stage by a multi-view CNN, which classifies 2D patches cut through each candidate along several planes.

Anatomical landmarks (points on the bones of a hand radiograph, cephalometric points of the skull, vertebra centres in spine CT) are single points. They are detected by heatmap regression: the network predicts, for landmark \ell, a map close to the Gaussian G_\ell(x) = \exp\big(-\|x - \ell\|^2 / (2\sigma^2)\big), and the prediction is the location of its maximum (Payer et al., 2019).

The soft-argmax lets the network be trained directly on the coordinate error, in millimetres.

3.6 Evaluation: FROC analysis

A scan contains a variable number of lesions, and a detector can raise any number of false alarms; specificity is undefined. Detection is evaluated per lesion and per scan.

The code computes IoU and non-maximum suppression on four boxes, then an FROC curve on simulated detections: 90 % of the nodules receive a candidate, and each scan has about 20 false candidates.

import numpy as np

def iou(a, b):
    """Boxes as (x0, y0, x1, y1); a: (n, 4), b: (m, 4) -> (n, m) IoU matrix."""
    x0 = np.maximum(a[:, None, 0], b[None, :, 0]); y0 = np.maximum(a[:, None, 1], b[None, :, 1])
    x1 = np.minimum(a[:, None, 2], b[None, :, 2]); y1 = np.minimum(a[:, None, 3], b[None, :, 3])
    inter = np.clip(x1 - x0, 0, None) * np.clip(y1 - y0, 0, None)
    area = lambda r: (r[:, 2] - r[:, 0]) * (r[:, 3] - r[:, 1])
    return inter / (area(a)[:, None] + area(b)[None, :] - inter)

def nms(boxes, scores, thr=0.5):
    order, keep = np.argsort(-scores), []
    while order.size:
        i = order[0]
        keep.append(i)
        order = order[1:][iou(boxes[[i]], boxes[order[1:]])[0] < thr]
    return np.array(keep)

boxes = np.array([[10, 10, 30, 30], [12, 11, 31, 32], [50, 50, 70, 72], [11, 9, 29, 31.0]])
scores = np.array([0.9, 0.8, 0.7, 0.6])
print("IoU(box0, box1) =", round(iou(boxes[[0]], boxes[[1]])[0, 0], 3))
print("kept after NMS:", nms(boxes, scores))

# FROC on synthetic scans: each scan has a few nodules; the detector outputs scored candidates
rng = np.random.default_rng(1)
n_scans, hits, all_scores, n_nodules = 100, [], [], 0
for s in range(n_scans):
    k = rng.integers(1, 4); n_nodules += k
    tp = rng.random(k) < 0.9                          # 90 % of nodules get a candidate
    all_scores += list(rng.normal(1.5, 1.0, tp.sum())); hits += [1] * tp.sum()
    nf = rng.poisson(20)                              # false candidates in this scan
    all_scores += list(rng.normal(0.0, 1.0, nf)); hits += [0] * nf
order = np.argsort(-np.array(all_scores)); hits = np.array(hits)[order]
sens = np.cumsum(hits) / n_nodules                    # sensitivity at each threshold
fp_per_scan = np.cumsum(1 - hits) / n_scans           # average false positives per scan
cpm = [sens[fp_per_scan <= r].max() for r in (1/8, 1/4, 1/2, 1, 2, 4, 8)]
print("sensitivity at 1/8 ... 8 FP/scan:", np.round(cpm, 3))
print(f"CPM (mean of the seven) = {np.mean(cpm):.3f}")

The IoU of the first two boxes is 342 / 457 \approx 0.748, and NMS keeps boxes 0 and 2. The FROC sensitivity rises from 0.13 to 0.79 between 1/8 and 8 false positives per scan, for a CPM of 0.43: with 20 false candidates per scan, false positive reduction is the real problem.

4. Segmentation

Segmentation labels every voxel: which voxels belong to the left ventricle, to the hippocampus, to the tumour? It gives volumes, shapes and boundaries, from which clinicians compute measurements such as the ejection fraction of the heart or the volume of a tumour. It is the task where deep learning has changed practice the most.

4.1 Classical segmentation

• Thresholding labels the voxels whose intensity lies in an interval; Otsu’s method chooses the threshold that maximises the between-class variance of the histogram. It works for bone in CT, rarely in MRI.
• Region growing starts from a seed and adds neighbouring voxels of similar intensity; it leaks through weak boundaries.
• Energy-based methods (active contours, level sets, graph cuts) minimise an energy that balances fidelity to the image and smoothness of the boundary. Graph cuts find the global minimum of such binary energies when the pairwise terms are submodular, for instance when they penalise label changes between neighbours.
• Atlas-based segmentation registers one or several labelled images (atlases) to the new image and transfers their labels; multi-atlas methods fuse the transferred labels by voting. It uses the registration of ch. 5 and was the state of the art for brain structures before deep learning.

4.2 Fully convolutional networks and U-Net

A fully convolutional network has no dense layer, so it maps an image of any size to a map of class scores of the same size. U-Net (Ronneberger, Fischer and Brox, 2015) became the reference architecture for biomedical segmentation.

The encoder provides context: after four downsamplings, a voxel’s features depend on a large region. The skip connections bring back the fine localisation that downsampling lost. Without them, boundaries are blurred.

4.3 Three-dimensional networks and patch-based training

The 3D U-Net (Çiçek et al., 2016) replaces every 2D operation by its 3D counterpart, and V-Net (Milletari, Navab and Ahmadi, 2016) adds residual blocks and the Dice loss below. A full volume rarely fits in GPU memory (ch. 0.2), so training and inference work on patches.

• Training: sample patches of, say, 128^3 voxels, with a fixed fraction centred on foreground voxels, so that rare structures are seen.
• Inference: slide the patch over the volume with overlap, and average the overlapping predictions, often with weights that decrease towards the patch border, where context is missing.

4.4 Losses and metrics

For a binary segmentation, let P be the predicted set of voxels and G the ground truth.

Because the true negatives are absent, the Dice loss is not swamped by the background, which can fill 99.9 % of a brain MRI with a small lesion. The cross-entropy, a mean over all voxels, is dominated by that background. It gives smoother gradients, however, and the sum of the two losses is a robust default.

The overlap measures miss the boundary: a thin false-positive island far away barely changes the Dice coefficient but can matter clinically.

4.5 nnU-Net and applications

Isensee et al. (2021) showed with nnU-Net that a carefully configured U-Net outperformed most specialised architectures on 23 public data sets. nnU-Net reads a fingerprint of the data set (spacings, image sizes, intensity statistics) and sets the preprocessing, patch size, network depth and batch size by fixed rules. It then trains 2D, 3D and cascaded variants by five-fold cross-validation and selects or ensembles them empirically.

• Brain tumours: the BraTS challenge (Menze et al., 2015) segments, from four MRI contrasts, the enhancing tumour, the tumour core and the whole tumour, three nested regions.
• Cardiac MRI: the ACDC challenge (Bernard et al., 2018) segments the left and right ventricular cavities and the myocardium on cine MRI. From the left-ventricle volumes at end-diastole (EDV) and end-systole (ESV) follows the ejection fraction (\mathrm{EDV} - \mathrm{ESV}) / \mathrm{EDV}: with 150 mL and 60 mL, it is 60 %.

The code trains a two-level U-Net on noisy ellipses with an intensity shading, with the sum of cross-entropy and soft Dice, and compares it with Gaussian smoothing followed by a threshold.

import numpy as np
import torch
import torch.nn as nn
import torch.nn.functional as F
from scipy.ndimage import gaussian_filter

torch.manual_seed(0)
torch.set_num_threads(1)
rng = np.random.default_rng(0)

def make(n, size=32):
    """Noisy images of one ellipse (the 'organ'), with its binary mask."""
    ii, jj = np.mgrid[:size, :size]
    x, m = np.zeros((n, 1, size, size)), np.zeros((n, 1, size, size))
    for k in range(n):
        ci, cj = rng.uniform(10, 22, 2); a, b = rng.uniform(3, 8, 2)
        m[k, 0] = ((ii - ci) / a) ** 2 + ((jj - cj) / b) ** 2 < 1
        x[k, 0] = 0.6 * m[k, 0] + 0.2 * np.sin(jj / 3.0 + rng.uniform(0, 6))  # contrast + shading
    x += 0.35 * rng.standard_normal(x.shape)
    return torch.tensor(x, dtype=torch.float32), torch.tensor(m, dtype=torch.float32)

def block(cin, cout):
    return nn.Sequential(nn.Conv2d(cin, cout, 3, padding=1), nn.ReLU(),
                         nn.Conv2d(cout, cout, 3, padding=1), nn.ReLU())

class TinyUNet(nn.Module):
    def __init__(self):
        super().__init__()
        self.enc1, self.enc2 = block(1, 8), block(8, 16)
        self.up = nn.ConvTranspose2d(16, 8, 2, stride=2)
        self.dec1 = block(16, 8)                       # 8 (skip) + 8 (upsampled) channels
        self.out = nn.Conv2d(8, 1, 1)
    def forward(self, x):
        e1 = self.enc1(x)
        e2 = self.enc2(F.max_pool2d(e1, 2))
        d1 = self.dec1(torch.cat([e1, self.up(e2)], dim=1))   # skip connection
        return self.out(d1)

def soft_dice_loss(logits, m, eps=1.0):
    p = torch.sigmoid(logits)
    inter = (p * m).sum(dim=(1, 2, 3))
    return 1 - ((2 * inter + eps) / (p.sum(dim=(1, 2, 3)) + m.sum(dim=(1, 2, 3)) + eps)).mean()

xtr, mtr = make(400)
xte, mte = make(100)
net = TinyUNet()
opt = torch.optim.Adam(net.parameters(), lr=3e-3)
for it in range(200):
    idx = torch.randint(0, len(xtr), (16,))
    logits = net(xtr[idx])
    loss = F.binary_cross_entropy_with_logits(logits, mtr[idx]) + soft_dice_loss(logits, mtr[idx])
    opt.zero_grad(); loss.backward(); opt.step()

with torch.no_grad():
    pred = (net(xte) > 0).float()
dice = 2 * (pred * mte).sum((1, 2, 3)) / (pred.sum((1, 2, 3)) + mte.sum((1, 2, 3)))
# classical baseline: Gaussian smoothing, then a global threshold
thr = torch.tensor(gaussian_filter(xte.numpy(), sigma=(0, 0, 1.5, 1.5)) > 0.3).float()
dice_thr = 2 * (thr * mte).sum((1, 2, 3)) / (thr.sum((1, 2, 3)) + mte.sum((1, 2, 3)))
print(f"U-Net Dice: mean {dice.mean():.3f}, worst case {dice.min():.3f}")
print(f"smoothing + threshold Dice: mean {dice_thr.mean():.3f}")

After 200 iterations the U-Net reaches a mean Dice of about 0.89, against 0.79 for the smoothed threshold, whose parameters were even tuned on the test images. The worst case, near 0.58, is a small ellipse, as Example 4.5 predicts: report the distribution of per-case scores, not only the mean.

5. Registration

Registration finds the spatial transformation that aligns two images: two scans of the same patient taken months apart, an MRI and a CT, or a patient and an atlas. It is needed to compare images voxel by voxel, to transfer labels and to measure change.

5.1 The registration problem

The warped image M \circ \varphi is evaluated at non-integer positions, by linear or B-spline interpolation. Writing \varphi(x) = x + u(x) with a displacement field u, a typical regulariser is the smoothness term \int \|\nabla u\|^2.

5.2 Transformation models

• Rigid: a rotation and a translation, 6 parameters in 3D; enough between two scans of the same head.
• Affine: a linear map and a translation, 12 parameters; it absorbs global differences of size and shear between subjects.
• Free-form deformations (Rueckert et al., 1999): a displacement interpolated by cubic B-splines from a coarse grid of control points.
• Dense displacements: one vector per voxel, regularised by \mathcal{R}.

A transformation is physically plausible if it does not tear or fold tissue: it should be invertible, with positive local volume change.

Diffeomorphisms (smooth bijections with smooth inverse) are obtained by integrating a velocity field. With a stationary velocity field v, the transformation is the time-one flow of the equation \frac{d}{dt} \phi_t = v \circ \phi_t, with \phi_0 = \mathrm{id}, written \exp(v) (Arsigny et al., 2006; Ashburner, 2007).

In practice \phi_{2^{-N}} is replaced by \mathrm{id} + 2^{-N} v, then squared N times by composing it with itself (about 7 compositions). Each composition is one interpolation, so the whole computation is differentiable.

5.3 Similarity measures

• Sum of squared differences, \int (F - M \circ \varphi)^2, for images of the same modality with the same intensities.
• Normalised cross-correlation, often computed in local windows, for the same modality with different contrast or bias field.
• Mutual information (Viola and Wells, 1997; Maes et al., 1997), for different modalities such as MRI and CT: \mathrm{MI}(F, M) = \sum_{a, b} p(a, b) \log \frac{p(a, b)}{p(a)\, p(b)}, where p(a, b) is the joint histogram of the intensities of F and M \circ \varphi at the same voxels.

Mutual information therefore only asks that intensities correspond consistently, not that they be equal: bone is bright in CT and dark in a T_1-weighted MRI.

5.4 Classical optimisation

Classical methods minimise E for each pair of images, by gradient descent over the parameters, from coarse to fine resolution to avoid local minima. Demons algorithms (Thirion, 1998; Vercauteren et al., 2009) alternate a force step and a Gaussian smoothing of the field; SyN (Avants et al., 2008) optimises symmetric diffeomorphisms. These methods are accurate but take minutes to hours per pair of 3D images.

5.5 Learning registration

A network g_\theta takes the pair (F, M) and outputs a displacement field u = g_\theta(F, M). The warped image is computed by a spatial transformer (Jaderberg et al., 2015), a differentiable interpolation layer.

VoxelMorph (Balakrishnan et al., 2019) trains g_\theta without ground-truth deformations, on pairs of images, by minimising the energy of Definition 5.1: \mathcal{L}(\theta) = \mathbb{E}_{(F, M)}\Big[\mathcal{D}\big(F, M \circ (\mathrm{id} + g_\theta(F, M))\big) + \lambda \|\nabla g_\theta(F, M)\|^2\Big].

The optimisation is amortised: once trained, the network registers a new pair in one forward pass, in about a second on a GPU. When segmentations are available for training, a Dice term between warped and fixed labels improves anatomical accuracy. Dalca et al. (2019) make the output a stationary velocity field integrated by scaling and squaring layers (Theorem 5.3), which gives diffeomorphisms.

Supervised alternatives regress known deformations, either synthetic ones applied to real images or ones computed by a classical method; they inherit the errors of the method they imitate.

Evaluation. A good image match is easy and meaningless: an unregularised field can make any image look like any other. Registration is evaluated on what it should preserve: the Dice overlap of anatomical labels transferred by the transformation, the target registration error between corresponding landmarks in millimetres, and the fraction of voxels with J_\varphi \leq 0.

The code registers an image to an affinely transformed copy of itself by gradient descent through grid_sample, then compares a direct displacement with its scaling-and-squaring version.

import numpy as np
import torch
import torch.nn.functional as F
from scipy.ndimage import gaussian_filter

torch.manual_seed(0)
torch.set_num_threads(1)
rng = np.random.default_rng(0)
n = 64
img = gaussian_filter((rng.random((n, n)) > 0.5).astype(float), 3)      # a blobby 'anatomy'
fixed = torch.tensor(img, dtype=torch.float32)[None, None]

def warp(im, theta):                     # theta: 2 x 3 affine in normalised [-1, 1] coordinates
    grid = F.affine_grid(theta[None], im.shape, align_corners=False)
    return F.grid_sample(im, grid, align_corners=False, padding_mode="border")

a = np.deg2rad(10)
true = torch.tensor([[1.05 * np.cos(a), -np.sin(a), 0.10], [np.sin(a), 0.95 * np.cos(a), -0.05]],
                    dtype=torch.float32)
moving = warp(fixed, true)               # the moving image is a warped copy of the fixed one

# 1) affine registration: find theta such that moving o T_theta ~ fixed  (SSD + gradient descent)
theta = torch.tensor([[1.0, 0, 0], [0, 1.0, 0]], requires_grad=True)
opt = torch.optim.Adam([theta], lr=0.01)
for it in range(300):
    loss = ((warp(moving, theta) - fixed) ** 2).mean()
    opt.zero_grad(); loss.backward(); opt.step()
A = torch.cat([theta.detach(), torch.tensor([[0, 0, 1.0]])])
B = torch.cat([true, torch.tensor([[0, 0, 1.0]])])
print(f"final SSD {loss.item():.2e}; composed transform (should be ~identity):")
print(np.round((B @ A).numpy()[:2], 3))

# 2) folding: Jacobian determinant of x -> x + u(x) for a large random displacement u
def jac_det(u):                           # u: (2, n, n) displacement in pixels, u[a] along axis a
    d = lambda f, axis: np.gradient(f, axis=axis)
    return (1 + d(u[0], 0)) * (1 + d(u[1], 1)) - d(u[0], 1) * d(u[1], 0)

v = np.stack([gaussian_filter(rng.standard_normal((n, n)), 4) for _ in range(2)])
v *= 8 / np.abs(v).max()                  # smooth field, up to 8 pixels
print(f"u = v directly:      {100 * (jac_det(v) <= 0).mean():.1f} % folded pixels")

# 3) scaling and squaring: phi = exp(v) = (id + v / 2^N) composed 2^N times
def compose(u, w):                        # displacement of (id + u) o (id + w)
    ii, jj = np.mgrid[:n, :n].astype(float)
    coords = torch.tensor(np.stack([jj + w[1], ii + w[0]]) / (n - 1) * 2 - 1, dtype=torch.float32)
    grid = coords.permute(1, 2, 0)[None]
    u_at = F.grid_sample(torch.tensor(u, dtype=torch.float32)[None], grid,
                         align_corners=True, padding_mode="border")[0].numpy()
    return w + u_at

N = 7
u = v / 2 ** N
for _ in range(N):
    u = compose(u, u)
print(f"u = exp(v) - id:     {100 * (jac_det(u) <= 0).mean():.1f} % folded pixels,"
      f" min det {jac_det(u).min():.2f}")

The recovered transform composed with the true one is the identity to within 0.01. Used directly as a displacement, the smooth field folds 3.7 % of the pixels; its exponential folds none, with a minimal Jacobian determinant of 0.15.

6. Denoising and reconstruction

Scanners do not measure images: CT measures projections and MRI measures Fourier coefficients, both corrupted by noise. Reducing the radiation dose or the acquisition time degrades these measurements further. This chapter writes the forward models, recalls the classical reconstructions, and shows how networks denoise and reconstruct while staying consistent with the data.

6.1 Forward models

All reconstruction problems of this chapter have the form y = \mathcal{A} x + n, where x is the unknown image, \mathcal{A} the known acquisition operator and n the noise.

In MRI, \mathcal{A} = M \mathcal{F}: the Fourier transform followed by a sampling mask M that keeps only some k-space locations. Acquiring one line in R shortens the scan by a factor R, the acceleration.

6.2 Classical reconstruction

Writing the inverse 2D Fourier transform in polar coordinates gives filtered backprojection: f(x) = \int_0^\pi \int_{\mathbb{R}} \hat p_\alpha(\omega)\, |\omega|\, e^{2i\pi \omega \langle x, \theta \rangle}\, d\omega\, d\alpha .

Each projection is filtered by the ramp |\omega|, which amplifies high-frequency noise, then smeared back along its rays. With few angles or a low dose, the image shows streaks and noise.

For undersampled MRI, the zero-filled reconstruction \mathcal{F}^{-1}(y) shows aliasing. Compressed sensing (Lustig, Donoho and Pauly, 2007) samples k-space randomly, so that aliasing looks like noise, and solves \hat x = \arg\min_x \tfrac12 \|M \mathcal{F} x - y\|^2 + \lambda \|\Psi x\|_1, with \Psi a wavelet or finite-difference transform in which images are sparse.

6.3 Noise models

• CT: photon counts are Poisson; at low dose, the noise of the log-transformed data grows where few photons arrive, and it becomes correlated after reconstruction.
• PET: Poisson counts, at much lower numbers.
• MRI: the complex k-space data have white complex Gaussian noise; the displayed magnitude image does not.

The background of a magnitude image is therefore not zero, and low-signal regions are biased upwards; a denoiser trained with Gaussian noise on magnitude images learns the wrong model (Gudbjartsson and Patz, 1995).

6.4 Learned denoising

A network trained on pairs (noisy, clean) with the squared loss estimates the conditional mean of the clean image (Proposition 0.2 and its analogue for the squared loss). For low-dose CT, Chen et al. (2017) trained a residual encoder–decoder on pairs of simulated low-dose and normal-dose images; Jin et al. (2017) post-processed sparse-view filtered backprojections with a U-Net.

Clean targets are rare in medicine: no patient is scanned twice at full dose to make a training set. Noise2Noise removes the need for them.

Two independent acquisitions, or two halves of the raw data, are enough. Noise2Self (Batson and Royer, 2019) goes further: it trains on single noisy images by predicting each pixel from the others, which is valid when the noise is independent between pixels.

6.5 Unrolled reconstruction networks

The most reliable learned reconstructions alternate a network and the physics. Starting from the zero-filled image, each iteration applies a learned denoiser D_\theta, then a data-consistency step that brings the image back towards the measurements.

As \lambda \to \infty, the measured coefficients are simply put back. Cascades of such blocks were trained end to end by Schlemper et al. (2018); the variational network of Hammernik et al. (2018) and MoDL (Aggarwal, Mani and Jacob, 2019) unroll gradient or conjugate-gradient steps in the same way. The fastMRI data set (Zbontar et al., 2018) provides raw k-space for training and benchmarking.

Learned reconstructions can be unstable to tiny perturbations of the data, and can hallucinate or erase small structures (Antun et al., 2020). Radiologists’ reading of the reconstructed images, not only PSNR, is the final test.

The code simulates a fourfold-accelerated Cartesian acquisition of a head phantom, then alternates a median-filter denoiser and the exact data-consistency step, a hand-made version of an unrolled network. It also checks Proposition 6.3.

import numpy as np
from scipy.ndimage import median_filter

rng = np.random.default_rng(0)
n = 128
yy, xx = (np.mgrid[:n, :n] - n / 2) / (n / 2)
def ellipse(cx, cy, a, b):
    return ((xx - cx) / a) ** 2 + ((yy - cy) / b) ** 2 < 1
# a simple head phantom: skull, brain, ventricles, a small lesion
x = 1.0 * ellipse(0, 0, 0.72, 0.92) - 0.7 * ellipse(0, 0, 0.66, 0.86)
x += 0.2 * ellipse(-0.18, -0.1, 0.1, 0.3) + 0.2 * ellipse(0.18, -0.1, 0.1, 0.3)
x += 0.3 * ellipse(0.3, 0.45, 0.06, 0.06)

F = lambda u: np.fft.fft2(u, norm="ortho")              # unitary DFT: k-space
Fi = lambda k: np.fft.ifft2(k, norm="ortho")
# Cartesian undersampling: keep the 9 central phase-encoding lines + about 20 % of the others
rows = rng.random(n) < 0.2
rows[np.abs(np.fft.fftfreq(n) * n) < 5] = True
M = np.repeat(rows[:, None], n, axis=1)
sigma = 0.01
y = M * (F(x) + sigma * (rng.standard_normal((n, n)) + 1j * rng.standard_normal((n, n))))
psnr = lambda u: 10 * np.log10(x.max() ** 2 / np.mean((np.abs(u) - x) ** 2))
print(f"acceleration R = {n / rows.sum():.1f}")

z = Fi(y)                                                 # zero-filled reconstruction
print(f"zero-filled PSNR: {psnr(z):.2f} dB")

# alternate a denoiser D and an exact data-consistency step (Proposition 6.5, lambda -> infinity)
u = z
for k in range(50):
    d = median_filter(u.real, size=3) + 1j * median_filter(u.imag, size=3)   # D
    u = Fi(M * y + (1 - M) * F(d))                                           # DC
print(f"denoise + data consistency, 50 iterations: {psnr(u):.2f} dB")

# Rician bias: magnitude of pure complex Gaussian noise (zero signal)
noise = sigma * (rng.standard_normal(10**6) + 1j * rng.standard_normal(10**6))
print(f"mean magnitude of noise: {np.abs(noise).mean():.5f}, "
      f"sigma*sqrt(pi/2) = {sigma * np.sqrt(np.pi / 2):.5f}")

The zero-filled image reaches 16.6 dB; fifty iterations of median filtering and data consistency bring it to 21.7 dB, with no training at all. A learned denoiser in place of the median filter, trained end to end through the iterations, is exactly a cascade network.

7. Generative models: autoencoders and GANs

Generative models learn the distribution of images. In medical imaging they serve to detect anomalies as departures from normal anatomy, to translate one modality into another (MRI to CT), and to synthesise training data. This chapter derives the objectives of autoencoders, variational autoencoders and GANs, and shows what can go wrong.

7.1 Autoencoders and PCA

Non-linear convolutional autoencoders generalise PCA to curved families of images. Their reconstructions are blurred, because the squared loss averages what the code does not determine (Proposition 0.2).

7.2 Variational autoencoders

A variational autoencoder (Kingma and Welling, 2014) is a latent-variable model p_\theta(x) = \int p_\theta(x \mid z)\, p(z)\, dz, with p(z) = \mathcal{N}(0, I), trained together with an encoder q_\phi(z \mid x) that approximates the posterior.

With a Gaussian encoder \mathcal{N}(\mu, \mathrm{diag}(s^2)), the KL term is \frac12 \sum_j (\mu_j^2 + s_j^2 - \log s_j^2 - 1). The expectation is estimated by sampling z = \mu + s \odot \varepsilon with \varepsilon \sim \mathcal{N}(0, I), the reparameterisation that lets gradients flow to \mu and s.

7.3 Anomaly detection

Annotating every possible pathology is impossible, but healthy scans are plentiful. A generative model trained on healthy images only reconstructs healthy anatomy well; a lesion is badly reconstructed, so the residual |x - D(E(x))| highlights it. Schlegl et al. (2017) did this with a GAN on retinal OCT images, and Baur et al. (2021) compared many autoencoder variants on brain MRI.

The code uses the optimal linear autoencoder of Theorem 7.2, fitted on 2000 synthetic healthy “brain slices”, to detect faint lesions in 100 of 200 test slices.

import numpy as np
from scipy.ndimage import gaussian_filter
from sklearn.decomposition import PCA
from sklearn.metrics import roc_auc_score

rng = np.random.default_rng(0)
n = 32
yy, xx = (np.mgrid[:n, :n] - n / 2) / (n / 2)

def healthy(m):
    """'Brain slices': an ellipse whose size, ventricle width and contrast vary between subjects."""
    out = np.empty((m, n, n))
    for k in range(m):
        a, b, w, c = rng.uniform(0.7, 0.9), rng.uniform(0.8, 0.95), rng.uniform(0.1, 0.25), rng.uniform(0.8, 1.2)
        brain = (xx / a) ** 2 + (yy / b) ** 2 < 1
        vent = (xx / w) ** 2 + (yy / 0.4) ** 2 < 1
        out[k] = gaussian_filter(c * brain - 0.5 * c * vent, 1.0)
    return out + 0.15 * rng.standard_normal(out.shape)

train = healthy(2000).reshape(2000, -1)
pca = PCA(n_components=20, random_state=0).fit(train) # the optimal linear autoencoder (Thm 7.2)

test = healthy(200)
lesion = np.zeros_like(test, dtype=bool)
for k in range(100):                                  # half of the test subjects get a lesion
    ci, cj = rng.uniform(-0.4, 0.4, 2)
    lesion[k] = ((xx - cj) ** 2 + (yy - ci) ** 2) < 0.12 ** 2
test = test + 0.25 * lesion
flat = test.reshape(200, -1)
recon = pca.inverse_transform(pca.transform(flat))
resid = gaussian_filter(np.abs(flat - recon).reshape(200, n, n), (0, 1, 1))   # anomaly map

print(f"pixel-level AUC of the residual map: {roc_auc_score(lesion.ravel(), resid.ravel()):.3f}")
is_patho = np.arange(200) < 100
print(f"subject-level AUC of max residual:   {roc_auc_score(is_patho, resid.max(axis=(1, 2))):.3f}")
print(f"explained variance of 20 components: {pca.explained_variance_ratio_.sum():.3f}")

The residual map separates lesion pixels well (AUC about 0.94), but the maximum of the map separates patients less well (about 0.77): noise produces high residuals in healthy slices too. Real anatomy varies far more than this toy model, which is the central difficulty of the approach.

7.4 Generative adversarial networks

A GAN (Goodfellow et al., 2014) trains a generator G, which maps noise z to images, against a discriminator D, which tells real images from generated ones: \min_G \max_D\ V(D, G) = \mathbb{E}_{x \sim p_{\text{data}}} \log D(x) + \mathbb{E}_{z} \log\big(1 - D(G(z))\big).

In practice the generator maximises \log D(G(z)) instead, which gives stronger gradients early in training. GAN samples are sharp, but training is unstable and can collapse onto a few modes.

7.5 Image-to-image translation

• Paired data: pix2pix (Isola et al., 2017) conditions both networks on the source image s, and adds to the adversarial loss the term \lambda\, \mathbb{E}\|G(s) - t\|_1 over pairs (s, t) of source and target images.
• Unpaired data: CycleGAN (Zhu et al., 2017) trains two generators G : S \to T and F : T \to S with two discriminators, and adds a cycle-consistency term \mathbb{E}\|F(G(s)) - s\|_1 + \mathbb{E}\|G(F(t)) - t\|_1 .

Wolterink et al. (2017) synthesised CT from brain MRI with unpaired data, useful for radiotherapy planning from MRI alone. Frid-Adar et al. (2018) augmented a liver-lesion classifier with GAN-generated lesions.

The danger is specific to medicine. Cohen, Luck and Honari (2018) showed that distribution-matching losses can hallucinate or remove tumours when the training distributions of the two domains differ: the output looks realistic for the target domain, whatever the input says. A translated image must never be read as a measurement.

8. Validation, interpretation and reproducibility

A medical imaging model is only as good as the evidence that it works on new patients. This chapter shows how evaluations go wrong (leakage, small test sets, unpaired comparisons), how to explain a prediction without fooling oneself, and how to make an experiment reproducible.

8.1 Splitting the data without leakage

Data leakage is any path by which information about the test set reaches training. The most common forms in medical imaging are:

• slices, patches or augmented copies of the same subject in both training and test sets;
• several visits of the same patient split across sets in a longitudinal study;
• preprocessing, feature selection or normalisation statistics computed on all the data before splitting;
• hyperparameters or the number of epochs chosen on the test set.

Wen et al. (2020) reviewed convolutional networks for Alzheimer’s disease classification from MRI and found that a substantial share of the published studies may have suffered from leakage, often through slice-level splits. They showed that such splits inflate accuracy, and released an open-source framework for reproducible evaluation.

Splits must be made at the level of the subject: every image of a patient goes to the same set. Hyperparameters are chosen by cross-validation inside the training set (nested cross-validation), and the test set is used once. An external test set, from other hospitals and scanners, measures the robustness to changes of scanner, protocol and population (the distribution shifts surveyed in ch. 9.4).

In the code below, each subject contributes 10 slices with a strong individual “fingerprint” and a weak disease effect. A one-nearest-neighbour classifier recognises the subject, not the disease.

import numpy as np
from sklearn.model_selection import KFold, GroupKFold, cross_val_predict
from sklearn.neighbors import KNeighborsClassifier
from sklearn.metrics import roc_auc_score, accuracy_score

rng = np.random.default_rng(0)
n_subj, per_subj, d = 100, 10, 50
label = rng.integers(0, 2, n_subj)                       # one diagnosis per subject
finger = rng.standard_normal((n_subj, d)) * 2.0          # subject-specific anatomy
signal = np.zeros(d); signal[:5] = 0.4                   # weak disease effect on 5 features
X = np.concatenate([finger[s] + label[s] * signal + 0.5 * rng.standard_normal((per_subj, d))
                    for s in range(n_subj)])
y = np.repeat(label, per_subj)
groups = np.repeat(np.arange(n_subj), per_subj)          # subject identifier of each slice

clf = KNeighborsClassifier(n_neighbors=1)
splits = {"slice-level split (leaky)": KFold(5, shuffle=True, random_state=0).split(X),
          "subject-level split": GroupKFold(5).split(X, y, groups)}
for name, cv in splits.items():
    pred = cross_val_predict(clf, X, y, cv=list(cv))
    print(f"{name:28s} accuracy = {accuracy_score(y, pred):.3f}")

# bootstrap confidence interval of a test AUC, resampling subjects (n = 100)
score = X[:: per_subj, :5].sum(axis=1)                   # one slice per subject, a simple score
auc = roc_auc_score(label, score)
boot = []
for b in range(2000):
    i = rng.integers(0, n_subj, n_subj)
    if label[i].min() != label[i].max():
        boot.append(roc_auc_score(label[i], score[i]))
lo, hi = np.percentile(boot, [2.5, 97.5])
print(f"AUC = {auc:.3f}, 95 % bootstrap interval [{lo:.3f}, {hi:.3f}]")

The leaky split gives an accuracy of 1.000, the subject-level split 0.527: chance level for this classifier. The second part of the code is used in ch. 8.2.

8.2 The uncertainty of a performance estimate

With n = 100 and \hat p = 0.8, the interval is 0.8 \pm 0.078; with n = 1000, 0.8 \pm 0.025. Varoquaux (2018) showed that, with the sample sizes common in neuroimaging, cross-validation error bars are of the order of 10 points, larger than many published improvements.

For metrics without a simple formula (AUC, Dice, sensitivity at a fixed specificity), the bootstrap resamples the test subjects with replacement many times, recomputes the metric, and takes the 2.5 % and 97.5 % percentiles. In the code above, an AUC of 0.71 on 100 subjects has the interval [0.60, 0.81].

8.3 Comparing models, choosing metrics, calibration

Two models evaluated on the same test cases must be compared with paired tests, which remove the variability shared by both.

For segmentation, the per-case Dice coefficients of two methods are compared by a paired Wilcoxon test or a paired bootstrap. When many methods are compared, the tests must be corrected for multiplicity. Maier-Hein et al. (2018) showed that challenge rankings can change with the metric, the aggregation and the test cases, and the Metrics Reloaded recommendations (Maier-Hein et al., 2024) guide the choice of a metric from the clinical question.

Probabilities matter when a clinician reads them. A model is calibrated if, among the cases given probability 0.8, 80 % are positive. The expected calibration error groups predictions into B bins of confidence and averages the gap: \mathrm{ECE} = \sum_{b=1}^{B} \frac{|I_b|}{n}\, \big|\mathrm{acc}(I_b) - \mathrm{conf}(I_b)\big| .

Modern networks are often over-confident; temperature scaling divides the logits by a scalar fitted on validation data, which improves calibration without changing the ranking (Guo et al., 2017).

8.4 Interpretation

An explanation method assigns to each input voxel a relevance for a given prediction f_c(x) of class c.

• Gradient saliency (Simonyan et al., 2014) uses |\partial f_c / \partial x|.
• Occlusion masks regions and records the drop of the score.
• Grad-CAM (Selvaraju et al., 2017) weights the feature maps A^k of the last convolutional layer by their averaged gradients.

Explanations must themselves be validated. Adebayo et al. (2018) proposed sanity checks: randomise the weights of the trained network, or train it on random labels, and verify that the explanation changes. Some popular methods produce almost the same maps, which then reflect edges of the image rather than the model.

Explanations are most useful for detecting shortcut learning. Zech et al. (2018) found that pneumonia detectors could identify the hospital of a radiograph, and that their performance was often lower on radiographs from other hospitals; DeGrave, Janizek and Lee (2021) found that COVID-19 detectors relied on shortcuts such as laterality markers and features near the image borders, rather than on the lungs. Performance must also be reported by subgroup (sex, age, site): Larrazabal et al. (2020) showed that gender imbalance in training data degrades chest X-ray diagnosis for the under-represented group.

8.5 Reproducibility

A result is reproducible if another team can obtain it again from the same data and code. In practice:

• fix and report the random seeds, and run several seeds, since GPU convolutions are not always deterministic;
• version the code, the data splits (as lists of subject identifiers) and the preprocessing pipeline, and record library versions;
• report the data set, inclusion criteria, labels and their source, following a checklist such as CLAIM (Mongan, Moy and Kahn, 2020);
• release the code and, when allowed, the trained weights; use public data for at least one experiment.

9. Research directions and open problems

This chapter is for self-study. It surveys active research topics that extend the methods of ch. 2–8, with entry points into the literature; the MICCAI, IPMI, ISBI and MIDL conferences named in the official references are where most of this work appears.

9.1 Foundation models

Large models pretrained on broad data and adapted to many tasks are reaching medical imaging. The Segment Anything Model (Kirillov et al., 2023) segments an object from a prompt (a point or a box); work by Ma and coauthors (2024) fine-tuned it on more than a million medical image–mask pairs across many modalities to obtain MedSAM. Moor et al. (2023) describe the broader programme of “generalist” medical models that combine images, text and other data.

Open questions: do such models beat a task-specific nnU-Net when a few hundred annotated cases exist? How should they be validated, given that their pretraining data are rarely documented well enough to exclude overlap with test sets?

9.2 Learning with fewer labels

Self-supervised pretraining learns representations from unlabelled scans. Models Genesis (Zhou et al., 2021) learns to restore deliberately corrupted 3D CT patches; Azizi et al. (2021) showed that contrastive pretraining on unlabelled medical images improved dermatology and chest X-ray classifiers. Radiology reports are a free source of weak labels: Tiu et al. (2022) trained an image–text model on chest radiographs and their reports that detects pathologies without any explicit annotation.

The survey of Cheplygina, de Bruijne and Pluim (2019) organises semi-supervised, multiple-instance and transfer learning in medical imaging. Annotation noise and disagreement between experts remain poorly handled by losses that assume a single true label.

9.3 Privacy and federated learning

Medical data rarely leave the hospital. In federated learning, each site trains on its own data and only model updates are aggregated; Rieke et al. (2020) and Kaissis et al. (2020) review the approach and its privacy limits. Model updates can still leak information about the training data, which motivates secure aggregation and differential privacy, at a cost in accuracy that is not yet well understood for 3D segmentation.

9.4 Distribution shift, causality and fairness

Models degrade when the scanner, protocol, hospital or population changes. Castro, Walker and Glocker (2020) argue that the causal structure of the data-generating process (does the image cause the label, or the reverse?) determines which shifts a model can survive. Oakden-Rayner et al. (2020) showed that average performance can hide clinically important subgroups where models fail (“hidden stratification”).

Gichoya et al. (2022) found that deep networks can predict patients’ self-reported race from medical images, even when clinicians cannot, which raises the question of how such information is used in diagnostic models. Measuring and correcting such effects is an open problem.

9.5 Uncertainty quantification

A clinician needs to know when not to trust a prediction. Deep ensembles (Lakshminarayanan, Pritzel and Blundell, 2017) and Monte Carlo dropout (Gal and Ghahramani, 2016) estimate predictive uncertainty empirically. Conformal prediction gives distribution-free guarantees (Angelopoulos and Bates, 2023).

The guarantee is marginal, averaged over patients, and it fails under the distribution shifts of ch. 9.4: conditional and shift-robust versions are an active topic.

9.6 Diffusion models and new architectures

Diffusion models, which generate images by learning to reverse a gradual noising process, now dominate image synthesis. Chung and Ye (2022) used score-based diffusion models as priors for accelerated MRI reconstruction, combining them with the data-consistency step of ch. 6.5. Wolleb et al. (2022) used them for anomaly detection, and Pinaya et al. (2022) generated synthetic 3D brain MRI with latent diffusion models. They inherit the hallucination risks of ch. 7.5.

Transformers have been combined with U-Nets (TransUNet, Chen et al., 2021; UNETR, Hatamizadeh et al., 2022). Isensee et al. (2024) argued, after re-evaluating many recent architectures under rigorous validation, that well-configured convolutional U-Nets remain competitive: claims of progress depend on the strength of the baseline (ch. 8).

9.7 Open problems

• Validation that predicts clinical benefit: prospective studies and the comparison with clinicians on the decisions that matter, not only on benchmark metrics.
• Robustness to new scanners and hospitals without retraining, and the detection of inputs outside the training distribution.
• Learning from heterogeneous labels: different annotation protocols, partial labels, disagreement between experts.
• Reconstruction guarantees: when can a learned reconstruction be trusted not to invent or erase a lesion?
• Longitudinal and multimodal models that combine images over time with clinical data, as clinicians do.

Exercises

The exercises are sorted by level: ★ applies the course, ★★ combines results or proves a short statement, ★★★ is open-ended. Exercises 10, 12 and 14 are computational.

1. ★ Hounsfield units and windows. Take \mu_{\text{water}} = 0.19 cm⁻¹ and \mu_{\text{air}} = 0. (a) Convert \mu = 0.20 and \mu = 0.38 cm⁻¹ to Hounsfield units. (b) Map these values and fat (-100 HU) through the soft-tissue window of level 40 and width 400. (c) Which tissues can no longer be told apart in this window?

   Solution: (a) By Definition 1.2, 1000 \times 0.01 / 0.19 \approx 53 HU and 1000 \times 0.19 / 0.19 = 1000 HU.

   2. The window maps [-160, 240] HU linearly to [0, 1]: 53 HU gives 213 / 400 \approx 0.53, fat gives 0.15, and bone saturates at 1.

   3. Everything above 240 HU (all bone, contrast agent, metal) maps to 1, and everything below -160 HU (lung, air) to 0. Lungs are read with another window, for instance level -600 and width 1500.

2. ★ Voxels, millimetres and millilitres. A brain MRI of size 182 \times 218 \times 182 has the affine matrix with diagonal (-1, 1, 1, 1) and translation (90, -126, -72). (a) Which voxel is at the world origin? (b) Which voxel contains the point (10, -20, 30) mm? (c) A hippocampus mask has 3500 voxels; give its volume, and its expected number of voxels after resampling to 2 mm.

   Solution: (a) By Definition 1.3, x = 90 - i, y = j - 126, z = k - 72, so the origin is voxel (90, 126, 72).

   2. i = 90 - 10 = 80, j = 106, k = 102.

   3. The voxels are 1 mm³, so the volume is 3.5 mL. A 2 mm voxel has 8 mm³, hence about 3500 / 8 \approx 438 voxels; partial-volume effects at the boundary make the count depend on the interpolation.

3. ★ Predictive values and prior shift. A test has sensitivity 0.95 and specificity 0.85. (a) Compute its positive and negative predictive values at a prevalence of 5 %. (b) A network trained on balanced classes outputs 0.7 for a patient from this population. What probability should be reported?

   Solution: (a) Per 1000 people: 47.5 true positives, 142.5 false positives, 807.5 true negatives and 2.5 false negatives. Hence \mathrm{PPV} = \frac{47.5}{190} = 0.25, \qquad \mathrm{NPV} = \frac{807.5}{810} \approx 0.997 .

   2. By Proposition 2.4 with \pi_0 = 0.5 and \pi = 0.05, the numerator is 0.1 \times 0.7 = 0.07 and the second denominator term is 1.9 \times 0.3 = 0.57. The probability is 0.07 / 0.64 \approx 0.11.

4. ★ IoU and non-maximum suppression. A ground-truth box is (0, 0, 10, 10). A detector outputs A = (1, 1, 11, 11) with score 0.9, B = (0, 2, 10, 12) with score 0.8 and C = (20, 20, 30, 30) with score 0.7. (a) Compute the IoU of A and B with the ground truth and with each other. (b) Apply NMS at threshold 0.5. (c) Give the precision and recall at IoU 0.5, with and without NMS.

   Solution: (a) A and the ground truth intersect on 9 \times 9 = 81, so the IoU is 81 / 119 \approx 0.68.

   B gives 80 / 120 \approx 0.67, and A with B gives 81 / 119 \approx 0.68.

   2. A is kept, B is suppressed (IoU above 0.5 with A), C is kept.

   3. With NMS, A is a true positive and C a false positive: precision 1/2, recall 1. Without NMS, B is a duplicate, counted as a false positive: precision 1/3.

5. ★ Focal loss and FROC. (a) Compute the cross-entropy and the focal loss with \gamma = 2 for p_t = 0.9, 0.5 and 0.1. By how much does the ratio between the hardest and the easiest example change? (b) On 50 scans with 100 nodules, a detector at some threshold hits 85 nodules and raises 60 false positives. Give the point of the FROC curve.

   Solution: (a) The cross-entropy is 0.105, 0.693 and 2.303; the focal loss is 0.00105, 0.173 and 1.865 (Definition 3.3). The hard-to-easy ratio goes from about 22 to about 1800.

   2. Sensitivity 0.85 at 60 / 50 = 1.2 false positives per scan (Definition 3.7).

6. ★★ The AUC by hand. Positives have scores 0.9, 0.7, 0.4 and negatives 0.8, 0.4, 0.3, 0.1. (a) Compute the AUC with the Mann–Whitney count, ties counting one half. (b) Show that the AUC does not change if all scores are replaced by g(\text{score}) with g strictly increasing. (c) Why can a model with a higher AUC still be worse at the operating threshold used in the clinic?

   Solution: (a) There are 12 pairs. The positive 0.9 beats 4 negatives, 0.7 beats 3, and 0.4 beats 2 and ties 1. The AUC is (4 + 3 + 2.5) / 12 \approx 0.792 (Theorem 2.2).

   2. The count only depends on the order of the scores, which g preserves.

   3. Two ROC curves can cross: the AUC averages over all thresholds, and a model can be better on average yet worse in the high-sensitivity region that screening uses.

7. ★★ The gradient of the Dice loss. With \varepsilon = 0, write I = \sum_i p_i g_i, P = \sum_i p_i and G = \sum_i g_i.

   1. Compute \partial L_{\text{Dice}} / \partial p_i. (b) Evaluate it at p = 0 for an image with 10^6 voxels and G = 100 foreground voxels. (c) Compare the accuracy, the mean cross-entropy and the Dice loss of the “all background” prediction p_i = 0.01 for every voxel.

   Solution: (a) Since L = 1 - 2I / (P + G), \frac{\partial L}{\partial p_i} = -\frac{2 g_i (P + G) - 2I}{(P + G)^2} .

   2. At p = 0, I = P = 0: the gradient is -2 / G = -0.02 on foreground voxels and 0 on background voxels, whatever their number.

   3. At threshold 0.5 the prediction has accuracy 99.99\,\%. Its mean cross-entropy is about 0.0105, of which the 100 missed foreground voxels contribute only 100 \times \log(100) / 10^6 \approx 0.00046.

   The Dice loss is 1 - 2 / 10\,100 \approx 0.9998, almost its maximum: it sees the failure that accuracy and the averaged cross-entropy hide (Proposition 4.3).

8. ★★ Folding and diffeomorphisms. (a) In 1D, for which a is \varphi(x) = x + a \sin x invertible with a positive derivative?

   2. Let v(x) = a \sin x be a velocity field. Is \exp(v) a diffeomorphism for a = 2? (c) Compute the Jacobian determinant of \varphi(x, y) = (x + c\, y,\ y + c\, x) for c = 0.5 and c = 1.2.

   Solution: (a) \varphi'(x) = 1 + a \cos x is positive everywhere if and only if |a| < 1; for |a| > 1 it changes sign and \varphi folds (Definition 5.2).

   2. Yes. The flow of \dot x = 2 \sin x never crosses the fixed points k\pi and preserves order, so by Theorem 5.3 its time-one map is a diffeomorphism, for every a.

   3. The determinant is 1 - c^2: 0.75 for c = 0.5, and -0.44 for c = 1.2, an orientation reversal everywhere.

9. ★★ Data consistency as a MAP estimate. (a) Show that the problem of Proposition 6.5 is the MAP estimate of x when the prior is \mathcal{N}(z, \tau^2 I) and the measured k-space values have complex Gaussian noise of variance \sigma^2; express \lambda.

   2. Compute k_j for (\mathcal{F} z)_j = 1 + i, y_j = 2 and \lambda = 3.

   Solution: (a) The negative log-posterior is, up to constants, \frac{1}{2\sigma^2} \|M \mathcal{F} x - y\|^2 + \frac{1}{2\tau^2} \|x - z\|^2 . Multiplying by 2\tau^2 gives Proposition 6.5 with \lambda = \tau^2 / \sigma^2: trust the data more when the noise is small or the denoiser uncertain.

   2. k_j = (1 + i + 6) / 4 = 1.75 + 0.25\, i.

10. ★★ Rician bias, numerically. (a) Simulate the magnitude |a + n_1 + i n_2| with \sigma = 1 for a = 0, 1, 2, 5, and compare its mean with a and with \sqrt{a^2 + \sigma^2}. (b) Estimate \sigma from a background region of 5000 voxels in two ways.

    Hints and expected result: with 10^6 samples, the means are about 1.253, 1.549, 2.273 and 5.102. The bias is large at low signal-to-noise ratio, and the mean approaches \sqrt{a^2 + \sigma^2} at high signal-to-noise ratio.

    2. Proposition 6.3 gives \hat\sigma = \overline{R} / \sqrt{\pi/2}; the second moment gives \hat\sigma = \sqrt{\overline{R^2} / 2}. Both return about 1.00.

11. ★★ VAE and GAN objectives. (a) Compute the Kullback–Leibler divergence of \mathcal{N}(\mu, s^2) from \mathcal{N}(0, 1) for \mu = 1 and s = 0.5.

    2. Let p_{\text{data}} be uniform on [0, 1] and p_g uniform on [2, 3]. Give D^\star, the Jensen–Shannon divergence and V(D^\star, G).

    What happens when p_g moves to [1.5, 2.5], and why is it a problem for training?

    Solution: (a) By the formula of ch. 7.2, \tfrac12 (1 + 0.25 - \log 0.25 - 1) \approx 0.818.

    2. By Proposition 7.4, D^\star = 1 on [0, 1] and 0 on [2, 3]. The supports are disjoint, so \mathrm{JS} = \log 2 and V = -\log 4 + 2 \log 2 = 0.

    Moving p_g to [1.5, 2.5] leaves the supports disjoint: the divergence stays \log 2, so it gives no signal about how close the generator is. This motivates other divergences, such as the Wasserstein distance.

12. ★★ Confidence intervals and a paired test, numerically. A model is correct on 80 of 100 test patients. (a) Compute the normal-approximation interval of Proposition 8.1, a percentile-bootstrap interval, and the Hoeffding radius of Theorem 0.4. (b) A second model is correct on 85; on the discordant cases, the first is right 3 times and the second 8 times. Is the difference significant at 5 %?

    Hints and expected result: (a) The normal interval is [0.722, 0.878], the bootstrap gives about [0.72, 0.88], and Hoeffding’s radius is 0.136.

    2. McNemar’s statistic (Proposition 8.2) is 25 / 11 \approx 2.27, with a p-value of 0.13; the exact binomial test gives 0.23. Five points of accuracy on 100 patients are not significant.

13. ★★ Integrated gradients for a logistic model. Let f(x) = \sigma(w^\top x + b), with \sigma the logistic function, and the baseline x' = 0.

    1. Compute \mathrm{IG}_i(x) in closed form, for w^\top x \neq 0, and check Proposition 8.4. (b) Compare with the gradient saliency when |w^\top x| is large. (c) Propose a sanity check, in the sense of Adebayo et al. (2018), that the attribution x_i alone (the input itself) would fail.

    Solution: (a) Along the path, \partial_i f = w_i\, \sigma'(t\, w^\top x + b). Integrating in t, \mathrm{IG}_i(x) = x_i w_i\, \frac{\sigma(w^\top x + b) - \sigma(b)}{w^\top x}. Summing over i gives \sigma(w^\top x + b) - \sigma(b) = f(x) - f(0).

    2. The gradient w_i\, \sigma'(w^\top x + b) vanishes when the logistic saturates, even though the features drove the decision; IG does not.

    3. Randomise w: IG changes with w, whereas the input-only map does not, so it explains the image, not the model.

14. ★★★ Leakage with a convolutional network. Modify the code of ch. 2.5 so that the data come from 150 “subjects” of 8 slices each: every subject has its own background texture, shared by all its slices, and its own label. Train the CNN with a slice-level split and with a subject-level split, and compare the test AUCs. How does the gap change with the strength of the subject texture and of the lesion?

    Hints and expected result: draw one texture per subject and add small slice-specific noise. Split with GroupShuffleSplit for the subject level.

    The slice-level split gives a higher AUC, because the network can recognise the subject’s texture seen in training. The gap grows with the subject-specific variability and shrinks when the lesion signal is strong (ch. 8.1).

15. ★★★ Designing a multi-centre validation study. A team wants to show that a hippocampus segmentation network trained on one public data set works in three hospitals. Design the study: data and splits, reference segmentations, metrics, confidence intervals, subgroups, comparison with a classical method, and what result would count as success.

    Hints and expected result: train and tune only on the public data, with subject-level cross-validation (ch. 8.1). Hold out the three hospitals entirely as external test sets, with reference segmentations by two raters to measure inter-rater variability.

    Report Dice and HD95 in millimetres per hospital with bootstrap intervals (ch. 4.4, 8.2), and compare with a multi-atlas method by paired tests (ch. 8.3). Analyse age, disease stage and scanner as subgroups (ch. 9.4).

    Success: a per-hospital performance within inter-rater variability, with intervals reported, and a pre-registered analysis plan.

Projects and assessment

The official assessment (master-mva.com) is an article presentation and a project; the Université Paris-Saclay table describes it as partial continuous assessment, and no weights are displayed. The solo subjects below are sized for the course project, and the paper studies among them double as article presentations. The group subjects are for study groups or a longer project, and say how the work divides. The written exam is a proposal for self-assessment, not an official component.

Solo projects

1. [Solo · ≈ 25 h] Learned registration of brain MRI slices (paper study). Reimplement VoxelMorph (Balakrishnan et al., 2019) in 2D on slices of the OASIS brain MRI data (Marcus et al., 2007), as preprocessed and labelled for the Learn2Reg challenge (Hering et al., 2023). Compare with affine registration and with per-pair iterative optimisation of a dense field, then add scaling and squaring.
   • Deliverables: code, a 4–6 page report, a 10-minute talk.
   • Chapters: ch. 5, 8.
   • Success criterion: a mean label Dice at least 5 points above the affine baseline on held-out pairs, under 1 % of folded pixels with scaling and squaring, and a run-time comparison with the iterative method.
2. [Solo · ≈ 20 h] Noise2Noise for magnitude MRI (paper study). Reproduce the Noise2Noise principle (Lehtinen et al., 2018) on slices of the public IXI brain MRI data set: add complex Gaussian noise, train a residual CNN denoiser with clean targets and with independent noisy targets, and compare. Then lower the signal-to-noise ratio and show the Rician bias of Proposition 6.3.
   • Deliverables: code, a report, a talk.
   • Chapters: ch. 6, 8.
   • Success criterion: noisy-target training within 0.5 dB of clean-target training at high SNR, and a measured bias in the background at low SNR, removed by training on complex data.
3. [Solo · ≈ 20 h] A calibrated and explained chest X-ray classifier. Train a small CNN on PneumoniaMNIST from MedMNIST v2 (Yang et al., 2023). Report the AUC with a bootstrap interval, calibrate with temperature scaling (Guo et al., 2017), and compute Grad-CAM maps with the randomisation sanity check of Adebayo et al. (2018).
   • Deliverables: code, a report with reliability diagrams and saliency figures, a talk.
   • Chapters: ch. 2, 8.
   • Success criterion: a lower expected calibration error at an unchanged AUC, and a figure showing whether the maps survive weight randomisation.

Group projects

1. [Group of 2–3 · ≈ 30 h each] Cardiac segmentation and ejection fraction on ACDC. On the ACDC data (Bernard et al., 2018), public after registration, compare a 2D U-Net, a 3D U-Net and, if a free GPU tier allows, nnU-Net, with cross-entropy, Dice and combined losses. Compute the ejection fraction from the segmented volumes.
   • Split: one architecture per member; shared preprocessing and evaluation code.
   • Measure: Dice and HD95 per structure and phase, with bootstrap intervals; Bland–Altman analysis of the ejection fraction.
   • Chapters: ch. 4, 8. Deliverables: code, an 8-page report, a 20-minute talk.
   • Success criterion: a Dice above 0.90 on the left-ventricular cavity, and paired tests between the variants.
2. [Group of 2–3 · ≈ 30 h each] Accelerated MRI reconstruction. On a subset (about 100 training volumes) of the single-coil knee data of fastMRI (Zbontar et al., 2018), available under a free data-use agreement, compare zero filling, compressed sensing with total variation, a U-Net post-processing network and a cascade with data-consistency layers (Schlemper et al., 2018) at acceleration 4 and 8.
   • Split: compressed sensing and evaluation; U-Net; cascade.
   • Measure: PSNR and SSIM, and a hallucination test: insert a small synthetic lesion in k-space and check whether each method preserves it (in the spirit of Antun et al., 2020).
   • Chapters: ch. 6, 8. Deliverables: code, a report, a talk.
   • Success criterion: the cascade beats the post-processing U-Net at equal parameters, and the lesion test is reported for every method.
3. [Group of 2–3 · ≈ 30 h each] Unsupervised anomaly detection in brain MRI (paper study). Following the comparative study of Baur et al. (2021), train an autoencoder and a variational autoencoder on healthy T2-weighted slices of IXI, and detect tumours on BraTS (Menze et al., 2015) after registration to a common template and intensity normalisation.
   • Split: preprocessing and registration; models; evaluation.
   • Measure: pixel-level area under the precision–recall curve and the best achievable Dice, against a simple intensity-threshold baseline.
   • Chapters: ch. 1, 5, 7, 8. Deliverables: code, a report, a talk.
   • Success criterion: reproduce the qualitative conclusions of the paper, and quantify how much of the detection comes from differences between the two data sets rather than from the lesions.

Written exam (proposal)

Format: 2 hours, one handwritten A4 sheet allowed, no electronic device; four exercises, 20 points in total.

1. Question (5 points). Course question. Define data leakage, give three forms of it in medical imaging, and describe a protocol that avoids them. A model has an accuracy of 0.85 on 200 independent test patients: give an approximate 95 % confidence interval.

   Expected: ch. 8.1 (subject-level splits, nested cross-validation, preprocessing inside the training folds, external test set). By Proposition 8.1, the standard error is \sqrt{0.85 \times 0.15 / 200} \approx 0.025, hence 0.85 \pm 0.049.

2. Question (5 points). Computation. A tumour segmentation has 800 true positive, 200 false positive and 100 false negative voxels, with voxels of 1 \times 1 \times 2 mm. Compute the Dice coefficient, the Jaccard index, the predicted and true volumes, and check the relation between Dice and Jaccard. Explain why a Dice of 0.84 means less for a small lesion than for the liver.

   Expected: Dice = 1600 / 1900 \approx 0.842, Jaccard = 800 / 1100 \approx 0.727, and 2J/(1 + J) \approx 0.842 (Proposition 4.3). The volumes are 2.0 mL predicted and 1.8 mL true. Example 4.5: boundary errors weigh more on small structures.

3. Question (5 points). Short proof. State and prove the Noise2Noise theorem. Then explain why training on two noisy magnitude MRI images of the same slice gives a biased denoiser at low signal-to-noise ratio.

   Expected: Theorem 6.4 and its proof. For magnitude images, the noise of the target is not centred given the clean image (Proposition 6.3): the hypothesis fails, and the network learns the mean of the Rician magnitude, which is above the true signal.

4. Question (5 points). Registration. Write the VoxelMorph training loss and explain each term. Why is a low image dissimilarity not evidence of a good registration? Compute the Jacobian determinant of \varphi(x, y) = (x + a \sin x,\ y) and say for which a it folds; how does scaling and squaring avoid folding?

   Expected: the loss of ch. 5.5 (dissimilarity of the warped image, smoothness of the field). An unregularised field can match intensities while tearing anatomy; evaluation uses labels, landmarks and folding. The determinant is 1 + a \cos x: folding for |a| > 1. Theorem 5.3: exponentials of smooth velocity fields are diffeomorphisms.

Grading rubric (proposal)

References

Official references (master-mva.com)

• Goodfellow, I., Bengio, Y. and Courville, A. (2016). Deep Learning. MIT Press.
• Curated lists of public medical imaging data sets named on the course page: https://github.com/beamandrew/medical-data, https://sites.google.com/site/aacruzr/image-datasets, https://github.com/sfikas/medical-imaging-datasets, and the challenges of https://grand-challenge.org/challenges/.
• Software tools: PyTorch — Paszke, A. et al. (2019). PyTorch: an imperative style, high-performance deep learning library. Advances in Neural Information Processing Systems (NeurIPS); TensorFlow — Abadi, M. et al. (2016). TensorFlow: a system for large-scale machine learning. USENIX Symposium on Operating Systems Design and Implementation (OSDI).
• Related conferences: MICCAI (Medical Image Computing and Computer-Assisted Intervention), IPMI (Information Processing in Medical Imaging), ISBI (International Symposium on Biomedical Imaging), MIDL (Medical Imaging with Deep Learning).

Surveys, textbooks and software

• Bushberg, J. T., Seibert, J. A., Leidholdt, E. M. and Boone, J. M. (2012). The Essential Physics of Medical Imaging, 3rd edition. Lippincott Williams & Wilkins.
• Cardoso, M. J. et al. (2022). MONAI: an open-source framework for deep learning in healthcare. arXiv:2211.02701.
• Colliot, O. (ed.) (2023). Machine Learning for Brain Disorders. Neuromethods 197, Humana Press (Springer), open access.
• Kak, A. C. and Slaney, M. (1988). Principles of Computerized Tomographic Imaging. IEEE Press.
• Litjens, G., Kooi, T., Bejnordi, B. E., Setio, A. A. A., Ciompi, F., Ghafoorian, M., van der Laak, J. A. W. M., van Ginneken, B. and Sánchez, C. I. (2017). A survey on deep learning in medical image analysis. Medical Image Analysis, 42, 60–88.
• Pérez-García, F., Sparks, R. and Ourselin, S. (2021). TorchIO: a Python library for efficient loading, preprocessing, augmentation and patch-based sampling of medical images in deep learning. Computer Methods and Programs in Biomedicine, 208, 106236.
• Shen, D., Wu, G. and Suk, H.-I. (2017). Deep learning in medical image analysis. Annual Review of Biomedical Engineering, 19, 221–248.
• Tustison, N. J., Avants, B. B., Cook, P. A., Zheng, Y., Egan, A., Yushkevich, P. A. and Gee, J. C. (2010). N4ITK: improved N3 bias correction. IEEE Transactions on Medical Imaging, 29(6), 1310–1320.
• Zhou, S. K., Greenspan, H. and Shen, D. (eds.) (2017). Deep Learning for Medical Image Analysis. Academic Press.

Classification and detection

• Campanella, G. et al. (2019). Clinical-grade computational pathology using weakly supervised deep learning on whole slide images. Nature Medicine, 25, 1301–1309.
• Ilse, M., Tomczak, J. M. and Welling, M. (2018). Attention-based deep multiple instance learning. International Conference on Machine Learning (ICML).
• Jaeger, P. F., Kohl, S. A. A., Bickelhaupt, S., Isensee, F., Kuder, T. A., Schlemmer, H.-P. and Maier-Hein, K. H. (2020). Retina U-Net: embarrassingly simple exploitation of segmentation supervision for medical object detection. Machine Learning for Health Workshop (ML4H) at NeurIPS 2019, PMLR 116.
• Lambin, P. et al. (2012). Radiomics: extracting more information from medical images using advanced feature analysis. European Journal of Cancer, 48(4), 441–446.
• Lin, T.-Y., Goyal, P., Girshick, R., He, K. and Dollár, P. (2017). Focal loss for dense object detection. IEEE International Conference on Computer Vision (ICCV).
• Payer, C., Štern, D., Bischof, H. and Urschler, M. (2019). Integrating spatial configuration into heatmap regression based CNNs for landmark localization. Medical Image Analysis, 54, 207–219.
• Raghu, M., Zhang, C., Kleinberg, J. and Bengio, S. (2019). Transfusion: understanding transfer learning for medical imaging. Advances in Neural Information Processing Systems (NeurIPS).
• Redmon, J., Divvala, S., Girshick, R. and Farhadi, A. (2016). You only look once: unified, real-time object detection. IEEE Conference on Computer Vision and Pattern Recognition (CVPR).
• Ren, S., He, K., Girshick, R. and Sun, J. (2015). Faster R-CNN: towards real-time object detection with region proposal networks. Advances in Neural Information Processing Systems (NIPS).
• Setio, A. A. A. et al. (2016). Pulmonary nodule detection in CT images: false positive reduction using multi-view convolutional networks. IEEE Transactions on Medical Imaging, 35(5), 1160–1169.
• Setio, A. A. A. et al. (2017). Validation, comparison, and combination of algorithms for automatic detection of pulmonary nodules in computed tomography images: the LUNA16 challenge. Medical Image Analysis, 42, 1–13.

Segmentation

• Bernard, O. et al. (2018). Deep learning techniques for automatic MRI cardiac multi-structures segmentation and diagnosis: is the problem solved? IEEE Transactions on Medical Imaging, 37(11), 2514–2525.
• Çiçek, Ö., Abdulkadir, A., Lienkamp, S. S., Brox, T. and Ronneberger, O. (2016). 3D U-Net: learning dense volumetric segmentation from sparse annotation. Medical Image Computing and Computer-Assisted Intervention (MICCAI).
• Iglesias, J. E. and Sabuncu, M. R. (2015). Multi-atlas segmentation of biomedical images: a survey. Medical Image Analysis, 24(1), 205–219.
• Isensee, F., Jaeger, P. F., Kohl, S. A. A., Petersen, J. and Maier-Hein, K. H. (2021). nnU-Net: a self-configuring method for deep learning-based biomedical image segmentation. Nature Methods, 18, 203–211.
• Long, J., Shelhamer, E. and Darrell, T. (2015). Fully convolutional networks for semantic segmentation. IEEE Conference on Computer Vision and Pattern Recognition (CVPR).
• Menze, B. H. et al. (2015). The multimodal brain tumor image segmentation benchmark (BRATS). IEEE Transactions on Medical Imaging, 34(10), 1993–2024.
• Milletari, F., Navab, N. and Ahmadi, S.-A. (2016). V-Net: fully convolutional neural networks for volumetric medical image segmentation. International Conference on 3D Vision (3DV).
• Ronneberger, O., Fischer, P. and Brox, T. (2015). U-Net: convolutional networks for biomedical image segmentation. Medical Image Computing and Computer-Assisted Intervention (MICCAI).
• Taha, A. A. and Hanbury, A. (2015). Metrics for evaluating 3D medical image segmentation: analysis, selection, and tool. BMC Medical Imaging, 15, 29.

Registration

• Arsigny, V., Commowick, O., Pennec, X. and Ayache, N. (2006). A log-Euclidean framework for statistics on diffeomorphisms. Medical Image Computing and Computer-Assisted Intervention (MICCAI).
• Ashburner, J. (2007). A fast diffeomorphic image registration algorithm. NeuroImage, 38(1), 95–113.
• Avants, B. B., Epstein, C. L., Grossman, M. and Gee, J. C. (2008). Symmetric diffeomorphic image registration with cross-correlation: evaluating automated labeling of elderly and neurodegenerative brain. Medical Image Analysis, 12(1), 26–41.
• Balakrishnan, G., Zhao, A., Sabuncu, M. R., Guttag, J. and Dalca, A. V. (2019). VoxelMorph: a learning framework for deformable medical image registration. IEEE Transactions on Medical Imaging, 38(8), 1788–1800.
• Dalca, A. V., Balakrishnan, G., Guttag, J. and Sabuncu, M. R. (2019). Unsupervised learning of probabilistic diffeomorphic registration for images and surfaces. Medical Image Analysis, 57, 226–236.
• Hering, A. et al. (2023). Learn2Reg: comprehensive multi-task medical image registration challenge, dataset and evaluation in the era of deep learning. IEEE Transactions on Medical Imaging, 42(3).
• Jaderberg, M., Simonyan, K., Zisserman, A. and Kavukcuoglu, K. (2015). Spatial transformer networks. Advances in Neural Information Processing Systems (NIPS).
• Maes, F., Collignon, A., Vandermeulen, D., Marchal, G. and Suetens, P. (1997). Multimodality image registration by maximization of mutual information. IEEE Transactions on Medical Imaging, 16(2), 187–198.
• Rueckert, D., Sonoda, L. I., Hayes, C., Hill, D. L. G., Leach, M. O. and Hawkes, D. J. (1999). Nonrigid registration using free-form deformations: application to breast MR images. IEEE Transactions on Medical Imaging, 18(8), 712–721.
• Sotiras, A., Davatzikos, C. and Paragios, N. (2013). Deformable medical image registration: a survey. IEEE Transactions on Medical Imaging, 32(7), 1153–1190.
• Thirion, J.-P. (1998). Image matching as a diffusion process: an analogy with Maxwell’s demons. Medical Image Analysis, 2(3), 243–260.
• Vercauteren, T., Pennec, X., Perchant, A. and Ayache, N. (2009). Diffeomorphic demons: efficient non-parametric image registration. NeuroImage, 45(1), S61–S72.
• Viola, P. and Wells, W. M. (1997). Alignment by maximization of mutual information. International Journal of Computer Vision, 24(2), 137–154.

Reconstruction and denoising

• Aggarwal, H. K., Mani, M. P. and Jacob, M. (2019). MoDL: model-based deep learning architecture for inverse problems. IEEE Transactions on Medical Imaging, 38(2), 394–405.
• Antun, V., Renna, F., Poon, C., Adcock, B. and Hansen, A. C. (2020). On instabilities of deep learning in image reconstruction and the potential costs of AI. Proceedings of the National Academy of Sciences, 117(48), 30088–30095.
• Batson, J. and Royer, L. (2019). Noise2Self: blind denoising by self-supervision. International Conference on Machine Learning (ICML).
• Chen, H., Zhang, Y., Kalra, M. K., Lin, F., Chen, Y., Liao, P., Zhou, J. and Wang, G. (2017). Low-dose CT with a residual encoder-decoder convolutional neural network. IEEE Transactions on Medical Imaging, 36(12), 2524–2535.
• Gudbjartsson, H. and Patz, S. (1995). The Rician distribution of noisy MRI data. Magnetic Resonance in Medicine, 34(6), 910–914.
• Hammernik, K., Klatzer, T., Kobler, E., Recht, M. P., Sodickson, D. K., Pock, T. and Knoll, F. (2018). Learning a variational network for reconstruction of accelerated MRI data. Magnetic Resonance in Medicine, 79(6), 3055–3071.
• Jin, K. H., McCann, M. T., Froustey, E. and Unser, M. (2017). Deep convolutional neural network for inverse problems in imaging. IEEE Transactions on Image Processing, 26(9), 4509–4522.
• Lehtinen, J., Munkberg, J., Hasselgren, J., Laine, S., Karras, T., Aittala, M. and Aila, T. (2018). Noise2Noise: learning image restoration without clean data. International Conference on Machine Learning (ICML).
• Lustig, M., Donoho, D. and Pauly, J. M. (2007). Sparse MRI: the application of compressed sensing for rapid MR imaging. Magnetic Resonance in Medicine, 58(6), 1182–1195.
• Schlemper, J., Caballero, J., Hajnal, J. V., Price, A. N. and Rueckert, D. (2018). A deep cascade of convolutional neural networks for dynamic MR image reconstruction. IEEE Transactions on Medical Imaging, 37(2), 491–503.
• Zbontar, J. et al. (2018). fastMRI: an open dataset and benchmarks for accelerated MRI. arXiv:1811.08839.

Generative models

• Baur, C., Denner, S., Wiestler, B., Navab, N. and Albarqouni, S. (2021). Autoencoders for unsupervised anomaly segmentation in brain MR images: a comparative study. Medical Image Analysis, 69, 101952.
• Cohen, J. P., Luck, M. and Honari, S. (2018). Distribution matching losses can hallucinate features in medical image translation. Medical Image Computing and Computer-Assisted Intervention (MICCAI).
• Frid-Adar, M., Diamant, I., Klang, E., Amitai, M., Goldberger, J. and Greenspan, H. (2018). GAN-based synthetic medical image augmentation for increased CNN performance in liver lesion classification. Neurocomputing, 321, 321–331.
• Goodfellow, I. et al. (2014). Generative adversarial nets. Advances in Neural Information Processing Systems (NIPS).
• Isola, P., Zhu, J.-Y., Zhou, T. and Efros, A. A. (2017). Image-to-image translation with conditional adversarial networks. IEEE Conference on Computer Vision and Pattern Recognition (CVPR).
• Kingma, D. P. and Welling, M. (2014). Auto-encoding variational Bayes. International Conference on Learning Representations (ICLR).
• Schlegl, T., Seeböck, P., Waldstein, S. M., Schmidt-Erfurth, U. and Langs, G. (2017). Unsupervised anomaly detection with generative adversarial networks to guide marker discovery. Information Processing in Medical Imaging (IPMI).
• Wolterink, J. M., Dinkla, A. M., Savenije, M. H. F., Seevinck, P. R., van den Berg, C. A. T. and Išgum, I. (2017). Deep MR to CT synthesis using unpaired data. Simulation and Synthesis in Medical Imaging (SASHIMI), MICCAI workshop.
• Zhu, J.-Y., Park, T., Isola, P. and Efros, A. A. (2017). Unpaired image-to-image translation using cycle-consistent adversarial networks. IEEE International Conference on Computer Vision (ICCV).

Validation, interpretation and reproducibility

• Adebayo, J., Gilmer, J., Muelly, M., Goodfellow, I., Hardt, M. and Kim, B. (2018). Sanity checks for saliency maps. Advances in Neural Information Processing Systems (NeurIPS).
• DeGrave, A. J., Janizek, J. D. and Lee, S.-I. (2021). AI for radiographic COVID-19 detection selects shortcuts over signal. Nature Machine Intelligence, 3, 610–619.
• Efron, B. and Tibshirani, R. J. (1993). An Introduction to the Bootstrap. Chapman & Hall.
• Guo, C., Pleiss, G., Sun, Y. and Weinberger, K. Q. (2017). On calibration of modern neural networks. International Conference on Machine Learning (ICML).
• Larrazabal, A. J., Nieto, N., Peterson, V., Milone, D. H. and Ferrante, E. (2020). Gender imbalance in medical imaging datasets produces biased classifiers for computer-aided diagnosis. Proceedings of the National Academy of Sciences, 117(23), 12592–12594.
• Maier-Hein, L. et al. (2018). Why rankings of biomedical image analysis competitions should be interpreted with care. Nature Communications, 9, 5217.
• Maier-Hein, L. et al. (2024). Metrics reloaded: recommendations for image analysis validation. Nature Methods, 21, 195–212.
• McNemar, Q. (1947). Note on the sampling error of the difference between correlated proportions or percentages. Psychometrika, 12(2), 153–157.
• Mongan, J., Moy, L. and Kahn, C. E. (2020). Checklist for Artificial Intelligence in Medical Imaging (CLAIM): a guide for authors and reviewers. Radiology: Artificial Intelligence, 2(2), e200029.
• Selvaraju, R. R., Cogswell, M., Das, A., Vedantam, R., Parikh, D. and Batra, D. (2017). Grad-CAM: visual explanations from deep networks via gradient-based localization. IEEE International Conference on Computer Vision (ICCV).
• Simonyan, K., Vedaldi, A. and Zisserman, A. (2014). Deep inside convolutional networks: visualising image classification models and saliency maps. ICLR Workshop.
• Sundararajan, M., Taly, A. and Yan, Q. (2017). Axiomatic attribution for deep networks. International Conference on Machine Learning (ICML).
• Thibeau-Sutre, E., Collin, S., Burgos, N. and Colliot, O. (2023). Interpretability of machine learning methods applied to neuroimaging. In O. Colliot (ed.), Machine Learning for Brain Disorders, Neuromethods 197, Humana Press.
• Varoquaux, G. (2018). Cross-validation failure: small sample sizes lead to large error bars. NeuroImage, 180, 68–77.
• Varoquaux, G. and Cheplygina, V. (2022). Machine learning for medical imaging: methodological failures and recommendations for the future. npj Digital Medicine, 5, 48.
• Varoquaux, G. and Colliot, O. (2023). Evaluating machine learning models and their diagnostic value. In O. Colliot (ed.), Machine Learning for Brain Disorders, Neuromethods 197, Humana Press.
• Wen, J., Thibeau-Sutre, E., Diaz-Melo, M., Samper-González, J., Routier, A., Bottani, S., Dormont, D., Durrleman, S., Burgos, N. and Colliot, O. (2020). Convolutional neural networks for classification of Alzheimer’s disease: overview and reproducible evaluation. Medical Image Analysis, 63, 101694.
• Zech, J. R., Badgeley, M. A., Liu, M., Costa, A. B., Titano, J. J. and Oermann, E. K. (2018). Variable generalization performance of a deep learning model to detect pneumonia in chest radiographs: a cross-sectional study. PLOS Medicine, 15(11), e1002683.

Research directions

• Angelopoulos, A. N. and Bates, S. (2023). Conformal prediction: a gentle introduction. Foundations and Trends in Machine Learning, 16(4).
• Azizi, S. et al. (2021). Big self-supervised models advance medical image classification. IEEE/CVF International Conference on Computer Vision (ICCV).
• Castro, D. C., Walker, I. and Glocker, B. (2020). Causality matters in medical imaging. Nature Communications, 11, 3673.
• Chen, J. et al. (2021). TransUNet: transformers make strong encoders for medical image segmentation. arXiv:2102.04306.
• Cheplygina, V., de Bruijne, M. and Pluim, J. P. W. (2019). Not-so-supervised: a survey of semi-supervised, multi-instance, and transfer learning in medical image analysis. Medical Image Analysis, 54, 280–296.
• Chung, H. and Ye, J. C. (2022). Score-based diffusion models for accelerated MRI. Medical Image Analysis, 80, 102479.
• Gal, Y. and Ghahramani, Z. (2016). Dropout as a Bayesian approximation: representing model uncertainty in deep learning. International Conference on Machine Learning (ICML).
• Gichoya, J. W. et al. (2022). AI recognition of patient race in medical imaging: a modelling study. The Lancet Digital Health, 4(6), e406–e414.
• Hatamizadeh, A. et al. (2022). UNETR: transformers for 3D medical image segmentation. IEEE/CVF Winter Conference on Applications of Computer Vision (WACV).
• Isensee, F. et al. (2024). nnU-Net revisited: a call for rigorous validation in 3D medical image segmentation. Medical Image Computing and Computer-Assisted Intervention (MICCAI).
• Kaissis, G. A., Makowski, M. R., Rückert, D. and Braren, R. F. (2020). Secure, privacy-preserving and federated machine learning in medical imaging. Nature Machine Intelligence, 2, 305–311.
• Kirillov, A. et al. (2023). Segment anything. IEEE/CVF International Conference on Computer Vision (ICCV).
• Lakshminarayanan, B., Pritzel, A. and Blundell, C. (2017). Simple and scalable predictive uncertainty estimation using deep ensembles. Advances in Neural Information Processing Systems (NIPS).
• Ma, J., He, Y., Li, F., Han, L., You, C. and Wang, B. (2024). Segment anything in medical images. Nature Communications, 15, 654.
• Moor, M. et al. (2023). Foundation models for generalist medical artificial intelligence. Nature, 616, 259–265.
• Oakden-Rayner, L., Dunnmon, J., Carneiro, G. and Ré, C. (2020). Hidden stratification causes clinically meaningful failures in machine learning for medical imaging. ACM Conference on Health, Inference, and Learning (CHIL).
• Pinaya, W. H. L. et al. (2022). Brain imaging generation with latent diffusion models. Deep Generative Models workshop (DGM4MICCAI), MICCAI.
• Rieke, N. et al. (2020). The future of digital health with federated learning. npj Digital Medicine, 3, 119.
• Tiu, E., Talius, E., Patel, P., Langlotz, C. P., Ng, A. Y. and Rajpurkar, P. (2022). Expert-level detection of pathologies from unannotated chest X-ray images via self-supervised learning. Nature Biomedical Engineering, 6, 1399–1406.
• Wolleb, J., Bieder, F., Sandkühler, R. and Cattin, P. C. (2022). Diffusion models for medical anomaly detection. Medical Image Computing and Computer-Assisted Intervention (MICCAI).
• Zhou, Z., Sodha, V., Pang, J., Gotway, M. B. and Liang, J. (2021). Models Genesis. Medical Image Analysis, 67, 101840.

Data sets

• Marcus, D. S., Wang, T. H., Parker, J., Csernansky, J. G., Morris, J. C. and Buckner, R. L. (2007). Open Access Series of Imaging Studies (OASIS): cross-sectional MRI data in young, middle aged, nondemented, and demented older adults. Journal of Cognitive Neuroscience, 19(9), 1498–1507.
• Yang, J., Shi, R., Wei, D., Liu, Z., Zhao, L., Ke, B., Pfister, H. and Ni, B. (2023). MedMNIST v2: a large-scale lightweight benchmark for 2D and 3D biomedical image classification. Scientific Data, 10, 41.
• IXI data set: about 600 brain MRI scans of healthy subjects (T1, T2, PD, MRA, diffusion), distributed by Imperial College London at https://brain-development.org/ixi-dataset/.