evaluate continuous outputs as probabilities: the *expected calibration error* (ECE) and the *Brier score*. The ECE measures whether, on samples predicted with a score *s*, the error rate is indeed *s*, in which case the classifier is said to be calibrated. The *Brier score* is minimal when the classifier score is the true probability of the class given the data for an individual, for instance, the probability of the presence of a tumor given the image. These two notions are similar, but it is important to understand that ECE controls average error rates, while Brier score controls individual probabilities, which is much more stringent and more useful to the practitioner [5]. Accurate probabilities of individual predictions can be used for optimal decision-making, e.g., opting for brain surgery only for individuals for which a diagnostic model predicts cancer with high confidence.

A given value of ECE is easy to interpret, as it qualifies probabilities mostly independently of prediction performance. On the other hand, the Brier score accounts for both the quality of probabilities and corresponding binary decisions as a low Brier score captures the ability to give good probabilistic prediction of the output. For any classification problem, there exist many classifiers with 0 expected calibration errors, including some with very poor predictions. On the other hand, even the best possible prediction has a non-zero Brier score, unless the output is a deterministic function of the data. The Brier skill score, a variant of the Brier score, is often used to assess how far a predictor is from the best possible prediction, more independent of the intrinsic uncertainty in the data. The Brier skill score is a rescaled version of the Brier score taking as a reference a reasonable baseline: 1 is a perfect prediction, while negative values mean predictions worse than guessing from class prevalence.

## To Conclude

When assessing a classifier:

* Always look at all the individual metrics: false positives and false negatives are seldom equivalent. Understand the medical problem to know the right trade-off [4].

* Never trust a single summary metric (accuracy, balanced accuracy, ROC AUC, ...).

* Consider the prevalence in your target population. It may be that the prevalence in your testing sample is not representative of that of the target population. In that case, aside from LR+ and LR-, performance metrics computed from the testing sample will not be representative of those in the target population.

## 2.1.2 Multi-class Classification

When there are multiple classes to distinguish, the main difference with two-class classification is that the problem can no longer be separated into a positive class (typically individuals with the medical condition of interest) and a negative class (individuals without). As a consequence, sensitivity and specificity no longer have a meaning
for the whole data, nor do $F_1$ score, or the ROC or precision-recall curves. Accuracy is still defined and easy to compute, but still suffers from its common drawbacks, in particular that it may not be straightforward to interpret in the face of class imbalance.

A classic approach is to aggregate metrics for binary settings considering successively each class as the positive instances and all the others as the negatives, in a form of “one versus all.” There are different approaches to averaging the results for each class. *Macro-averaging* computes the metric, for instance, the ROC AUC, for each class, and then averages the results. One drawback is that it may put too much emphasis on classes that are more infrequent. *Weighted or micro-averaging* combines the results of the different classes weighed by the number of instances of each class. The difference between the two is that *weighted averaging* computes the average of the metric weighted by the number of true instances for each class, while *micro-averaging* computes the metric by adding the number of TPs (resp., TNs, FPs, FNs) across all classes.

Inspecting the confusion matrix extended to multi-class settings gives an interesting tool to understand errors: it displays how many times a given true class is predicted as another (Fig. 5). A perfect prediction has non-zero entries only on the diagonal. The confusion matrix may be interesting to reveal which classes are commonly confused, as its name suggests. In our example, instances that are actually of class C2 are often predicted as of class C3.

Multilabel Classification

Multilabel settings are when the multiple classes are not mutually exclusive: for instance, if an individual can have multiple pathologies. The problem is then to detect the presence or absence of each label for an individual. In terms of evaluation, multilabel settings can be understood as several binary classification problems, and thus the corresponding metrics can be used on each label. As in the multi-class settings, there are different ways to average the results for each label—macro, micro—that put more or less emphasis on the rare labels.

<table><caption>**Fig. 5** Multi-class confusion matrix, for a 3-class problem, C1, C2, C3. Each entry gives the number of instances predicted of a given class, knowing the actual class. A perfect prediction would give non-zero entries only on the diagonal</caption><thead><tr><th rowspan="2"></th><th rowspan="2"></th><th colspan="3">Predicted</th></tr><tr><th>C1</th><th>C2</th><th>C3</th></tr></thead><tbody><tr><th rowspan="3">True</th><th>C1</th><td>133</td><td>0</td><td>0</td></tr><tr><th>C2</th><td>0</td><td>107</td><td>36</td></tr><tr><th>C3</th><td>0</td><td>0</td><td>92</td></tr></tbody></table>
**2.2 Metrics for Regression**

In regression settings, the outcome to predict y is continuous, for instance, an individual's age, cognitive scores, or glucose level. Corresponding error metrics gauge how far the prediction ŷ is from the observed y.

**R2 Score.** The go-to metric here is typically the R2 score, sometimes called explained variance—however, the term R2 score should be preferred, as some authors define explained variance as ignoring bias. Mathematically, the R2 score is the fraction of variance of the outcome y explained by the prediction ŷ, relative to the variance explained by the mean ŷ on the test set:

$$R2 = 1 - \frac{SS(y - \hat{y})}{SS(y - \bar{y})},$$

where SS is the sum of squares on the test data. A strong benefit of this metric is that it comes with a natural scale: an R2 of 1 implies perfect prediction, while an R2 of zero implies a trivial and not very useful prediction. Note that chance-level predictions (as obtained for instance by learning on permuted y) yield slightly negative predictions: indeed, even when the data do not support a prediction of y—as in chance settings—it is impossible to estimate the mean y perfectly and predictions will be worse than the actual mean. In this respect, the R2 score has a different behavior in machine learning settings compared to inferential statistics settings not focused on prediction: *in-sample* (for inferential statistics) versus *out-of-sample* settings (for machine learning). Indeed, when the mean of y is computed on the same data as the model, the R2 score is positive and is the square of the correlation between y and ŷ. This is not the case in predictive settings, and the correlation between y and ŷ should not be used to judge the quality of a prediction [6], because it discards errors on the mean and the scale of the prediction, which are important in practice.

**Absolute Error Measures.** Reporting only the R2 score is not sufficient to characterize well a predictive model. Indeed, the R2 score depends on the variance of the outcome y in the study population and thus does not enable comparing predictive models on different samples. For this purpose, it is important to report also an absolute error measure. The root mean square error (RMSE) and the mean absolute error (MAE) are two of such measures that give an error in the scale of the outcome: if the outcome y is an age in years, the error is also in years. The mean absolute error is easier to interpret. Compared to the root mean square error, the mean absolute error will put much less weight on some rare large deviations. For instance, consider the following prediction error (on 11 observations):
<figure><img src="image_225.png" /><figcaption>**Fig. 6** Visualizing prediction errors—plotting the predicted outcome as a function of the observed one enables to detect structure in the error beyond summary metric. Here the error increases for large values of y, for which there is also a systematic undershoot</figcaption></figure>

$$ \begin{array}{l} \text{error} = [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 100] \\ \text{MAE} = 10 \\ \text{RMSE} \approx 30.17. \end{array} $$

Note that if the error was uniformly equal to the same value (10, for instance), both measures would give the same result.

**Assessing the Distribution of Errors.** The difference between the mean absolute error and the root mean square error arises from the fact that both measures account differently for the tails of the distribution of errors. It is often useful to visualize these errors, to understand how they are structured. Figure 6 shows such visualization: predicted y as a function of observed y. It reveals that for large values of y, the predictive model has a larger prediction error, but also that it tends to undershoot: predict a value that underestimates the observed value. This aspect of the prediction error is not well captured by the summary metrics because there are comparatively much less observations with large y.

**Concluding Remarks on Performance Metrics.** Whether it is in regression or in classification, a single metric is not enough to capture all aspects of prediction performance that are important for applications. Heterogeneity of the error, as we have just seen in our last example, can be present not only as a function of prediction target, but of any aspect of the problem, for instance, the sex of the individuals. Problems related to *fairness*, where some groups (e.g., demographic, geographic, socioeconomic groups) suffer more errors than others, can lead to loss of trust or amplification of inequalities [7]. For these reasons, it may be important to also report error metrics on relevant subgroups, following common medical research practice of stratification.
# 3 Evaluation Strategies

The previous section detailed metrics for assessing the performance of a ML model. We now focus on how to estimate the expected prediction performance of the model with these metrics. Importantly, we draw the difference between evaluating a learning procedure, or learner, and a learned model. While these two questions are often conflated in the literature, the first one must account for uncontrolled fluctuations in the learning procedure, while the second one controls a given model on a target external population. The first question is typically of interest to the methods researcher, to conclude on learning procedures, while the second is central to the medical research, to conclude on the clinical application of a model.

Additional information on validation strategies, seen from the perspective of regulatory science, can be found in Subheading 3 of Chap. 23. We focus here on an accessible discussion of the main concepts to have in mind concerning model evaluation strategies, and Raschka [8] gives a more mathematically detailed coverage of related topics.

## 3.1 Evaluating a Learning Procedure

We first focus on assessing the expected performance of a learning procedure on data drawn from a given population. Here, the model is validated on data with similar characteristics to the one used for training, a validation sometimes called *internal validation*. Most importantly, performance should not be evaluated using the same data that were used for training [6]. Therefore, the first step is to split the data into a training set and a testing set. This should be done before starting any work on the data, be it training a ML model or even doing simple statistics for identifying interesting features. Splitting the data can be done using `sklearn.model_selection train_test_split` or `sklearn.model_selection.ShuffleSplit(n Splits=1)` from scikit-learn. When one simply performs a single split of the data into training and testing set, the validation method is called “hold-out.” One should nevertheless check that the training and testing sets have similar characteristics. More precisely, we want the output variable distribution to be approximately the same in the training and testing sets. This is called stratification. For instance, for classification, the proportion of diseased individuals should approximately be the same in the two sets. To that purpose, use Stratified-ShuffleSplit(n Splits=1). In medical applications, it is recommended to control not only for the disease status but also for other variables, such as sociodemographic information (age, sex, ...) or some relevant clinical variables. It will often be difficult (and it is not even necessary) to obtain almost identical distributions between training and testing sets. In practice, it is often
sufficient to have similar means and variances for continuous vari-
ables and similar proportions for categorical variables. The first two
rows of Fig. 7 illustrate the concepts of “hold-out” and
stratification.

**Non-independent Samples.** Prediction may be performed across non-independent data points, for instance, different points in a time series, or repeated measures of the same individual. In such case, it is important that samples in the train and test sets are independent, which may require selecting separated time windows. Also, the cross-validation should mimic the intended usage of the predictor. For instance, a diagnostic model intended to be applied to new individuals should be evaluated making sure that there are no shared individuals between the train and test sets.

3.1.1 Cross-validation

The split between train and test sets is arbitrary. With the same machine learning algorithm, two different data splits will lead to two different observed performances, both of which are noisy estimates of the expected generalization performance of prediction models built with this learning procedure. A common strategy to obtain better estimates consists in performing multiple splits of the whole dataset into training and testing sets: a so-called *cross-validation* loop. For each split, a model is trained using the training set, and the performances are computed using the testing set. The performances over all the testing sets are then aggregated. Figure 7 displays different cross-validation methods. k-fold cross-validation consists in splitting the data into k sets (called folds) of approximately equal size. It ensures that each sample in the dataset is used exactly once for testing. For classification, `sklearn.model_selection.StratifiedKFold` performs stratified k-fold cross-validation.

In each split, ideally, one would want to have a large training
set, because it usually allows training better performing models,
and a large testing set, because it allows a more accurate estimation
of the performance. But the dataset size is not infinite. Splitting out
10–20% for the test set is a good trade-off [9], which amounts to
k=5 or 10 in a k-fold. With small datasets, to maximize the amount
of train data, it may be tempting to leave out only one observation,
in a so-called *leave-one-out* cross-validation. However, such deple-
tion of the test set gives overall worse estimates of the generaliza-
tion performance. Increasing the number of splits is, however,
useful, and thus another strategy consists in performing a large
number of random splits of the data, breaking from the regularity
of the k-fold. If the number of splits is sufficiently large, all samples
will be approximately used the same number of times for training
and testing. This strategy can be done using `sklearn.model-
selection.StratifiedShuffleSplit(n_splits)` and is
called “Repeated hold-out” or “Monte-Carlo cross-validation.”
<figure><img src="image_226.png" /></figure>

**Hold out**

<figure><img src="image_227.png" /></figure>

**Stratification**

<figure><img src="image_228.png" /></figure>

**k-fold cross-validation (here k=5)**

<figure><img src="image_229.png" /></figure>

**Repeated hold-out**

<figure><img src="image_230.png" /><figcaption>**Fig. 7** Different validation methods, from top to bottom. The first method, called "hold-out," involves a single split of the dataset into training and testing sets. It is thus not a cross-validation method. Stratification is the procedure that controls that the output variable (for instance, disease vs. healthy) has approximately the same distribution in the training and testing sets. k-fold cross-validation consists in splitting the data into k sets (called folds) of approximately equal size. Repeated hold-out consists in performing a large number of random splits of the data</figcaption></figure>
Beyond giving a good estimate of the generalization performance, an important benefit of this strategy is that it enables to study the variability of the performances. However, running many splits may be computationally expensive with models that are slow to train.

### 3.1.2 The Need of an Additional Validation Set

Often, it is useful to make choices on the model to maximize prediction performance: make changes on the architecture, tune hyper-parameters, perform early stopping,.... As the test set performance is our best estimate of prediction performance, it would be natural to run cross-validation and pick the best model. However, in such a situation, the performances reported on the testing set will have an optimistic bias: a data-dependent choice has been made on this test set. There are two main solutions to this issue. The first one is usually applied when the model training is fast and the dataset is of small size. It is called nested cross-validation. It consists in running two loops of cross-validation, one nested into the other. The inner loop serves for hyper-parameter tuning or model selection, while the outer loop is used to evaluate the performance. The second solution is to separate from the whole dataset the test set, which will only be used to evaluate the performances. Then, the remainder of the dataset can be further split into training data and data used to make modeling choices, called the validation set.<fnref n="1" /> Such a procedure is illustrated in Fig. 8. Commonly, the training and validation sets will be used in a cross-validation manner. They can then be used to experiment with different models, tune parameters, .... It is absolutely crucial that the test set is isolated at the very beginning, before any experiment is done. It should be left untouched and used only at the end of the study to report the performances. As for the split between training and validation sets, it is desirable that stratification is done when isolating the test set.

If the dataset is very small, nested cross-validation should be preferred as it gives better testing power than hold-out: all the data are used alternatively for model testing. If the dataset feels too small to split into train, validation, test, it may be too small to conduct a trustworthy machine learning study [10].

### 3.1.3 Various Sources of Data Leakage

Data leakage denotes cases where some information from the training set has “leaked” into the test set. As a consequence, the estimation of the performances is likely to be optimistic. Data leakage can be introduced in many ways, some of which are particularly insidi-

<footnote n="1">¹ In Chapter 23, the validation set is called the tuning set, as it is the standard practice in regulatory science and because it insists on the fact that it should not be used to evaluate the final performance, which should be done on an independent test set. In the present chapter, we use the term validation set as it is the most common in the academic setting.</footnote>
Training, validation and test set

<figure><img src="image_231.png" /></figure>

Use to train the model, experiment with different architectures...

**Fig. 8** A standard approach consists in splitting the whole dataset into training, validation, and test sets. The test set must be isolated from the very beginning, left untouched until the end of the study and only be used to evaluate the performance. The training and validation sets are often used in a cross-validation manner. They can be used to experiment with different architectures and tune parameters

ous and may not be obvious to a researcher that is not familiar with
a specific application field. Below, we describe some common causes
of data leakage. A summary can be found in Box 3.

**Box 3: Some Common Causes of Data Leakage**

- Perform feature selection using the whole dataset.

- Perform dimensionality reduction using the whole dataset.

- Perform parameter selection using the whole dataset or the test set.

- Perform model or architecture search using the whole dataset or the test set.

- Report the performance obtained on the validation set that was used to decide when to stop training (in deep learning).

- For a given patient, put some of its visits in the training set and some in the validation set.

- For a given 3D medical image, put some 2D slices in the training set and some in the validation set.
A first basic cause of data leakage is to use the whole dataset for performing various operations on the data. A very common example is to perform feature selection using the whole dataset and then to use the selected features for model training. A similar situation is when dimensionality reduction is performed on the whole dataset. If this is done in an unsupervised manner (for example, using principal component analysis), it is likely to introduce less bias in the performance estimation because the target is not used. It nevertheless remains, in principle, a bad practice. A common practice in deep learning is to perform early stopping, i.e., use the validation set to determine when to stop training. If this is the case, the validation performances can be overoptimistic, and a separate test dataset should be used to report performance. Another cause of data leakage is when there are multiple longitudinal visits (i.e., the patient is evaluated at several time points) or multiple modalities for a given patient. In such a case, one should never put data from the same patient in both the training and validation sets. For instance, one should not, for a given patient, put the visit at month 0 in the training set and the visit at month 6 in the validation set. Similarly, one should not use the magnetic resonance imaging (MRI) data of a given patient for training and the positron emission tomography (PET) image for validation. A similar situation arises when dealing with 3D medical image. It is absolutely mandatory to avoid putting some of the 2D slices of a given patient in the training set and the rest of the slices in the validation set. More generally, in medical applications, the split between training and test sets should always be done at the patient level. Unfortunately, data leakage is still prevalent in many machine learning studies on brain disorders. For instance, a literature review identified that up to 40% of the studies on convolutional neural networks for automatic classification of Alzheimer's disease from T1-weighted MRI potentially suffered from data leakage [11].

### 3.1.4 Statistical Testing

#### Sources of Variance

Train-test splits, cross-validation, and the like seek to estimate the expected generalization performance of a learning procedure. Keeping test data rigorously independent from algorithm development minimizes the bias of this estimation. However, there are multiple sources of arbitrary variations in these estimates. The most obvious one is the intrinsic randomness of certain aspects of learning procedures, such as the random initial weights in deep learning. Indeed, while fixing the seed of the random number generator may remove the randomness on given train data, this stability is misleading given this choice is arbitrary and not representative of the overall behavior of the machine learning algorithm on the data distribution of interest [12]. A systematic study of machine learning benchmarks [13] shows that their most important sources of variance are:

**Choice of test data/split.**

A given test set is an arbitrary sample of the actual population that we are trying to generalize to. As a result, the corresponding
measure of performance is an imperfect estimate of the actual expected performance. Subheading 3.2, below, gives the resulting confidence intervals for a fixed test set. Using multiple splits, and thus multiple test sets, improves the estimation [13], though it makes computing confidence intervals hard [14].

**Hyper-parameter optimization.**

The choice of hyper-parameters is imperfect, for instance, because of limited resources to tune these hyper-parameters. Another attempt to tune hyper-parameter would lead to a slightly different choice. Thus benchmarks do not give an absolute characterization of a learning procedure but are muddied by imperfect hyper-parameters.

**Random seeds.**

As mentioned above, random choices in a learning procedure—initial weights, random drop-out for neural networks, or boot-straps in bagging—lead to uncontrolled fluctuations in benchmarking results that do not characterize the procedure's ability to generalize to new data.

Conclusions Must Account for Benchmarking Variance

With all these sources of arbitrary variance, the question is: given benchmarks of a learning procedure performance, or improvement, is it likely to generalize reliably to new data or rather to be due to benchmarking fluctuations? Considering, for instance, the performance metrics in Table 1, it seems a safe bet to say that the convolutional neural network outperforms the two others but what about the difference between the two other models? From an application perspective, the question is whether this observed difference is likely to generalize to new data.

To answer this question, we must account for estimation error for the expected generalization performance from the different sources of uncontrolled variance in the benchmarks, as listed above. The first source of error comes from the limited sample size to test the predictions of the different learning procedures. Indeed, suppose that the testing set was composed of 100 samples. In that case, if only 3 more samples had been misclassified by the
<table><caption>**Table 1**
Accuracies obtained by different ML models on a binary classification task. Which model performs best? While it is quite likely that the convolutional neural network outperforms the two other models, it is less clear for the two other models. It seems that the support vector machine results in a slightly higher accuracy but is it due to random fluctuations in the benchmarks? Will the difference carry over to new data?</caption><thead><tr><td>Model</td><td>Accuracy</td></tr></thead><tbody><tr><td>Logistic regression</td><td>0.72</td></tr><tr><td>Support vector machine</td><td>0.75</td></tr><tr><td>Convolutional neural network</td><td>0.95</td></tr></tbody></table>

support vector machine, the two models would have had the same performance. A difference of 3 out of 100 could be easily due to having drawn 3 samples not representative of the population. Other sources of variance are due to how stable the learning pipeline is: sensitivity to hyper-parameters, random initialization, etc.

**Box 4: Statistical Procedure to Characterize a Learner**

1. Perform *k* runs of:

(a) Randomly splitting out a test set

(b) Training the learning procedure on the train set

(c) Measuring the performance *p* on the test set

Choose different values of arbitrary parameters (such as random seeds) on each run, and if enough computing power, run hyper-parameter optimization each time. This results in a set of performance measures $\mathcal{M} = \{m_1, ..., m_k\}$.

2. Use all the values {$m_1, ..., m_k$} to conclude on the performance of the learner:

**Confidence intervals** are given by percentiles of $\mathcal{M}$.

**Standard deviation** of $\mathcal{M}$ can be used to gauge typical variance of performance, as it requires performing a smaller number of runs *k* than percentiles. Standard error should not be used (see text).

**Learner comparison** can be done by comparing two such set of values $\mathcal{M}$ and $\mathcal{M}'$, typically counting the fractions of values in $\mathcal{M}$ that outperform $\mathcal{M}'$ (without any pairing). Statistical procedures such as t-test should not be used (see text).
Training and testing a prediction pipeline multiple times are needed to estimate the variability of the performance measure. The simplest solution is to do this several times while varying the arbitrary factors, such as split between the train and the test or random initialization (see Box 4). The resulting set of performance measures is similar to bootstrap samples and can be used to draw conclusions on the distribution of performances in a test set. Confidence intervals can be computed using percentiles of this distribution. Two learning procedures can be compared by counting the number of times that one outperforms the other: outperforming 75% of the times is typically considered as a reliable improvement [13]. If the available computing power enables training learning procedures only a few times, empirical standard deviations should be used, as they require less runs to estimate. The improvements brought by a learning procedure can then be compared to these standard deviations.

Note these procedures do not perform classic null-hypothesis significance testing, which is difficult here. In particular, the standard error across the various runs should not be used instead of the standard deviation: the standard error is the standard deviation divided by the number of runs. The number of runs can be made arbitrarily large given enough compute power, thus making the standard error arbitrarily small. But in no way does the uncertainty due to the limited test data vanish. This uncertainty can be quantified for a fixed test set—see Subheading 3.2, but in repeated splits or cross-validation, it is difficult to derive confidence intervals because the runs are not independent [14, 15]. In particular, it is invalid to use a standard hypothesis test—such as a T-test—across the different folds of a cross-validation. There are some valid options to perform hypothesis testing in a cross-validation setting [14, 16], but they must be implemented with care.

Another reason not to rely on null-hypothesis testing is that their statistical significance only asserts that the expected performance—or improvement—is non-zero over a test population of infinite size. From a practical perspective, we care about meaningful improvements on test sets of finite size, which is related to the notion of acceptance tests —as opposed to significance—in the Neyman-Pearson framework of statistical testing [17]. Unlike null-hypothesis significance testing, it requires choosing a non-zero difference considered as acceptable, for instance as implicitly set by considering that a new learning procedure should improve upon an existing one 75% of the times—far from chance, which lies at 50%.

**3.2 Generalization to an External Population**

**The Importance of External Validation**

The procedures described above characterize the expected error of a learning procedure applied on a given population. A related, but different, question is that of characterizing the error
of a given predictive model, typically output by a training machine
learning procedure on a study population. That second question,
related to the notion of *external validity*, is important for two
reasons. First, it characterizes the specific predictive model that
will be used in practice, "in production." Indeed, variance in the
learning procedure will lead to arbitrary variation in model perfor-
mance as large as typical improvements achieved by developing
better models [13]. Second, characterizing the model on the *target*
population may be important, as it may differ markedly from the
study population. Indeed, the techniques in the previous section
rely on splitting the initial dataset in training and testing
(or validation) sets; hence, these different sets are by construction
drawn from the same population and have similar characteristics
(data coming from the same hospital/centers/countries, similar
age/sex,...). They only demonstrate the ability of the model to
generalize to new but similar data. To better assess model utility,
guidelines on evaluating clinical prediction models insist on exter-
nal validation using data collected later in time, or in a different
geographical area [18].

Testing whether a prediction model can generalize to dissimilar
data is important as it is all too frequent that the study sample, on
which the model was developed, does not represent the target
population [19]. The target data may, for instance, come from
different hospitals and different countries, be acquired with differ-
ent acquisition devices and protocols or with different sociodemo-
graphic or clinical characteristics than those of the training data. For
instance, it has been shown that the type of MRI scanner can have a
substantial impact on the generalization ability of ML models. To
assess such generalization ability, a common practice is to use one
or several additional datasets for testing, these datasets being
acquired using different protocols and at different sites (Fig. 9).
Most often, these datasets come from other research studies (dif-
ferent from the one used for training). However, research studies
do not usually reflect well clinical routine data. Indeed, in research
studies, the acquisition protocols are often standardized and rigor-
ous data quality control is applied. Moreover, participants may not
be representative of the target population. This can be due to
inclusion/exclusion criteria (for instance, excluding patients with
vascular abnormalities in a study on Alzheimer's disease) or due to
uncontrolled biases. For instance, participants to research studies
tend to have a higher socioeconomic status than the general popu-
lation. Therefore, it is highly valuable to also perform validation on
clinical routine data, whenever possible, as it is more likely to reflect
"real-life" situations. One should nevertheless be aware that a given
clinical routine dataset may come with specificities that may not
generalize to all settings. For instance, data collected within a
specialized center of a university hospital may substantially differ
from that seen by a general practitioner.
<figure><img src="image_232.png" /><figcaption>**Fig. 9** In order to assess the generalization ability of a model under different conditions (such as data coming from different hospitals/countries, acquired with different devices and protocols...), a common practice is to use one or several additional datasets that come from other studies than the one used for training</figcaption></figure>

## Testing Procedures for External Validation

External validation of a predictive model relies on an independent test set and not cross-validation. Statistical testing thus amounts to derive confidence intervals or null-hypothesis significance testing for the metric of interest on this test set, exactly as when characterizing a diagnostic test [20].

For simple metrics that rely on counting successes, such as accuracy, sensitivity, PPV, NPV, the sampling distribution can be deduced from a binomial law. Table 2 gives such confidence intervals for a different set of the test set and different values of the ground-truth accuracy. These can be easily adapted to other counts of errors as follows:

<table><tr><td>Accuracy</td><td>N is the size of the test set</td></tr><tr><td>Sensitivity</td><td>N is the number of negative samples in the test set</td></tr><tr><td>Specificity</td><td>N is the number of positive samples in the test set</td></tr><tr><td>PPV</td><td>N is the number of positively classified test samples</td></tr><tr><td>NPV</td><td>N is the number of negatively classified test samples</td></tr></table>

We believe it is very important to have in mind the typical orders of magnitude reported in Table 2. It is not uncommon to find medical classification studies where the test set size is about a hundred or less. In such a situation, the uncertainty on the estimation of the performance is very high.

These parametric confidence intervals are easy to compute and refer to. But actual confidence intervals may be wider if the samples are not i.i.d. In addition, some interesting metrics, such as AUC ROC, do not come with such parametric confidence interval. A
<table><caption>**Table 2**
Binomial confidence intervals on accuracy (95% CI) for different values of ground-truth accuracy</caption><thead><tr><th>N</th><th>65%</th><th>80%</th><th>90%</th><th>95%</th></tr></thead><tbody><tr><td>100</td><td>[-9.0% 9.0%]</td><td>[-8.0% 8.0%]</td><td>[-6.0% 5.0%]</td><td>[-5.0% 4.0%]</td></tr><tr><td>1000</td><td>[-3.0% 2.9%]</td><td>[-2.5% 2.4%]</td><td>[-1.9% 1.8%]</td><td>[-1.4% 1.3%]</td></tr><tr><td>10,000</td><td>[-0.9% 0.9%]</td><td>[-0.8% 0.8%]</td><td>[-0.6% 0.6%]</td><td>[-0.4% 0.4%]</td></tr><tr><td>100,000</td><td>[-0.3% 0.3%]</td><td>[-0.2% 0.2%]</td><td>[-0.2% 0.2%]</td><td>[-0.1% 0.1%]</td></tr></tbody></table>

general and good option, applicable to all situations, is to approximate the sampling distribution of the metric of interest by bootstrapping the test set [8].

Finally, note that all these confidence intervals assume that the available labels are the ground truth. In practice, medical truth is difficult to establish, and label error may bias the estimation of error rates.

When comparing two classifiers, a McNemar's test is useful to test whether the observed difference in errors can be explained solely by sampling noise [21, 22]. The test is based on the number of samples misclassified by one classifier and not the other, $n_{01}$ and vice versa $n_{10}$. The test statistics is then written $(|n_{01} - n_{10}| - 1)^2 / (n_{01} + n_{10})$; it is distributed under the null as a $\chi^2$ with 1 degree of freedom. To compare classifiers scanning the trade-off between specificity and sensitivity without choosing a specific threshold on their score, one option is to compare areas under the curve of the ROC, using the DeLong test [23] or a permutation scheme to define the null [24].

# 4 Conclusion

Evaluating machine learning models is crucial. Can we claim that a new model outperforms an existing one? Is a given model trustworthy enough to be “deployed,” making decisions in actual clinical settings? A good answer to these questions requires model evaluation experiments adapted to the application settings. There is no one-size-fits-all solution. Multiple performance metrics are often important, chosen to reflect target population and cost-benefit trade-offs of decisions, as discussed in Subheading 2. The prediction model must *always* be evaluated on unseen “test” data, but different evaluation goals lead to procedures to choose these test data. Evaluating a “learner”—a model construction algorithm—leads to cross-validation, while evaluating the fitness
of a given prediction rule—as output by model fitting—calls for left-out data representative of the target population. In all settings, accounting for uncertainty or variance of the performance estimate is important, for instance, to avoid investing in models that bring no reliable improvements.

## Acknowledgements

This work was supported by the French government under management of Agence Nationale de la Recherche as part of the “Investissements d’avenir” program, reference ANR-19-P3IA-0001 (PRAIRIE 3IA Institute), ANR-10-IAIHU-06 (Agence Nationale de la Recherche-10-IA Institut Hospitalo-Universitaire-6), ANR-20-CHIA-0026 (LearnI). We thank Sebastian Raschka for detailed feedback.

## Appendix

**Odds Ratio and Diagnostic Tests Evaluation**

Odds and odds ratio are frequently used in biostatistics and epidemiology, but less in machine learning. Here we give a quick introduction to these topics.

**Odds**

Odds are a measure of likelihood of an outcome: the ratio of the number of events that produce that outcome to the number that do not. The odds $O(a)$ of an outcome $a$ are simply related to the probability $P(a)$ of this outcome:

$$ \text{Odds of } a \qquad O(a) = \frac{P(a)}{1 - P(a)}. \qquad (1) $$

In other words, $O(a)$ is the number of times the event $a$ would occur for each occurrence of the opposite event. This intuitive explanation has led odds to be often used for sports gambling. For instance, if the odds are 3 (or more specifically in gambling terminology 3 : 1) for FC Barcelona vs. Real Madrid, it means that FC Barcelona has a probability of winning against Real Madrid of 75% ($P(a) = \frac{O(a)}{O(a)+1}$). Coming back to diseases, supposing that only a minority of the population is affected, if the odds of the disease are 1%, which can be written as 1 : 100, this means that for every diseased person in the population, there are 100 persons without it. The prevalence is thus $\frac{1}{101} = 0.99\% \approx 1\%$. One can see that when the prevalence is low, it is close to the odds, which is not the case when prevalence gets higher. This is true in general of probabilities and odds: when the probability is low, it is close to the odds.
The odds ratio measures the association between two events, *a* and *b*, which we can arbitrarily call respectively outcome and property. The odds ratio is defined as the ratio of the odds of the outcome in the group where the property holds to that in the group where the property does not hold:

$$ \text{Odds ratio between } a \text{ and } b \qquad OR(a, b) = \frac{O(a|b=+)}{O(a|b=-)}. \qquad (2) $$

To compute the odds ratio, the problem is fully specified by the counts in the following contingency table:

$$ \begin{array}{c|c|c} & \multicolumn{2}{c}{\text{Outcome } a} \\ \cline{2-3} \text{Property } & a+ & a- \\ \cline{2-3} b & & \\ \cline{2-3} b- & n_{++} & n_{-+} \\ & & \\ & n_{+-} & n_{--} \end{array} \qquad (3) $$

The odds are written: $O(a|b=+)=\frac{n_{++}}{n_{-+}}$ and $O(a|b=-)=\frac{n_{+-}}{n_{--}}$; hence, the odds ratio reads

$$ OR(a, b) = \frac{n_{++} n_{--}}{n_{-+} n_{+-}}. \qquad (4) $$

Note that this expression is unchanged swapping the role of *a* and *b*; the odds ratio is symmetric, $OR(a, b) = OR(b, a)$.

## Invariance to Class Sampling

Suppose we have sampled the population selecting with a frequency *f* on the outcome *a*+, for instance, to oversample the positive outcome or the positive property.<fnref n="2" /> In Eq. 4, *n*<sub>++</sub> is replaced by *f* *n*<sub>++</sub> and *n*<sub>+-</sub> by *f* *n*<sub>+-</sub>; however, the factor *f* cancels out and the overall expression of the odds ratio is unchanged. This is a central property of the odds ratio:

The odds ratio is unchanged by sample selection bias on one of the variables (*a* or *b*).

This property is one reason why odds and odds ratio are so central to biostatistics and epidemiology: sampling or recruitment bias is an important concern in these fields. For instance, a case-control study has a very different prevalence as the target population, where the frequency of the disease is typically very low.

<footnote n="2">² Indeed, thankfully, many diseases have a prevalence much lower than 50%, e.g., 1%, which is already considered a frequent disease. Therefore, in order to have a sufficient number of diseased individuals in the sample without dramatically increasing the cost of the study, diseased participants will be oversampled. One extreme example, but very common in medical research, is a case-control study where the number of diseased and healthy individuals is equal.</footnote>
Confusion with Risk Ratio

The odds ratio is often wrongly interpreted as a risk ratio—or
relative risk—which is more easily understood.

The risk ratio is the ratio of the probability of an outcome in a
group where the property holds to the probability of this outcome
in a group where this property does not hold. The risk ratio thus
differs from the odds ratio in that it is expressed for probabilities
and not odds. Even though the values for odds ratio and risk ratio
are often close because, in most diseases being diseased in much less
likely than not, they are fundamentally different because the odds
ratio does not depend on sampling whereas the risk ratio does.

Likelihood Ratio of
Diagnostic Tests or
Classifiers

The likelihood ratio used to characterize diagnostic tests or classi-
fiers is strongly related to the odds ratio introduced above, though
it is not strictly speaking an odds ratio. It is defined as

$$
\text{LR} + = \frac{P(T+ | D+)}{P(T+ | D-)} \qquad (5)
$$

Using the expressions in Box 1 and the fact that P(T+|D+) = 1 - P(T-D+), the LR+ can be written as

$$
\mathrm{LR} + = \frac{\text{Sensitivity}}{1 - \text{Specificity}}. \qquad (6)
$$

Link to Pre-test and Post-test Odds

We can write this in terms of the contingency table in Eq. 3 (the link
to the confusion matrix in Fig. 1 is given by a = D, b = T and thus
n<sub>++</sub> = TP, n<sub>+-</sub> = FP, n<sub>+-</sub> = FN, n<sub>--</sub> = TN):

$$
\mathrm{LR} + = \frac{n_{++}}{n_{++} + n_{+-}} \frac{n_{--} + n_{--}}{n_{-+}} \qquad (7)
$$

$$
= \underbrace{\frac{n_{++}}{n_{-+}}}_{\frac{P(D+|T+)}{P(D-|T+)} = O(D+|T+)} \underbrace{\frac{n_{-+} + n_{--}}{n_{++} + n_{+-}}}_{\frac{P(D-)}{P(D+)} = \frac{1}{O(D+)}} \quad (8)
$$

$$
\mathrm{LR} + = \frac{O(D + | T +)}{O(D +)}. \qquad (9)
$$

Indeed, $O(D + | T +) = \frac{P(D + | T +)}{1 - P(D + | T +)} = \frac{P(D + | T +)}{P(D - | T +)}$ and
$O(D +) = \frac{P(D +)}{1 - P(D +)} = \frac{P(D +)}{P(D -)}$.

$O(D+)$ is called the pre-test odds (the odds of having the
disease in the absence of test information). $O(D+|T+)$ is called
the post-test odds (the odds of having the disease once the test
result is known).
Equation 9 shows how the LR+ relates pre- and post-test odds, an important aspect of its practical interpretation.

## Invariance to Prevalence

If the prevalence of the population changes, the quantities are changed as follows: $n_{++} \to f n_{++}$, $n_{+-} \to f n_{+-}$, $n_{-+} \to (1-f) n_{-+}$, $n_{--} \to (1-f) n_{--}$, affecting LR+ as follows:

$$ LR + = \frac{f n_{++} (1-f) n_{-+} + (1-f) n_{--}}{f n_{++} + f n_{+-}}. \quad (10) $$

The factors f and (1-f) cancel out, and thus the expression of LR + is unchanged for a change of the pre-test frequency of the label (prevalence of the test population). This is alike odds ratio, though the likelihood ratio is not an odds ratio (and does not share all properties; for instance, it is not symmetric).

## References

1. Pedregosa F, et al (2011) Scikit-learn: machine learning in python. J Mach Learn Res 12(85): 2825-2830

2. Powers D (2011) Evaluation: from precision, recall and f-measure to ROC, informedness, markedness & correlation. J Mach Learn Technol 2(1):37-63

3. Naeini MP, Cooper G, Hauskrecht M (2015) Obtaining well calibrated probabilities using Bayesian binning. In: Twenty-Ninth AAAI Conference on Artificial Intelligence

4. Vickers AJ, Van Calster B, Steyerberg EW (2016) Net benefit approaches to the evaluation of prediction models, molecular markers, and diagnostic tests. BMJ 352:i6

5. Perez-Lebel A, Morvan ML, Varoquaux G (2023) Beyond calibration: estimating the grouping loss of modern neural networks. In: ICLR 2023 Conference

6. Poldrack RA, Huckins G, Varoquaux G (2020) Establishment of best practices for evidence for prediction: a review. JAMA Psychiatry 77(5): 534-540

7. Barocas S, Hardt M, Narayanan A (2019) Fairness and machine learning. http://www.fairmlbook.org

8. Raschka S (2018) Model evaluation, model selection, and algorithm selection in machine learning. Preprint arXiv:181112808

9. Varoquaux G, Raamana PR, Engemann DA, Hoyos-Idrobo A, Schwartz Y, Thirion B (2017) Assessing and tuning brain decoders:

cross-validation, caveats, and guidelines. NeuroImage 145:166-179

10. Varoquaux G (2018) Cross-validation failure: small sample sizes lead to large error bars. NeuroImage 180:68-77

11. Wen J, Thibeau-Sutre E, Diaz-Melo M, Samper-González J, Routier A, Bottani S, Dormont D, Durrleman S, Burgos N, Colliot O, et al (2020) Convolutional neural networks for classification of Alzheimer's disease: overview and reproducible evaluation. Med Image Anal 63:101694

12. Bouthillier X, Laurent C, Vincent P (2019) Unreproducible research is reproducible. In: International Conference on Machine Learning, pp 725-734

13. Bouthillier X, Delaunay P, Bronzi M, Trofimov A, Nichyporuk B, Szeto J, Mohammadi Sepahvand N, Raff E, Madan K, Voleti V, et al (2021) Accounting for variance in machine learning benchmarks. Proc Mach Learn Syst 3:747-769

14. Bates S, Hastie T, Tibshirani R (2021) Cross-validation: what does it estimate and how well does it do it? Preprint arXiv:210400673

15. Bengio Y, Grandvalet Y (2004) No unbiased estimator of the variance of k-fold cross-validation. J Mach Learn Res 5(Sep):1089-1105

16. Nadeau C, Bengio Y (2003) Inference for the generalization error. Mach Learn 52(3): 239-281
17. Perezgonzalez JD (2015) Fisher, Neyman-Pearson or NHST? A tutorial for teaching data testing. Front Psychol 6:223

18. Moons KG, Altman DG, Reitsma JB, Ioannidis JP, Macaskill P, Steyerberg EW, Vickers AJ, Ransohoff DF, Collins GS (2015) Transparent reporting of a multivariable prediction model for individual prognosis or diagnosis (tripod): explanation and elaboration. Ann Int Med 162(1):W1–W73

19. Dockès J, Varoquaux G, Poline JB (2021) Pre-venting dataset shift from breaking machine-learning biomarkers. GigaScience 10(9): giab055

20. Shapiro DE (1999) The interpretation of diagnostic tests. Statist Methods Med Res 8(2): 113–134

21. Leisenring W, Pepe MS, Longton G (1997) A marginal regression modelling framework for evaluating medical diagnostic tests. Statist Med 16(11):1263–1281

22. Dietterich TG (1998) Approximate statistical tests for comparing supervised classification learning algorithms. Neural Comput 10(7): 1895–1923

23. DeLong ER, DeLong DM, Clarke-Pearson DL (1988) Comparing the areas under two or more correlated receiver operating characteristic curves: a nonparametric approach. Biometrics 44:837–845

24. Bandos AI, Rockette HE, Gur D (2005) A permutation test sensitive to differences in areas for comparing ROC curves from a paired design. Statist Med 24(18):2873–2893

**Open Access** This chapter is licensed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license and indicate if changes were made.

The images or other third party material in this chapter are included in the chapter's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the chapter's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.
**Reproducibility in Machine Learning for Medical Imaging**

Olivier Colliot, Elina Thibeau-Sutre, and Ninon Burgos

**Abstract**

Reproducibility is a cornerstone of science, as the replication of findings is the process through which they become knowledge. It is widely considered that many fields of science are undergoing a reproducibility crisis. This has led to the publications of various guidelines in order to improve research reproducibility.

This didactic chapter intends at being an introduction to reproducibility for researchers in the field of machine learning for medical imaging. We first distinguish between different types of reproducibility. For each of them, we aim at defining it, at describing the requirements to achieve it, and at discussing its utility. The chapter ends with a discussion on the benefits of reproducibility and with a plea for a nondogmatic approach to this concept and its implementation in research practice.

**Key words** Reproducibility, Replicability, Reliability, Repeatability, Open science, Machine learning, Artificial intelligence, Deep learning, Medical imaging

# 1 Introduction

Reproducibility is at the core of the scientific method. In its general and most common meaning, it corresponds to the ability to reproduce the findings of a given experimental study. This is a necessary (but not sufficient) condition for a scientific statement to become accepted as new knowledge. Let's illustrate this with a simple example, considering the following statement: "the volume of the hippocampus is, on average, smaller in patients with Alzheimer's disease (AD) than in healthy people of comparable age." Such statement was the conclusion of studies which measured such volume from magnetic resonance images (MRI). To the best of our knowledge, the first study to assert this was that of Seab et al [1]. This was later reproduced by many other studies (e.g., [2, 3]). It is now widely accepted, which would not have been the case if the study had proven impossible to reproduce. Note that, as stated above, this is a *necessary* but not a *sufficient* condition. Indeed, there could be other reasons for such statement not to be considered as knowledge. For instance, let's imagine that some other
researchers discover that there is an artifact that is systematically present in the MRI of patients with AD and which leads to erroneous volume estimation. Then, the statement could not be considered new knowledge even though it had been reproduced several times.

Machine learning (ML) is, in part, an experimental science. This is not the case of the entirety of the discipline, part of which is theoretical (for instance, mathematical proofs of convergence or of approximation capabilities of different classes of models) or methodological (the invention of a new approach). Nevertheless, since ML ultimately aims at solving practical problems, its experimental component is essential. Typically, one would want to be able to make statements of the type described above from an experimental study. Here is an example of such statement: “this ML model (for instance, a specific convolutional neural network [CNN] architecture), using MRI data as input, is capable of classifying AD patients and healthy controls with an accuracy superior to 80%.” In order to end an article with such a statement, one needs to conduct an experimental study. For such findings to become knowledge, it needs to be subsequently reproduced. Of course, this statement is unlikely to be universal, and one would want to know under which conditions it holds: for instance, is it restricted to a specific class of MRI scanners, to specific disease stages, to specific age ranges?

**Box 1: Glossary**

The readers will find the definition of the terms we used in the present document.

* **Reproducibility, replicability, repeatability.** In the present document, these will be used as synonyms of reproducibility.

* **Original study.** Study that first showed a finding.

* **Replication study.** Study that subsequently aimed at replicating an original study, with the hope to support its findings.

* **Research artifact.** Any output of scientific research: papers, code, data, protocols. . . . Not to be confused with imaging artifacts which are defects of imaging data.

* **Claims.** The conclusions of a study. Basically a set of statements describing the results and a set of limitations which delineate the boundaries within which the claims are stated (the term “claim” is here used in the broad scientific sense not with the specific meaning it has in the context of regulation of medical devices although the two may be related).

* **Limitations.** A set of restrictions under which the claims may not hold (usually because the corresponding settings have not been explored).

(continued)
**Box 1** (continued)

* **Method.** The ML approach described in the paper, independently of its implementation.

* **Code.** The implementation of the method.

* **Software dependencies.** Other software packages that the main code relies on and which are necessary for its execution.

* **Public data.** Data that can be accessed by anybody with no or little restriction (for instance, the data hosted at https://openneuro.org).

* **Semi-public data.** Data which requires approval of a research project (for instance, the Alzheimer's Disease Neuroimaging Initiative [ADNI] http://www.adni-info.org). The researchers can then use the data only for the intended research purpose and cannot redistribute it.

* **Data split.** Separation into training, validation, and test sets.

* **Data leakage.** Faulty procedure which has led information from the training set to leak into the test set. See refs. 4, 5 for details.

* **Error margins.** A general term for providing the precision of the performance estimates (e.g., standard-error or confidence intervals).

* **Researcher degrees of freedom.** Number of different components (e.g., different architectures, hyperparameter values, subsamples...) which have been tried before arriving to the final method [6]. Too many degrees of freedom tend to produce methods that do not generalize.

* **p-hacking.** A bad practice that involves too many degrees of freedom and which consists in trying many different statistical procedures until a significant p-value is found.

* **Acquisition settings.** Factors that influence the scan of a given patient (imaging device, acquisition parameters, image quality).

* **Image artifacts.** Defects of a medical image, these may include noise, field heterogeneity, motion artifacts, and others.

* **Preregistration.** The deposit of the study protocol prior to performing the study. Limits degrees of freedom and increases likelihood of robust findings.

In the examples above, we have actually illustrated only one on the many possible meanings of reproducibility: the addition of new evidence to support a scientific finding of an *original* study through reproduction under different experimental conditions (*see* Box 1 for a glossary of some of the key concepts used). However, it is also used for very different meanings. In computational sciences, it is
often used for the ability to exactly reproduce the results (i.e., the
exact numbers) in a given study. In sciences which aim at providing
measurements (as is often the case in medical imaging), the word
may be used to describe the variability of a given measurement tool
under different acquisition settings. We shall provide more details
on these different meanings in Subheading 2. Finally, the topics of
reproducibility and open science are obviously related since the
latter favors the former. However, open science encompasses a
broader objective which is to make all research artifacts (code,
data, papers...) openly available for the benefit of the whole society.
Conversely, open research may still be unreproducible (e.g.,
because it has relied on faulty statistical procedures).

There has been increasing concern that science is undergoing a reproducibility crisis [7–10]. This is present in various fields from psychology [11] to preclinical oncology research [12]. ML [13–17], digital medicine [18], ML for healthcare [19, 20], and ML for medical imaging [21] are no exception. The concerns are multifaceted. In particular, they include two substantially different aspects: the report of failures to reproduce previous studies and the observation that many papers do not provide sufficient information for reproducing their results. It is important to have in mind that, while the two may be related, there is not a direct relationship between them: it may very well be that a paper seems to include all the necessary information for reproduction and that reproduction attempts fail (for instance, because the original study had too many degrees of freedom and led to a method that only works on a single dataset, see Subheading 4).

Various guidelines have been proposed to improve research reproducibility. Such guidelines may be general [10] or devoted to specific fields including brain imaging [22–25] and ML for healthcare and life sciences [26, 27]. Moreover many other papers, even though not strictly providing guidelines, provide very valuable pieces of advice for making research more trustworthy and in particular more reproducible (e.g., [14, 19, 28–32]).

This chapter is an introduction to the topic of reproducibility for researchers in the field of ML for medical imaging. It is not meant at providing a replacement for the aforementioned previously published guidelines that we strongly encourage the reader to refer to.

The remainder of the chapter is organized as follows. We first start by introducing different types of reproducibility (Subheading 2). For each of them, we attempt to clearly define it and describe what are the requirements to achieve it and the benefits it can provide (Subheadings 3, 4, 5, and 6). All this information is given with having the field of ML for medical imaging as a target, even though part of it may apply to other fields. Finally, we conclude with a discussion which both describes the benefits of reproducibility but also advocates for a nondogmatic point of view on the topic (Subheading 7).
## 2 The Polysemy of Reproducibility

The term “reproducibility” has been used with various meanings which may range from the exact reproduction of a study with the same material and methods, to the reproduction of a result using new experimental data to the support of a scientific idea using a completely different experimental setup [33, 34]. Moreover, various terms have been introduced including reproducibility, replicability, repeatability, reliability, robustness, generalizability. . . Some of these words, for instance, reproducibility vs replicability, have even been used by some authors with opposite meanings [33, 34]. We will not aim at assigning an unambiguous meaning to each of these words, as we find this of little interest, and will use the term “reproducibility,” “replicability,” and “repeatability” as synonyms. On the other hand, we believe, as many other authors [19, 23, 33, 35], that it is important to distinguish between different types of reproducibility. To that purpose, it is useful to have a *taxonomy* of reproducibility. Below, we describe such a taxonomy. We do not claim that it is novel, as it takes inspiration from other papers [14, 19, 23, 33, 35] nor that it should be universally adopted. Furthermore, boundaries between different types of reproducibility are partly fuzzy. We simply hope that it will be useful for the different concepts that we subsequently introduce and that it will be well adapted to the field of ML for medical imaging.

We distinguish between four main types of reproducibility: **exact reproducibility**, **statistical reproducibility**, **conceptual reproducibility**, and **measurement reproducibility**. We describe those four main types in the following sections. They are also summarized in Fig. 1. As will be explained below, the three first types have relationships with each other (this is why they have the same color in the figure) while the fourth is more separated.

## 3 Exact Reproducibility

**What Is It?** Exact reproducibility aims at reproducing strictly identical results as those of a previously published paper. Concretely, this amounts to being able to reproduce tables and figures as they appear in the original paper following the same procedures as the authors.

**What Does It Require?** Exact reproducibility requires to have access to all components that led to the results including, of course, data and code.
<figure><img src="image_233.png" /><figcaption>**Fig. 1** Different types of reproducibility. Note that, in the case of “statistical” and “conceptual” reproducibility, the terms come from [19] but the exact definition provided in each corresponding section may differ</figcaption></figure>

Access to data is obviously necessary [19, 22, 27]. Open data has been described (together with code and papers) as one of the pillars of open science [22]. It is widely accepted that scientific data should adhere to the FAIR (Findable, Accessible, Interoperable, Reusable) principles (please refer to https://www.go-fair.org/fair-principles
and [36] for more details). Among these principles, accessibility is
often the most difficult to adhere to for medical imaging data
(or healthcare data in general). It is very common in medical papers
that data is mentioned as available upon request. However, a study
has showed that, when data is subsequently requested, many
researchers actually do not comply with the data accessibility state-
ment [37]. This is worrisome, and more transparent ways of data
sharing would be welcome. However, as mentioned above, such
transparent sharing procedures may be difficult to put in place for
healthcare data. In particular, making the data public is often
difficult due to regulatory and privacy constraints [19]. Gorgo-
lewski and Poldrack [22] provide useful pieces of advice to facilitate
sharing, but there are cases where public sharing will remain impos-
sible. In particular, one must distinguish between research data
(acquired as part of a research protocol), which can often be made
*public* or *semi-public*<fnref n="1" /> provided that adequate measures have been
taken at data collection (e.g., adequate participant consent), and
routine clinical data (acquired as part of the routine clinical care of
the patients), which sharing can be much more complicated. It is
important that data is easily findable and that it is shared on a server
which has a long-term maintenance. General purpose data reposi-
tories such as Zenodo<fnref n="2" /> provide a good solution. Another important
aspect is to adhere to community standards for data organization,
so that it can easily be reused by researchers. For brain imaging, the
community standard is BIDS (Brain Imaging Data Structure) [38].<fnref n="3" />
This standard is already very mature and has already been extended
to incorporate other modalities such as microscopy images, for
instance (Microscopy-BIDS [39]). Note that there is an ongoing
proposal to extend it to other organs (MIDS – Medical Imaging
Data Structure [40]<fnref n="4" />). Finally, we would like to draw the attention
to an important point that is often overlooked. Even when a study
relies on public or semi-public data, it is absolutely necessary to
specify which samples (e.g., which participants and which scans)
have been used; otherwise, the study is not reproducible [41]. Ide-
ally, one would provide code to automatically make the data selec-
tion [42] in order to ease the replication.

Another key component is that the code is accessible [19, 22, 27]. Indeed, it would be delusional to think that exactly the same results could be obtained using a reimplementation based on information provided in the paper (even though it is good practice to provide as much detail as possible about the methods in the paper).

<footnote n="1">¹ See glossary Box 1.</footnote>

<footnote n="2">² https://zenodo.org/.</footnote>

<footnote n="3">³ https://bids.neuroimaging.io/.</footnote>

<footnote n="4">⁴ See BIDS extension proposal (BEP) number 25 (BEP025) https://bids.neuroimaging.io/get_involved.html#extending-the-bids-specification</footnote>
Theoretically, it does not mean that the code must come with an
open software license. However, doing so has many additional
benefits such as allowing other researchers to use the code or
parts of it for different purposes. The code should be made accord-
ing to good coding practices which include the use of a versioning
system and adequate documentation [20]. Furthermore, although
not strictly needed for reproducibility, the use of continuous inte-
gration makes the code more robust and eases its long-term main-
tenance. Besides, it is good to ease as much as possible the
installation of dependencies [27]. This can be done with pip<fnref n="5" />
when programming in Python. One can also use containers such
as Docker.<fnref n="6" /> One can find useful additional advice in the *Tips for
Publishing Research Code*.<fnref n="7" /> Note that we are not saying that all
these components should be present in any study or are prerequi-
sites for good research. They constitute an ideal
reproducibility goal.

Sharing well-curated notebooks is also a way to ease reproduc-
ibility of results by other researchers. This can be done through
standard means, but dedicated servers also exist. One can, for
instance, cite an interesting initiative called NeuroLibre<fnref n="8" /> which
provides a preprint server for reproducible data analysis in neuro-
science, in particular providing curated and reviewed Jupyter
notebooks [43].

In ML, sharing the code itself is not enough for exact repro-
ducibility. First, every element of the training procedure should be
stored: this includes the data splits and the criteria for model
selection. Moreover, there usually are non-deterministic compo-
nents so it is necessary to store random seeds [27]. Furthermore,
software/operating system versions, the GPU model/version, and
threading have been deemed necessary to obtain exact reproduc-
ibility [44]. The ClinicaDL software platform provides a framework
for easing exact reproducibility of deep learning for neuroimag-
ing [5].<fnref n="9" /> Although it is targeted at brain imaging, many of its
components and concepts are applicable to medical imaging in
general. Also in the field of brain imaging, NiLearn<fnref n="10" /> facilitates
the reproducibility of statistics and ML. One can also cite Pymia<fnref n="11" />
which provides data handling and validation functionalities for deep
learning in medical imaging [45].

<footnote n="5">⁵ https://pypi.org/project/pip/.</footnote>

<footnote n="6">⁶ https://www.docker.com/.</footnote>

<footnote n="7">⁷ https://github.com/paperswithcode/releasing-research-code.</footnote>

<footnote n="8">⁸ https://neurolibre.org/.</footnote>

<footnote n="9">⁹ https://clinicadl.readthedocs.io/en/latest/.</footnote>

<footnote n="10">¹⁰ https://nilearn.github.io/stable/index.html.</footnote>

<footnote n="11">¹¹ https://github.com/rundherum/pymia.</footnote>
Finally, it may seem obvious, but, even when the code is shared, the underlying theory of the method, all its components, and implementation details need to be present in the paper [14].

It is in principle possible to retrain models identically if the above elements are provided. It nevertheless remains a good practice to share trained models, in order to allow other researchers to check that retraining indeed led to the same results but also to save computational resources. However, models can be attacked to recover training data [27, 46]. This is not a problem when the training data is public. When it is privacy-sensitive, methods to preserve privacy exist [27, 47].

In medical imaging, preprocessing and feature extraction are often critical steps that will subsequently influence the ML results. It is thus necessary to also provide code for such parts. Several software initiatives including BIDSApps [48]<fnref n="12" /> and Clinica [49]<fnref n="13" /> provide ready-to-use tools for preprocessing and feature extraction for various brain imaging modalities. Applicable to many medical imaging modalities, the ITK [50, 51]<fnref n="14" /> framework provides a wide range of processing tools. It can ease the work of researchers who do not want to spend time on preprocessing and feature extraction pipelines and focus on the ML part of their work.

**Why Is It Useful?** It has been claimed that exact reproducibility is of little interest, that pursuing it is a waste of energy of the community, and that its only possible use would be the detection of outright fraud which is rare [52]. We disagree with that view. Let's start with fraud. It may be of low occurrence, although its exact prevalence is difficult to establish. Even so, it is of disastrous consequences as it leads to loss of trust by students, scientists, and the general public. In particular, a survey of 1,576 researchers indicated that 40% of them believe that fraud is a factor that “always/often” contributes to irreproducible research and that 70% of them think that it “sometimes” contributes [7]. Exact reproducibility can certainly contribute to reduce fraud as full transparency obviously makes fraud more difficult. Fraud remains possible (one could forge some data and share it), but it is more difficult to achieve under transparency constraints. Fraud may be rare but errors are much more common. The framework of exact reproducibility eases the detection of errors which is a service to science and even to the authors themselves. In particular, it may help discover “biases and artifacts in the data that were missed by the authors and that cannot be discovered if the data are never made available” [27]. Similarly, it can lead to the discovery of

<footnote n="12">¹² https://bids-apps.neuroimaging.io/.</footnote>

<footnote n="13">¹³ https://www.clinica.run/.</footnote>

<footnote n="14">¹⁴ https://itk.org/.</footnote>
wrong validation schemes, including data leakage or errors in
implementation that make it inconsistent with the methodology
presented in the paper. Overall, it may make progress slower, but it
will definitely make it steadier. However, this does not mean that
exact reproducibility should be aimed in all works or made a
requirement for all publications (see Subheading 7.4).

# 4 Statistical Reproducibility

**What Is It?** Statistical reproducibility aims at reproducing findings under statistically equivalent conditions.<fnref n="15" /> The specific definition may vary, but the following choices are often considered reasonable. The implementation of the method (the code) remains the same. Random components are left random. Regarding the data, the general idea is that the sample would be drawn from the same population. One could, for instance, use subsamples of the original data or another subsample of a larger source population. An interesting case is to study different data splits. A less restrictive view of statistical reproducibility would be to use another dataset whose characteristics are similar to those of the original dataset (for instance, similar age, sex, scanner distributions). Note that the boundary between statistical and conceptual reproducibility (defined in the next section) is fuzzy. We do not believe it is possible to draw exact frontiers that would delimit the statistical variations that are admissible in a statistical reproducibility study. Finally, it is important that those who conduct the statistical replication study clearly indicate which components of variability they study.

**What Does It Require?** Here one needs to distinguish between
two types of factors: those necessary to *attempt* reproducibility and
those that increase the likelihood of *successful* reproducibility.

Regarding the first type, most factors are common with those
for exact reproducibility. Code needs to be accessible so that varia-
tions coming from reimplementation do not impact the replication.
Random seeds, GPU model, or other software/execution para-
meters will not be set to be identical because the aim is precisely
to check if the findings of the study are statistically reproducible
under such variations. Knowing their value in the original study is
nevertheless useful in order to dissect potential reasons for failed
replication. Trained models are in a similar situation: they will
usually not be used for statistical replication (models will be
retrained) but shall prove useful to dissect potential failures. Data

<footnote n="15">¹⁵ We use the term of [19] although with a slightly different (more extensive) meaning.</footnote>
accessibility is also very valuable because it will allow studying
different data splits, or subsamples.

The abovementioned elements make it possible for other
researchers to attempt statistical replication of a given study. On
the other hand, there are features of the original study that will
make such replication more likely to be successful (equivalently, one
could say that the original findings are robust). One important
factor is that the original study reports error margins (reporting
the standard error or equivalently a confidence interval). It is
important in this specific context because statistical reproducibility
does not aim at obtaining (and cannot obtain) exactly the same
results. One wants the results to be *compatible* with original ones:
typically a successful replication would produce results which are
within the error margin of the original study. Beyond the topic of
statistical reproducibility, the report of error margins is of great
importance in general, in particular in the field of ML for medical
imaging, because it provides a precision on the estimates of the
performance. Unfortunately, this practice is still too uncommon in
the ML field as a whole [19]. Even worse, it is not uncommon to
find faulty interpretations of estimates. For instance, one should
never estimate standard errors (SE) from multiple runs of a cross-
validation, as the number of runs can be made arbitrarily large and
as a consequence the SE arbitrarily small (*see* [4]). A very common
example is papers which report empirical standard deviation
(SD) across *k*-folds (or more generally across splits). Unlike what
is quite widely believed, this value does not allow to gauge the
precision of the performance estimation. It provides some insight
on the variability of the learning procedure under variations of the
training and validation sets. Further, keep in mind that when *the
number of splits* is small, such gauge will be very rough. *When the
number of splits is* sufficiently large (and typically using random
splits rather than *k*-fold), it is possible to assess if a “learner” (i.e.,
an ML procedure to perform a task) is superior to another one by
counting the fraction of folds on which it obtains superior perfor-
mance (e.g., 75%) [53]. *See* Chap. 21 for more details on this
question. However, in no case can such procedures estimate the
precision of the performance of the trained model, in other words
the precision of the computed biomarker or computer-aided diag-
nosis tool. This requires an independent test set, from which SE
and confidence intervals can be computed.

**Why Is It Useful?** Statistical replication has many merits. First, by reassessing ML methods using different data splits, one can spot faulty procedures including data leakage which is prevalent in the field of medical imaging [54–57]. *See* refs. 4, 5 for more details on data leakage. Beyond procedures which are clearly wrong, it can also detect lack of robustness to different parameters. One would
consider that the procedure is not statistically replicable if it leads to substantially different results under different train/test data splits, different random seeds, or small changes in hyperparameters. Such an ML algorithm would display poor robustness and would be unlikely to be of future clinical use. Note that, regarding the use of different train/test data splits, these would need to preserve a distribution of metadata (for instance, age, sex, diagnosis...) between train and test that is similar to that of the original study. Most classically, if the original study has stratified the splits, the statistical replication study would also need to stratify the splits. Using different distributions (e.g., not stratified) is also interesting but, in our view, falls within conceptual rather than statistical reproducibility. Furthermore, it is very interesting to attempt replication on a different dataset with statistically equivalent characteristics: for instance, another subsample which has not been used in the original study (but comes from the same larger dataset) or a different dataset but with similar characteristics (e.g., same MRI scanners, similar age, similar disease stage...). Unsuccessful replication may be an indication of overfitting of the dataset of the original study through excessive experimentation with different architectures or hyperparameters which ended up with a method that would work only on this very specific dataset. This is referred to as the *researcher degrees of freedom* [6, 22]. This concept extends beyond the field of ML. It actually comes from experimental sciences where different statistical procedures are tried until a statistically significant result is found, a bad practice known as *p-hacking* [58]. It is important that researchers in our field have this problem in mind. Experimental sciences have proposed *preregistered* and *registered* studies as a potential solution to ban such bad practices. Preregistration means that the research plan is written down and made public before the study starts. It can, for example, be published on the Open Science Framework website.<fnref n="16" /> This mechanism reduces the researcher degrees of freedom and is thus likely to lead to more robust results. Registration goes one step further. The research plan is submitted to a journal and peer-reviewed. Thus (most of) the peer review is done before the results are known. It has the additional advantage of putting more focus on methodological soundness than on the groundbreaking nature of results (for instance, negative results will be published). More details about preregistration and registration can be found in [59]. Preregistration and registration are not yet widely used in ML for medical imaging. Such practices would certainly not fit all studies because they leave no room for methodological creativity. On the other hand, they should be very valuable to experimental

<footnote n="16">¹⁶ https://osf.io.</footnote>
studies aiming at validating ML methods. We believe that, as a
community, we should try to adapt such procedures to our field.

# 5 Conceptual Reproducibility

**What Is It?** Conceptual reproducibility can be seen as the ultimate goal: the one which lead to the consolidation of scientific knowledge. The general idea is to be able to validate the findings under conceptually similar conditions.<fnref n="17" /> Conceptually similar means that the method, the data, and the experiments are compatible with the claims of the original study but they are not identical. We will come back to the notion of claims of a study, and their relationships to generalizability and limitations, later in this section.

**What Does It Require?** Again, we may distinguish between factors that make it possible to *attempt* replication and those that will make it more likely to be *successful*.

In theory, nothing but the original paper should be strictly necessary. Nevertheless, this assumes that the original paper has adhered to the scientific gold standard of providing all details necessary for replication: not only a description of the methods which makes reimplementation possible but a detailed description of the datasets and experimental procedure. It is particularly worrisome that many medical imaging publications do not even report basic demographic statistics [30]. [14] argues that the replication should be independent of the implementation. We agree in principle but believe that such requirement would considerably lower the number of conceptual replication attempts, while more are needed to advance our field in a steadier manner. In practice, it is extremely useful to be able to access the code, not only to save a lot of time but also to make sure that an unsuccessful replication is not due to a faulty reimplementation. The same can be said for trained models. Access to the original data can be useful to dissect the potential reasons for differences in results. In summary, none of the elements of exact reproducibility are required, all of them are welcome.

There are several characteristics of an original study that make it less likely for it to be replicated. Low sample size not only means that it is less likely to find a true effect if it exists but also increases the odds that a positive finding is false [9]. This is not only true in ML but in experimental sciences in general. Ideally, the sample size should be justified by a previous power analysis [24]. Causes for failure of statistical reproducibility also apply here. In particular, too many researcher degrees of freedom increase the likelihood of

<footnote n="17">¹⁷ Again, we use the term of [19] although with a slightly different meaning.</footnote>
having built a method that is overly specific to a dataset. Another problem is that the datasets used in medical imaging ML papers are very often not representative of what would be found in the clinic [30]. Indeed, they often come from research datasets where the inclusion criteria are specific, the medical imaging protocols are harmonized, and the data quality is controlled. Thus, it is necessary to have more studies including clinical routine data (e.g., [60, 61]). Finally, it is very important to have in mind that most scientific findings will not universally replicate but that the replication will only succeed under specific conditions. This is why it is critical that scientific papers precisely define their claims and their limitations. For instance, a claim could be that a given algorithm can segment brain tumors with a Dice of 0.9 ± 0.02 when the MR images are acquired at 3 Tesla and have only minimal artifacts. The same paper would mention as limitations that it is unclear how the algorithm would perform at 1.5 Tesla or with data of lower quality. One can see that stating clear claims and limitations will allow defining the scope of conceptual replication studies. Studies outside that scope would aim at studying generalizability beyond original claims.

**Why Is It Useful?** As mentioned above, conceptual reproducibility is the ultimate goal, the one which, through accumulation of evidence, builds consensus about new scientific knowledge. Its utility in general is thus obvious. More specifically, it provides different benefits. In particular, in the field of ML for medical imaging, it allows studying the generalizability of a method. It is thus a step towards its applicability to the clinic. To that aim, the use of multiple datasets is of paramount importance. This will not only allow ruling out that a method is overly specific to a given dataset. It will allow defining which are the bounds within which the method applies. This includes the machine model, the acquisition parameters, and the data quality. It also includes factors which are unrelated to imaging such as population age, sex, geographic origin, disease severity, and others.

# 6 Measurement Reproducibility

**What Is It?** Measurement reproducibility is the study of the variability of a specific measurement under different acquisition conditions. We are aware that, at first sight, this concept does not fit ideally in our taxonomy (see Subheading 7.1 for a more detailed discussion). Nevertheless, we chose to present it as a separate entity because this is a very common meaning of the word reproducibility
in medical imaging<sup>18</sup> (e.g., [62–69]) and we thus believe that it deserves a special treatment. Here, we consider an algorithm that produces a measurement for each individual patient (for instance, the volume of an anatomical structure computed by a segmentation method). A prototypical example of measurement reproducibility is the test-retest reproducibility: how much does the measure vary when applied to two different scans of the same patient? One can then introduce different variations: scans on the same day or not, scans on the same or different machines, systematic addition of noise or artifacts to the data. ... Finally, some authors call inter-method reproducibility the comparison of different software packages for the measurement of the same anatomical entity [70]. We do not believe this falls within the topic of reproducibility but rather of methods’ comparison.<sup>19</sup>

**What Does It Require?** The code is necessary to make sure that variations do not depend on implementation and to ease the reproducibility study. The trained models are also very welcome to facilitate the process. It is then necessary to have access to test-retest data, meaning different acquisitions of the same patient. As mentioned above: the more varied these different acquisitions, the more extensive the study. Ideally, one would want to have access to scans performed on the same day [62, 67], on different days [65, 66], at different times during the day (e.g., before or after caffeine consumption, a factor which affects functional MRI measures [71]), on different imaging devices [63], and with different acquisition parameters [68]. ... It is unlikely to obtain that many scans for the same patients. A more feasible approach is to study these different factors of variations for different patients. Furthermore, starting with a given image, it is possible to simulate different types of alterations and defects by adding them to the original image. This can be very useful because it allows generating very large numbers of images easily and to control for specific imaging characteristics (such as, e.g., the level of noise or the strength of motion artifacts). Such simulations can involve completely synthetic images called phantoms [72] which mimic real images. It can also be done through the addition of defects to real images [73–75]. Ideally, measurement reproducibility should be performed in different populations of participants separately (for instance, a child with autism spectrum disorder or a patient with Parkinson’s disease is more likely to move during the acquisition, and the image is thus more likely to be affected by motion artifacts).

<footnote n="18"><sup>18</sup> Note that the word is used to evaluate reproducibility of automatic methods across different scans of the same subject but also when a human rater is involved (manual or semi-automated measurements), including intra-rater (measurement twice by the same rater) and inter-rater (two different raters) reproducibility from a single scan.</footnote>

<footnote n="19"><sup>19</sup> It is of interest to compare which of them is the most accurate or robust, with respect to a ground truth. However, as mentioned above, we do not believe it falls within the topic of reproducibility.</footnote>
**Why Is It Useful?** Measurement reproducibility is central for measurement sciences, and medical imaging is one of those. It is an extremely precious information to the user (for instance, the radiologist). Indeed, it provides, at the individual patient level, and ideally for different categories of patients, the precision that they may expect from the measurement tool. There is a wide tradition to perform such reproducibility studies in radiology journals. We believe that it would be very welcome that it becomes more commonplace in the ML for medical imaging community.

# 7 Discussion

## 7.1 About the
Different Types of
Reproducibility

We have presented different types of reproducibility. Our taxonomy is not original nor aims at being universal. The boundaries between types are partly fuzzy. For instance, to which degree replication with a different but similar dataset should be considered *statistical* or *conceptual* reproducibility? We do not believe such questions to be of great importance. Rather, it is fruitful, following Gundersen and Kjensmo [14] and Peng [76], to consider reproducibility as a spectrum. In particular, one can consider that the first three types provide increasing support for a finding: conceptual provides more support than statistical which in turns provides more support than exact. The amount of components necessary to perform them is in the reverse order: exact requires more than statistical which requires more than conceptual. Does it mean that only conceptual reproducibility matters? Absolutely not. As we mentioned, other types of reproducibility are necessary to dissect why a given replication has failed as well as to better specify the bounds within which a scientific claim is valid. Last but not least, exact reproducibility also helps build trust in science.

We must admit that measurement reproducibility does not fit very well in this landscape. Moreover, one could also argue that it is a type of conceptual reproducibility, which is partly true as it aims at studying the reproducibility when varying the input data. We nevertheless believe it deserves a special treatment, for several reasons. First, here reproducibility is studied at the individual (i.e., patient) level and not at the population level. Also, the emphasis is on the measurement rather than the finding. Even if it has its role in the building of scientific knowledge, it has specific practical implications for the user. Moreover, as mentioned above, this is actually the most widely used meaning of reproducibility in medical imaging, and it seemed important that the reader is acquainted with it.
## 7.2 The Many Benefits of Reproducibility

"Der Weg ist das Ziel" is a German saying which can be roughly translated as: "the path is the goal." Indeed, reproducibility allows researchers to discover many new places down the road before reaching the final destination. Even if this destination is never reached, the benefits of the travel are of major importance. Let us try to list some of them.

There are many individual benefits for researchers and labs. An important one is that aiming at reproducible research results in *reusable* research artifacts. How agreeable it is for a researcher to easily reuse an old code for a new project! How useful it is for a research lab to have data organized according to community standards making it easier to reuse and share! Moreover, papers that come with shared data [22, 77, 78] or code attract [79], on average, more citations. Thus aiming at reproducibility is also in the researchers' self-interest.

There are also considerable benefits for the scientific community as a whole. As mentioned before, reproducible research is often associated to open code, open data, and available trained models. This allows researchers not only to use them to perform replication studies but also to use these research artifacts for completely different purposes such as building new methods or conducting analysis on pooled datasets. In the specific case of ML for medical imaging, it also allows assessing independently the influence of preprocessing, feature extraction, and ML method. This is particularly important when claims of superiority of a new ML method are made, but the original paper uses overly specific preprocessing steps.

Of course, at the end of the path, the goal itself brings many benefits. These have already largely described in the previous sections so we will just mention them briefly. Conceptual replication studies are necessary for corroborating findings and thus building new scientific knowledge. Statistical replication allows ensuring that results are not due to cherry picking. Exact replication allows detecting errors and increases trust in science in general.

## 7.3 Awareness Is Rising

Throughout this chapter, we have referred to numerous papers, resources, and tools that demonstrate that awareness regarding reproducibility has strongly risen in the past years.

Various papers and studies have highlighted the lack of reproducibility in different fields (e.g., [11, 15, 19]). In machine learning for medical imaging, Simko et al. [21] have studied the reproducibility of methods (mainly code availability and usability thus restricted to exact reproducibility) published at the Medical Imaging with Deep Learning (MIDL) conference from 2018 to 2022
and found that about 20% of papers came with a repository that was deemed reproducible.

Various papers have been published providing advice and guidelines [10, 17, 22–27]. Some of the guidelines include reproducibility checklists. Some checklists are associated to a specific journal or conference and are provided to the reviewers so that they can take these aspects in consideration when evaluating papers. One can cite, for example, the MICCAI (Medical Image Computing and Computer-Assisted Intervention conference) reproducibility checklist.<fnref n="20" />,<fnref n="21" />

Finally, it is very important that reproducibility studies, assessing all aspects of reproducibility (exact, statistical, conceptual, measurement), are performed, published, and widely read. Unfortunately, it is still easier to publish in a high-impact journal a study that is not reproducible but describes exciting results than a replication study. The good news is that this is starting to change. Reproducibility challenges have been proposed in various fields including machine learning<fnref n="22" /> and medical image computing [80]. In the field of neuroimaging, the journal NeuroImage: Reports publishes Open Data Replication Reports.<fnref n="23" /> the Organization for Human Brain Mapping has a replication award.<fnref n="24" /> and the MRITogether workshop<fnref n="25" /> emphasizes reproducibility.

**7.4 One Size Does Not Fit All**

We hope the reader is now convinced of the benefits of aiming towards reproducible research. Does it mean that reproducibility requirements should be the same for all studies? We strongly believe the opposite. To take an extreme example, requiring all studies to be exactly reproducible with minimal efforts (like with running a single command) would be an awful idea. We believe, on the contrary, that reproducibility efforts should vary according to many factors including the type of study and the context in which it is performed. One would probably not have the same level of requirement for a methodological paper and for an extensive medical application with strong claims about clinical applicability. For the former, one may be satisfied with an experiment on a single or a few datasets. For the later, one would expect the study to include

<footnote n="20">²⁰ https://miccai2021.org/files/downloads/MICCAI2021-Reproducibility-Checklist.pdf.</footnote>

<footnote n="21">²¹ https://github.com/JunMa11/MICCAI-Reproducibility-Checklist.</footnote>

<footnote n="22">²² https://paperswithcode.com/rc2022.</footnote>

<footnote n="23">²³ https://www.journals.elsevier.com/neuroimage-reports/infographics/neuroimage-reports-presents-open-data-replication-reports?utm_campaign=STMJ_176479_SC&utm_medium=email&utm_acid=268008024&SID=&dgcid=STMJ_176479_SC&CMX_ID=&utm_in=DM292849&utm_source=AC_.</footnote>

<footnote n="24">²⁴ https://www.humanbrainmapping.org/i4a/pages/index.cfm?pageid=3731.</footnote>

<footnote n="25">²⁵ https://mritogether.esmrmb.org/.</footnote>
multiple datasets with varying characteristics and a comprehensive
assessment of generalizability under different factors such as imag-
ing devices and acquisition parameters. Also, there are some cases
where sharing the code is not desired (e.g., because an industrial
application is foreseen) or where the code will not adhere to best
development practices because it is just a prototype to test a new
methodology. Nevertheless, sharing weakly documented code is
always better than no sharing at all. Similarly, there are cases
where data sharing is difficult or even impossible due to regulatory
constraints. As mentioned above, reproducibility is a spectrum.
Where a given study should lie in this spectrum should depend on
the type of study and the constraints the researchers face.

We thus advocate for a nondogmatic approach to reproducibil-
ity. Guidelines are extremely useful, but they should not be carved
in stone. Also, we believe that the requirements should be assessed
by the reviewers on a case-by-case basis. Indeed, what matters is
that the reproducibility level matches the claims made in the paper.
Of course, it is a good thing that journals and conferences provide
requirements for reporting essential information. It is helpful to
researchers and makes the community progress towards better
science. Also, some bad practices such as data leakage or
p-hacking need to be banished. But we believe that very high
reproducibility requirements (e.g., requiring that exact reproduc-
ibility is feasible) at the level of a given journal or conference would
be counterproductive. Finally, we like the idea of a badging system
[27] which would tag papers according to their reproducibility
level. It remains to be seen how such system should be
implemented.

To conclude, we firmly believe that it is essential for researchers
and students in the field of ML for medical imaging to be trained to
the concepts and practice of reproducibility. It will be beneficial to
them as well as to the community in general. But this does not
mean that researchers should aim at perfect reproducibility in all
their studies. Diversity in research approaches and practices is also a
factor that drives science forward and which should be preserved.

Acknowledgements

The authors are grateful to G. Varoquaux for pointing them
towards useful references. This work was supported by the French
government under management of Agence Nationale de la
Recherche as part of the "Investissements d'avenir" program, ref-
erence ANR-19-P3IA-0001 (PRAIRIE 3IA Institute) and
ANR-10-IAIHU-06 (Agence Nationale de la Recherche-10-IA
Institut Hospitalo-Universitaire-6). ETS acknowledges funding from the 4TU Precision Medicine program supported by High Tech for a Sustainable Future, a framework commissioned by the four Universities of Technology of the Netherlands.

## References

1. Seab J, Jagust W, Wong S, Roos M, Reed BR, Budinger T (1988) Quantitative NMR measurements of hippocampal atrophy in Alzheimer's disease. Magn Reson Med 8(2):200-208

2. Lehericy S, Baulac M, Chiras J, Pierot L, Martin N, Pillon B, Deweer B, Dubois B, Marsault C (1994) Amygdalohippocampal MR volume measurements in the early stages of Alzheimer disease. Am J Neuroradiol 15(5): 929-937

3. Jack CR, Petersen RC, Xu YC, Waring SC, O'Brien PC, Tangalos EG, Smith GE, Ivnik RJ, Kokmen E (1997) Medial temporal atrophy on MRI in normal aging and very mild Alzheimer's disease. Neurology 49(3):786-794

4. Varoquaux G, Colliot O (2022) Evaluating machine learning models and their diagnostic value. HAL preprint hal-03682454. https://hal.archives-ouvertes.fr/hal-03682454/

5. Thibeau-Sutre E, Diaz M, Hassanaly R, Routier A, Dormont D, Colliot O, Burgos N (2022) ClinicaDL: an open-source deep learning software for reproducible neuroimaging processing. Comput Methods Prog Biomed 220:106818

6. Simmons JP, Nelson LD, Simonsohn U (2011) False-positive psychology: undisclosed flexibility in data collection and analysis allows presenting anything as significant. Psychol Sci 22: 1359-1366

7. Baker M (2016) 1,500 scientists lift the lid on reproducibility. Nature 533:452-454

8. Gundersen OE (2020) The reproducibility crisis is real. AI Mag 41(3):103-106

9. Ioannidis JP (2005) Why most published research findings are false. PLoS Med 2(8): e124

10. Begley CG, Ioannidis JP (2015) Reproducibility in science: improving the standard for basic and preclinical research. Circ Res 116(1): 116-126

11. Collaboration OS (2015) Estimating the reproducibility of psychological science. Science 349(6251):aac4716

12. Begley CG (2013) An unappreciated challenge to oncology drug discovery: pitfalls in preclinical research. Am Soc Clin Oncol Educ Book 33(1):466-468

13. Sonnenburg S, Braun ML, Ong CS, Bengio S, Bottou L, Holmes G, LeCunn Y, Muller KR, Pereira F, Rasmussen CE et al (2007) The need for open source software in machine learning. J Mach Learn Res 8(81):2443-2466

14. Gundersen OE, Kjensmo S (2018) State of the art: reproducibility in artificial intelligence. In: Proceedings of the AAAI Conference on Artificial Intelligence, vol 32

15. Hutson M (2018) Artificial intelligence faces reproducibility crisis. Science 359(6377): 725-726

16. Haibe-Kains B, Adam GA, Hosny A, Khodakarami F, Waldron L, Wang B, McIntosh C, Goldenberg A, Kundaje A, Greene CS, et al (2020) Transparency and reproducibility in artificial intelligence. Nature 586(7829):E14-E16

17. Pineau J, Vincent-Lamarre P, Sinha K, Larivière V, Beygelzimer A, d'Alché Buc F, Fox E, Larochelle H (2021) Improving reproducibility in machine learning research: a report from the neurips 2019 reproducibility program. J Mach Learn Res 22:1-20

18. Stupple A, Singerman D, Celi LA (2019) The reproducibility crisis in the age of digital medicine. NPJ Digit Med 2(1):1-3

19. McDermott M, Wang S, Marinsek N, Ranganath R, Ghassemi M, Foschini L (2019) Reproducibility in machine learning for health. arXiv preprint arXiv:190701463

20. Beam AL, Manrai AK, Ghassemi M (2020) Challenges to the reproducibility of machine learning models in health care. JAMA 323(4): 305-306

21. Simko A, Garpebring A, Jonsson J, Nyholm T, Löfstedt T (2022) Reproducibility of the methods in medical imaging with deep learning. arXiv preprint arXiv:221011146

22. Gorgolewski KJ, Poldrack RA (2016) A practical guide for improving transparency and reproducibility in neuroimaging research. PLoS Biol 14(7):e1002506

23. Nichols TE, Das S, Eickhoff SB, Evans AC, Glatard T, Hanke M, Kriegeskorte N, Milham MP, Poldrack RA, Poline JB et al (2017) Best practices in data analysis and sharing in neuroimaging using MRI. Nat Neurosci 20(3): 299-303
24. Poldrack RA, Baker CI, Durnez J, Gorgolewski KJ, Matthews PM, Munafò MR, Nichols TE, Poline JB, Vul E, Yarkoni T (2017) Scanning the horizon: towards transparent and reproducible neuroimaging research. Nat Rev Neurosci 18(2):115-126

25. Niso G, Botvinik-Nezer R, Appelhoff S, De La Vega A, Esteban O, Etzel JA, Finc K, Ganz M, Gau R, Halchenko YO et al (2022) Open and reproducible neuroimaging: from study inception to publication. Neuroimage 263:119623

26. Turkyilmaz-van der Velden Y, Dintzner N, Teperek M (2020) Reproducibility starts from you today. Patterns 1(6):100099

27. Heil BJ, Hoffman MM, Markowetz F, Lee SI, Greene CS, Hicks SC (2021) Reproducibility standards for machine learning in the life sciences. Nat Methods 18(10):1132-1135

28. Varoquaux G (2018) Cross-validation failure: small sample sizes lead to large error bars. Neuroimage 180:68-77

29. Button KS, Ioannidis J, Mokrysz C, Nosek BA, Flint J, Robinson ES, Munafò MR (2013) Power failure: why small sample size undermines the reliability of neuroscience. Nat Rev Neurosci 14(5):365-376

30. Varoquaux G, Cheplygina V (2022) Machine learning for medical imaging: methodological failures and recommendations for the future. NPJ Digit Med 5(1):1-8

31. Bouthillier X, Laurent C, Vincent P (2019) Unreproducible research is reproducible. In: International Conference on Machine Learning, PMLR, pp 725-734

32. Langer SG, Shih G, Nagy P, Landman BA (2018) Collaborative and reproducible research: goals, challenges, and strategies. J Digit Imaging 31(3):275-282

33. Goodman SN, Fanelli D, Ioannidis JP (2016) What does research reproducibility mean? Sci Transl Med 8(341):341ps12-341ps12

34. Plesser HE (2018) Reproducibility vs. replicability: a brief history of a confused terminology. Front Neuroinfor 11:76

35. McDermott MB, Wang S, Marinsek N, Ranganath R, Foschini L, Ghassemi M (2021) Reproducibility in machine learning for health research: still a ways to go. Sci Transl Med 13(586):eabb1655

36. Wilkinson MD, Dumontier M, Aalbersberg IJ, Appleton G, Axton M, Baak A, Blomberg N, Boiten JW, da Silva Santos LB, Bourne PE, et al (2016) The FAIR guiding principles for scientific data management and stewardship. Sci Data 3(1):1-9

37. Gabelica M, Bojčić R, Puljak L (2022) Many researchers were not compliant with their published data sharing statement: mixed-methods study. J Clin Epidemiol 150:33-41

38. Gorgolewski KJ, Auer T, Calhoun VD, Craddock RC, Das S, Duff EP, et al (2016) The brain imaging data structure, a format for organizing and describing outputs of neuroimaging experiments. Sci Data 3(1):1-9

39. Bourget MH, Kamentsky L, Ghosh SS, Mazzamuto G, Lazari A, Markiewicz CJ, Oostenveld R, Niso G, Halchenko YO, Lipp I, et al (2022) Microscopy-BIDS: an extension to the Brain imaging data structure for microscopy data. Front Neurosci 16: e871228

40. Saborit-Torres J, Saenz-Gamboa J, Montell J, Salinas J, Gómez J, Stefan I, Caparrós M, García-García F, Domenech J, Manión J, et al (2020) Medical imaging data structure extended to multiple modalities and anatomical regions. arXiv preprint arXiv:201000434

41. Cuingnet R, Gerardin E, Tessieras J, Auzias G, Lehéricy S, Habert MO, Chupin M, Benali H, Colliot O (2011) Automatic classification of patients with Alzheimer's disease from structural MRI: a comparison of ten methods using the ADNI database. Neuroimage 56(2): 766-781

42. Samper-González J, Burgos N, Bottani S, Fontanella S, Lu P, Marcoux A, Routier A, Guillon J, Bacci M, Wen J, et al (2018) Reproducible evaluation of classification methods in Alzheimer's disease: framework and application to MRI and PET data. Neuroimage 183:504-521

43. Karakuzu A, DuPre E, Tetrel L, Bermudez P, Boudreau M, Chin M, Poline JB, Das S, Bellec P, Stikov N (2022) Neurolibre: a preprint server for full-fledged reproducible neuroscience. OSF Preprints

44. Crane M (2018) Questionable answers in question answering research: reproducibility and variability of published results. Trans Assoc Comput Linguist 6:241-252

45. Jungo A, Scheidegger O, Reyes M, Balsiger F (2021) pymia: a python package for data handling and evaluation in deep learning-based medical image analysis. Comput Methods Program Biomed 198:105796. https://doi.org/10.1016/j.cmpb.2020.105796

46. Carlini N, Tramer F, Wallace E, Jagielski M, Herbert-Voss A, Lee K, Roberts A, Brown T, Song D, Erlingsson U, et al (2021) Extracting training data from large language models. In: 30th USENIX Security Symposium (USENIX Security 21), pp 2633-2650
47. Abadi M, Chu A, Goodfellow I, McMahan HB, Mironov I, Talwar K, Zhang L (2016) Deep learning with differential privacy. In: Proceedings of the 2016 ACM SIGSAC conference on computer and communications security, pp 308–318

48. Gorgolewski KJ, Alfaro-Almagro F, Auer T, Bellec P, Capotà M, Chakravarty MM, Churchill NW, Cohen AL, Craddock RC, Devenyi GA, et al (2017) BIDS apps: Improving ease of use, accessibility, and reproducibility of neuroimaging data analysis methods. PLoS Comput Biol 13(3):e1005209

49. Routier A, Burgos N, Díaz M, Bacci M, Bottani S, El-Rifai O, Fontanella S, Gori P, Guillon J, Guyot A, et al (2021) Clinica: an open-source software platform for reproducible clinical neuroscience studies. Front Neuroinform 15:e689675

50. McCormick M, Liu X, Jomier J, Marion C, Ibanez L (2014) ITK: enabling reproducible research and open science. Front Neuroinform 8:13

51. Yoo TS, Ackerman MJ, Lorensen WE, Schroeder W, Chalana V, Aylward S, Metaxas D, Whitaker R (2002) Engineering and algorithm design for an image processing API: a technical report on ITK-the insight toolkit. In: Medicine Meets Virtual Reality 02/10, IOS press, pp 586–592

52. Drummond C (2009) Replicability is not reproducibility: nor is it good science. In: Proceedings of the Evaluation Methods for Machine Learning Workshop at the 26th ICML, vol 1

53. Bouthillier X, Delaunay P, Bronzi M, Trofimov A, Nichyporuk B, Szeto J, Mohammadi Sepahvand N, Raff E, Madan K, Voleti V, et al (2021) Accounting for variance in machine learning benchmarks. Proc Mach Learn Syst 3:747–769

54. Wen J, Thibeau-Sutre E, Diaz-Melo M, Samper-González J, Routier A, Bottani S, Dormont D, Durrleman S, Burgos N, Colliot O (2020) Convolutional neural networks for classification of Alzheimer's disease: overview and reproducible evaluation. Med Image Anal 63:101694

55. Samala RK, Chan HP, Hadjiiski L, Koneru S (2020) Hazards of data leakage in machine learning: a study on classification of breast cancer using deep neural networks. In: Proceedings of SPIE Medical Imaging 2020: Computer-Aided Diagnosis, International Society for Optics and Photonics, vol 11314, p 1131416

56. Panwar H, Gupta PK, Siddiqui MK, Morales-Menendez R, Singh V (2020) Application of

deep learning for fast detection of COVID-19 in X-Rays using nCOVnet. Chaos, Solitons Fractals 138:109944

57. Bussola N, Marcolini A, Maggio V, Jurman G, Furlanello C (2021) AI slipping on tiles: Data leakage in Digital Pathology. In: Del Bimbo A, Cucchiara R, Sclaroff S, Farinella GM, Mei T, Bertini M, Escalante HJ, Vezzani R (eds) Pattern recognition. ICPR International Workshops and Challenges, Springer International Publishing, Cham. Lecture notes in computer science, pp 167–182. https://doi.org/10.1007/978-3-030-68763-2_13

58. Head ML, Holman L, Lanfear R, Kahn AT, Jennions MD (2015) The extent and consequences of p-hacking in science. PLoS Biol 13(3):e1002106

59. Henderson EL (2022) A guide to preregistration and registered reports. Preprint. https://osf.io/preprints/metaarxiv/x7aqr/download

60. Bottani S, Burgos N, Maire A, Wild A, Ströer S, Dormont D, Colliot O, Group AS et al (2022) Automatic quality control of brain T1-weighted magnetic resonance images for a clinical data warehouse. Med Image Anal 75: 102219

61. Perkuhn M, Stavrinou P, Thiele F, Shakirin G, Mohan M, Garmpis D, Kabbasch C, Borggrefe J (2018) Clinical evaluation of a multiparametric deep learning model for glioblastoma segmentation using heterogeneous magnetic resonance imaging data from clinical routine. Investig Radiol 53(11):647

62. Lukas C, Hahn HK, Bellenberg B, Rexilius J, Schmid G, Schimrigk SK, Przuntek H, Köster O, Peitgen HO (2004) Sensitivity and reproducibility of a new fast 3D segmentation technique for clinical MR-based brain volumetry in multiple sclerosis. Neuroradiology 46(11):906–915

63. Borga M, Ahlgren A, Romu T, Widholm P, Dahlqvist Leinhard O, West J (2020) Reproducibility and repeatability of MRI-based body composition analysis. Magn Reson Med 84(6): 3146–3156

64. Chard DT, Parker GJ, Griffin CM, Thompson AJ, Miller DH (2002) The reproducibility and sensitivity of brain tissue volume measurements derived from an SPM-based segmentation methodology. J Magn Reson Imaging: Off J Int Soc Magn Reson Med 15(3):259–267

65. de Boer R, Vrooman HA, Ikram MA, Vernooij MW, Breteler MM, van der Lugt A, Niessen WJ (2010) Accuracy and reproducibility study of automatic MRI brain tissue segmentation methods. Neuroimage 51(3):1047–1056
66. Lemieux L, Hagemann G, Krakow K, Woer-
mann FG (1999) Fast, accurate, and reproduc-
ible automatic segmentation of the brain in
T1-weighted volume MRI data. Magn
Reson Med: Off J Int Soc Magn Reson Med
42(1):127-135

67. Tudorascu DL, Karim HT, Maronge JM, Alhilali L, Fakhran S, Aizenstein HJ, Muschelli J, Crainiceanu CM (2016) Reproducibility and bias in healthy brain segmentation: comparison of two popular neuroimaging platforms. Front Neurosci 10:503

68. Yamashita R, Perrin T, Chakraborty J, Chou JF, Horvat N, Koszalka MA, Midya A, Gonen M, Allen P, Jarnagin WR et al (2020) Radiomic feature reproducibility in contrast-enhanced CT of the pancreas is affected by variabilities in scan parameters and manual segmentation. Euro Radiol 30(1):195-205

69. Poldrack RA, Whitaker K, Kennedy DN (2019) Introduction to the special issue on reproducibility in neuroimaging. Neuroimage 218:116357

70. Palumbo L, Bosco P, Fantacci M, Ferrari E, Oliva P, Spera G, Retico A (2019) Evaluation of the intra-and inter-method agreement of brain MRI segmentation software packages: a comparison between SPM12 and FreeSurfer v6.0. Phys Med 64:261-272

71. Laurienti PJ, Field AS, Burdette JH, Maldjian JA, Yen YF, Moody DM (2002) Dietary caffeine consumption modulates fMRI measures. Neuroimage 17(2):751-757

72. Collins DL, Zijdenbos AP, Kollokian V, Sled JG, Kabani NJ, Holmes CJ, Evans AC (1998) Design and construction of a realistic digital

brain phantom. IEEE Trans Med Imaging
17(3):463-468

73. Shaw R, Sudre C, Ourselin S, Cardoso MJ (2018) MRI K-space motion artefact augmentation: Model robustness and task-specific uncertainty. In: Medical Imaging with Deep Learning – MIDL 2018

74. Duffy BA, Zhang W, Tang H, Zhao L (2018) Retrospective correction of motion artifact affected structural MRI images using deep learning of simulated motion. In: Medical Imaging with Deep Learning – MIDL 2018

75. Loizillon S, Bottani S, Maire A, Ströer S, Dormont D, Colliot O, Burgos N (2023) Transfer learning from synthetic to routine clinical data for motion artefact detection in brain t1-weighted MRI. In: SPIE Medical Imaging 2023: Image Processing

76. Peng RD (2011) Reproducible research in computational science. Science 334(6060): 1226-1227

77. Piwowar HA, Day RS, Fridsma DB (2007) Sharing detailed research data is associated with increased citation rate. PLoS One 2(3): e308

78. Piwowar HA, Vision TJ (2013) Data reuse and the open data citation advantage. PeerJ 1:e175

79. Vandewalle P (2012) Code sharing is associated with research impact in image processing. Comput Sci Eng 14(4):42-47

80. Balsiger F, Jungo A, Chen J, Ezhov I, Liu S, Ma J, Paetzold JC, Sekuboyina A, Shit S, Suter Y et al (2021) The miccai hackathon on reproducibility, diversity, and selection of papers at the miccai conference. arXiv preprint arXiv:210305437
# Interpretability of Machine Learning Methods Applied to Neuroimaging

Elina Thibeau-Sutre, Sasha Collin, Ninon Burgos, and Olivier Colliot

## Abstract

Deep learning methods have become very popular for the processing of natural images and were then successfully adapted to the neuroimaging field. As these methods are non-transparent, interpretability methods are needed to validate them and ensure their reliability. Indeed, it has been shown that deep learning models may obtain high performance even when using irrelevant features, by exploiting biases in the training set. Such undesirable situations can potentially be detected by using interpretability methods. Recently, many methods have been proposed to interpret neural networks. However, this domain is not mature yet. Machine learning users face two major issues when aiming to interpret their models: which method to choose and how to assess its reliability. Here, we aim at providing answers to these questions by presenting the most common interpretability methods and metrics developed to assess their reliability, as well as their applications and benchmarks in the neuroimaging context. Note that this is not an exhaustive survey: we aimed to focus on the studies which we found to be the most representative and relevant.

**Key words** Interpretability, Saliency, Machine learning, Deep learning, Neuroimaging, Brain disorders

# 1 Introduction

## 1.1 Need for Interpretability

Many metrics have been developed to evaluate the performance of machine learning (ML) systems. In the case of supervised systems, these metrics compare the output of the algorithm to a ground truth, in order to evaluate its ability to reproduce a label given by a physician. However, the users (patients and clinicians) may want more information before relying on such systems. On which features is the model relying to compute the results? Are these features close to the way a clinician thinks? If not, why? This questioning coming from the actors of the medical field is justified, as errors in real life may lead to dramatic consequences. Trust into ML systems cannot be built only based on a set of metrics evaluating the performance of the system. Indeed, various examples of machine learning systems taking correct decisions for the wrong reasons
<figure><img src="image_234.png" /><figcaption>**Fig. 1** Example of an interpretability method highlighting why a network took the wrong decision. The explained classifier was trained on the binary task “Husky” vs “Wolf.” The pixels used by the model are actually in the background and highlight the snow. (Adapted from [1]. Permission to reuse was kindly granted by the authors)</figcaption></figure>

exist, e.g., [1–3]. Thus, even though their performance is high, they may be unreliable and, for instance, not generalize well to slightly different data sets. One can try to prevent this issue by interpreting the model with an appropriate method whose output will highlight the reasons why a model took its decision.

In [1], the authors show a now classical case of a system that correctly classifies images for wrong reasons. They purposely designed a biased data set in which wolves always are in a snowy environment whereas huskies are not. Then, they trained a classifier to differentiate wolves from huskies: this classifier had good accuracy but classified wolves as huskies when the background was snowy and huskies as wolves when there was no snow. Using an interpretability method, they further highlighted that the classifier was looking at the background and not at the animal (see Fig. 1).

Another study [2] detected a bias in ImageNet (a widely used data set of natural images) as the interpretation of images with the label “chocolate sauce” highlighted the importance of the spoon. Indeed, ImageNet “chocolate sauce” images often contained spoons, leading to a spurious correlation. There are also examples of similar problems in medical applications. For instance, a recent paper [3] showed with interpretability methods that some deep learning systems detecting COVID-19 from chest radiographs actually relied on confounding factors rather than on the actual pathological features. Indeed, their model focused on other regions than the lungs to evaluate the COVID-19 status (edges, diaphragm, and cardiac silhouette). Of note, their model was trained on public data sets which were used by many studies.
## 1.2 How to Interpret Models

According to [4], model interpretability can be broken down into two categories: transparency and post hoc explanations.

A model can be considered as transparent when it (or all parts of it) can be fully understood as such, or when the learning process is understandable. A natural and common candidate that fits, at first sight, these criteria is the linear regression algorithm, where coefficients are usually seen as the individual contributions of the input features. Another candidate is the decision tree approach where model predictions can be broken down into a series of understandable operations. One can reasonably consider these models as transparent: one can easily identify the features that were used to take the decision. However, one may need to be cautious not to push too far the medical interpretation. Indeed, the fact that a feature has not been used by the model does not mean that it is not associated with the target. It just means that the model did not need it to increase its performance. For instance, a classifier aiming at diagnosing Alzheimer's disease may need only a set of regions (for instance, from the medial temporal lobe of the brain) to achieve an optimal performance. This does not mean that other brain regions are not affected by the disease, just that they were not used by the model to take its decision. This is the case, for example, for sparse models like LASSO, but also standard multiple linear regressions. Moreover, features given as input to transparent models are often highly engineered, and choices made before the training step (preprocessing, feature selection) may also hurt the transparency of the whole framework. Nevertheless, in spite of these caveats, such models can reasonably be considered transparent, in particular when compared to deep neural networks which are intrinsically black boxes.

The second category of interpretability methods, post hoc interpretations, allows dealing with non-transparent models. Xie et al. [5] proposed a taxonomy in three categories: *visualization* methods consist in extracting an attribution map of the same size as the input whose intensities allow knowing where the algorithm focused its attention, *distillation* approaches consist in reproducing the behavior of a black box model with a transparent one, and *intrinsic* strategies include interpretability components within the framework, which are trained along with the main task (e.g., a classification). In the present work, we focus on this second category of methods (post hoc) and proposed a new taxonomy including other methods of interpretation (see Fig. 2). Post hoc interpretability is the most used category nowadays, as it allows interpreting deep learning methods that became the state of the art for many tasks in neuroimaging, as in other application fields.
<figure><img src="image_235.png" /><figcaption>Fig. 2 Taxonomy of the main interpretability methods</figcaption></figure>
## 1.3 Chapter Content and Outline

This chapter focuses on methods developed to interpret non-transparent machine learning systems, mainly deep learning systems, computing classification, or regression tasks from high-dimensional inputs. The interpretability of other frameworks (in particular generative models such as variational autoencoders or generative adversarial networks) is not covered as there are not enough studies addressing them. It may be because high-dimensional outputs (such as images) are easier to interpret “as such,” whereas small dimensional outputs (such as scalars) are less transparent.

Most interpretability methods presented in this chapter produce an attribution map: an array with the same dimensions as that of the input (up to a resizing) that can be overlaid on top of the input in order to exhibit an explanation of the model prediction. In the literature, many different terms may coexist to name this output such as saliency map, interpretation map, or heatmap. To avoid misunderstandings, in the following, we will only use the term “attribution map.”

The chapter is organized as follows. Subheading 2 presents the most commonly used interpretability methods proposed for computer vision, independently of medical applications. It also describes metrics developed to evaluate the reliability of interpretability methods. Then, Subheading 3 details their application to neuroimaging. Finally, Subheading 4 discusses current limitations of interpretability methods, presents benchmarks conducted in the neuroimaging field, and gives some advice to the readers who would like to interpret their own models.

Mathematical notations and abbreviations used during this chapter are summarized in Tables 1 and 2. A short reminder on neural network training procedure and a brief description of the diseases mentioned in the present chapter are provided in Appendices A and B.

## 2 Interpretability Methods

This section presents the main interpretability methods proposed in the domain of computer vision. We restrict ourselves to the methods that have been applied to the neuroimaging domain (the applications themselves being presented in Subheading 3). The outline of this section is largely inspired from the one proposed by Xie et al. [5]:

1. **Weight visualization** consists in directly visualizing weights learned by the model, which is natural for linear models but quite less informative for deep learning networks.
2. **Feature map visualization** consists in displaying intermediate results produced by a deep learning network to better understand its operation principle.

3. **Back-propagation methods** back-propagate a signal through the machine learning system from the output node of interest $o_c$ to the level of the input to produce an attribution map.

<table><caption>**Table 1**

**Mathematical notations**</caption><tr><td>$X_0$ is the input tensor given to the network, and $X$ refers to any input, sampled from the set $X$.</td></tr><tr><td>$y$ is a vector of target classes corresponding to the input.</td></tr><tr><td>$f$ is a network of $L$ layers. The first layer is the closest to the input; the last layer is the closest to the output. A layer is a function.</td></tr><tr><td>$g$ is a transparent function which aims at reproducing the behavior of $f$.</td></tr><tr><td>$w$ and $b$ are the weights and the bias associated to a linear function (e.g., in a fully connected layer).</td></tr><tr><td>$u$ and $v$ are locations (set of coordinates) corresponding to a node in a feature map. They belong respectively to the set $U$ and $V$.</td></tr><tr><td>$A_k^{(l)}(u)$ is the value of the feature map computed by layer $l$, of $K$ channels at channel $k$, at position $u$.</td></tr><tr><td>$R_k^{(l)}(u)$ is the value of a property back-propagated through the $l+1$, of $K$ channels at channel $k$, at position $u$. $R^{(l)}$ and $A^{(l)}$ have the same number of channels.</td></tr><tr><td>$o_c$ is the output node of interest (in a classification framework, it corresponds to the node of the class $c$).</td></tr><tr><td>$S_c$ is an attribution map corresponding to the output node $o_c$.</td></tr><tr><td>$m$ is a mask of perturbations. It can be applied to $X$ to compute its perturbed version $X^m$.</td></tr><tr><td>$\Phi$ is a function producing a perturbed version of an input $X$.</td></tr><tr><td>$\Gamma_c$ is the function computing the attribution map $S_c$ from the black-box function $f$ and an input $X_0$.</td></tr></table>

<table><caption>**Table 2**

**Abbreviations**</caption><tr><td>CAM</td><td>Class activation maps</td></tr><tr><td>CNN</td><td>Convolutional neural network</td></tr><tr><td>CT</td><td>Computed tomography</td></tr><tr><td>Grad-CAM</td><td>Gradient-weighted class activation mapping</td></tr><tr><td>LIME</td><td>Local interpretable model-agnostic explanations</td></tr><tr><td>LRP</td><td>Layer-wise relevance</td></tr><tr><td>MRI</td><td>Magnetic resonance imaging</td></tr><tr><td>SHAP</td><td>SHapley Additive exPlanations</td></tr><tr><td>T1w</td><td>T1-weighted [Magnetic Resonance Imaging]</td></tr></table>
4. Perturbation methods locally perturb the input and evaluate the difference in performance between using the original input and the perturbed version to infer which parts of the input are relevant for the machine learning system.

5. Distillation approximates the behavior of a black box model with a more transparent one and then draw conclusions from this new model.

6. Intrinsic includes the only methods of this chapter that are not post hoc explanations: in this case, interpretability is obtained thanks to components of the framework that are trained at the same time as the model.

Finally, for the methods producing an attribution map, a section is dedicated to the metrics used to evaluate different properties (e.g., reliability or human intelligibility) of the maps.

We caution readers that this taxonomy is not perfect: some methods may belong to several categories (e.g., LIME and SHAP could belong either to perturbation or distillation methods). Moreover, interpretability is still an active research field, and then some categories may (dis)appear or be fused in the future.

The interpretability methods were (most of the time) originally proposed in the context of a classification task. In this case, the network outputs an array of size C, corresponding to the number of different labels existing in the data set, and the goal is to know how the output node corresponding to a particular class *c* interacts with the input or with other parts of the network. However, these techniques can be extended to other tasks: for example, for a regression task, we will just have to consider the output node containing the continuous variable learned by the network. Moreover, some methods do not depend on the nature of the algorithm (e.g., standard perturbation or LIME) and can be applied to any machine learning algorithm.

## 2.1 Weight Visualization

At first sight, one of can be tempted to directly visualize the weights learned by the algorithm. This method is really simple, as it does not require further processing. However, even though it can make sense for linear models, it is not very informative for most networks unless they are specially designed for this interpretation.

This is the case for AlexNet [7], a convolutional neural network (CNN) trained on natural images (ImageNet). In this network the size of the kernels in the first layer is large enough (11 × 11) to distinguish patterns of interest. Moreover, as the three channels in the first layer correspond to the three color channels of the images (red, green, and blue), the values of the kernels can also be represented in terms of colors (this is not the case for hidden layers, in which the meaning of the channels is lost). The 96 kernels of the first layer were illustrated in the original article as in Fig. 3.
<figure><img src="image_236.png" /><figcaption>**Fig. 3** 96 convolutional kernels of size 3@11 × 11 learned by the first convolutional layer on the 3@224 × 224 input images by AlexNet. (Adapted from [7]. Permission to reuse was kindly granted by the authors)</figcaption></figure>

<figure><img src="image_237.png" /></figure>

<figure><img src="image_238.png" /><figcaption>**Fig. 4** The weights of small kernels in hidden layers (here 5 × 5) can be really difficult to interpret alone. Here some context allows better understanding how it modulates the interaction between concepts conveyed by the input and the output. (Adapted from [8] (CC BY 4.0))</figcaption></figure>

However, for hidden layers, this kind of interpretation may be misleading as nonlinearity activation layers are added between the convolutions and fully connected layers; this is why they only visualized the weights of the first layer.

To understand the weight visualization in hidden layers of a network, Voss et al. [8] proposed to add some context to the input and the output channels. This way they enriched the weight visualization with feature visualization methods able to generate an image corresponding to the input node and the output node (see Fig. 4). However, the feature visualization methods used to bring some context can also be difficult to interpret themselves, and then it only moves the interpretability problem from weights to features.

## 2.2 Feature Map Visualization

Feature maps are the results of intermediate computations done from the input and resulting in the output value. Then, it seems natural to visualize them or link them to concepts to understand how the input is successively transformed into the output.

Methods described in this section aim at highlighting which concepts a feature map (or part of it) A conveys.
### 2.2.1 Direct Interpretation

The output of a convolution has the same shape as its input: a 2D image processed by a convolution will become another 2D image (the size may vary). Then, it is possible to directly visualize these feature maps and compare them to the input to understand the operations performed by the network. However, the number of filters of convolutional layers (often a hundred) makes the interpretation difficult as a high number of images must be interpreted for a single input.

Instead of directly visualizing the feature map *A*, it is possible to study the latent space including all the values of the samples of a data set at the level of the feature map *A*. Then, it is possible to study the deformations of the input by drawing trajectories between samples in this latent space, or more simply to look at the distribution of some label in a manifold learned from the latent space. In such a way, it is possible to better understand which patterns were detected, or at which layer in the network classes begin to be separated (in the classification case). There is often no theoretical framework to illustrate these techniques, and then we referred to studies in the context of the medical application (see Subheading 3.2 for references).

### 2.2.2 Input Optimization

Olah et al. [9] proposed to compute an input that maximizes the value of a feature map *A* (see Fig. 5). However, this technique leads to unrealistic images that may be themselves difficult to interpret, particularly for neuroimaging data. To have a better insight of the behavior of layers or filters, another simple technique illustrated by the same authors consists in isolating the inputs that led to the highest activation of *A*. The combination of both methods, displayed in Fig. 6, allows a better understanding of the concepts conveyed by *A* of a GoogLeNet trained on natural images.

## 2.3 Back-Propagation Methods

The goal of these interpretability methods is to link the value of an output node of interest $o_c$ to the image $X_0$ given as input to a network. They do so by back-propagating a signal from $o_c$ to $X_0$:

<figure><img src="image_239.png" /><figcaption>**Fig. 5** Optimization of the input for different levels of feature maps. (Adapted from [9] (CC BY 4.0))</figcaption></figure>
**Dataset Examples** show us what neurons respond to in practice

<figure><img src="image_240.png" /></figure>

**Optimization** isolates the causes of behavior from mere correlations. A neuron may not be detecting what you initially thought.

<figure><img src="image_241.png" /><figcaption>**Fig. 7** Attribution map of an image found with gradients back-propagation. (Adapted from [10]. Permission to reuse was kindly granted by the authors)</figcaption></figure>

this process (backward pass) can be seen as the opposite operation than the one done when computing the output value from the input (forward pass).

Any property can be back-propagated as soon as its value at the level of a feature map *l*−1 can be computed from its value in the feature map *l*. In this section, the back-propagated properties are gradients or the relevance of a node *o*<sub>*c*</sub>.

### 2.3.1 Gradient Back-Propagation

During network training, gradients corresponding to each layer are computed according to the loss to update the weights. Then, we can see these gradients as the difference needed at the layer level to improve the final result: by adding this difference to the weights, the probability of the true class *y* increases.

In the same way, the gradients can be computed at the image level to find how the input should vary to change the value of *o*<sub>*c*</sub> (see example on Fig. 7. This gradient computation was proposed by [10], in which the attribution map *S*<sub>*c*</sub> corresponding to the input
image $X_0$ and the output node $o_c$ is computed according to the following equation:

$$S_c = \frac{\partial o_c}{\partial X} \bigg|_{X=X_0} \qquad (1)$$

Due to its simplicity, this method is the most commonly used to interpret deep learning networks. Its attribution map is often called a “saliency map”; however, this term is also used in some articles to talk about any attribution map, and this is why we chose to avoid this term in this chapter.

This method was modified to derive many similar methods based on gradient computation described in the following paragraphs.

**Gradient⊙Input** This method is the point-wise product of the gradient map described at the beginning of the section and the input. Evaluated in [11], it was presented as an improvement of the gradients method, though the original paper does not give strong arguments on the nature of this improvement.

**DeconvNet & Guided Back-Propagation** The key difference between this procedure and the standard back-propagation method is the way the gradients are back-propagated through the ReLU layer.

The ReLU layer is a commonly used activation function that sets to 0 the negative input values and does not affect positive input values. The derivative of this function in layer $l$ is the indicator function $\mathbb{1}_{A^{(l)} > 0}$: it outputs 1 (resp. 0) where the feature maps computed during the forward pass were positive (resp. negative).

Springenberg et al. [12] proposed to back-propagate the signal differently. Instead of applying the indicator function of the feature map $A^{(l)}$ computed during the forward pass, they directly applied ReLU to the back-propagated values $R^{(l+1)} = \frac{\partial o_c}{\partial A^{(l+1)}}$, which corresponds to multiplying it by the indicator function $\mathbb{1}_{R^{(l+1)} > 0}$. This “backward deconvnet” method allows back-propagating only the positive gradients, and, according to the authors, it results in a reconstructed image showing the part of the input image that is most strongly activating this neuron.

The guided back-propagation method (Eq. 4) combines the standard back-propagation (Eq. 2) with the backward deconvnet (Eq. 3): when back-propagating gradients through ReLU layers, a value is set to 0 if the corresponding top gradients or bottom data is negative. This adds an additional guidance to the standard back-propagation by preventing backward flow of negative gradients.
$$R^{(l)} = \mathbb{1}_{A^{(l)} > 0} * R^{(l+1)} \quad (2)$$

$$R^{(l)} = \mathbb{1}_{R^{(l+1)} > 0} * R^{(l+1)} \quad (3)$$

$$R^{(l)} = \mathbb{1}_{A^{(l)} > 0} * \mathbb{1}_{R^{(l+1)} > 0} * R^{(l+1)} \quad (4)$$

Any back-propagation procedure can be “guided,” as it only concerns the way ReLU functions are managed during back-propagation (this is the case, e.g., for guided Grad-CAM).

While it was initially adopted by the community, this method showed severe defects as discussed later in Subheading 4.

**CAM & Grad-CAM** In this setting, attribution maps are computed at the level of a feature map produced by a convolutional layer and then upsampled to be overlaid and compared with the input. The first method, class activation maps (CAM), was proposed by Zhou et al. [13] and can be only applied to CNNs with the following specific architecture:

1. A series of convolutions associated with activation functions and possibly pooling layers. These convolutions output a feature map A with N channels.

2. A global average pooling that extracts the mean value of each channel of the feature map produced by the convolutions.

3. A single fully connected layer

The CAM corresponding to $o_c$ will be the mean of the channels of the feature map produced by the convolutions, weighted by the weights $w_{kc}$ learned in the fully connected layer

$$S_c = \sum_{k=1}^{N} w_{kc} * A_k. \quad (5)$$

This map has the same size as $A_k$, which might be smaller than the input if the convolutional part performs downsampling operations (which is very often the case). Then, the map is upsampled to the size of the input to be overlaid on the input.

Selvaraju et al. [14] proposed an extension of CAM that can be applied to any architecture: Grad-CAM (illustrated on Fig. 8). As in CAM, the attribution map is a linear combination of the channels of a feature map computed by a convolutional layer. But, in this case, the weights of each channel are computed using gradient back-propagation

$$\alpha_{kc} = \frac{1}{|\mathcal{U}|} \sum_{u \in \mathcal{U}} \frac{\partial o_c}{\partial A_k(u)}. \quad (6)$$
<figure><img src="image_242.png" /><figcaption>**Fig. 8** Grad-CAM explanations highlighting two different objects in an image. (a) the original image, (b) the explanation based on the “dog” node, (c) the explanation based on the “cat” node. ©2017 IEEE. (Reprinted, with permission, from [14])</figcaption></figure>

The final map is then the linear combination of the feature maps weighted by the coefficients. A ReLU activation is then applied to the result to only keep the features that have a positive influence on class c

$$S_c = ReLU\left(\sum_{k=1}^{N} \alpha_{kc} * A_k\right). \qquad (7)$$

Similarly to CAM, this map is then upsampled to the input size.

Grad-CAM can be applied to any feature map produced by a convolution, but in practice the last convolutional layer is very often chosen. The authors argue that this layer is “the best compromise between high-level semantics and detailed spatial information” (the latter is lost in fully connected layers, as the feature maps are flattened).

Because of the upsampling step, CAM and Grad-CAM produce maps that are more human-friendly because they contain more connected zones, contrary to other attribution maps obtained with gradient back-propagation that can look very scattered. However, the smallest the feature maps $A_k$, the blurrier they are, leading to a possible loss of interpretability.

**2.3.2 Relevance Back-Propagation**

Instead of back-propagating gradients to the level of the input or of the last convolutional layer, Bach et al. [15] proposed to back-propagate the score obtained by a class c, which is called the relevance. This score corresponds to $o_c$ after some post-processing (e.g., softmax), as its value must be positive if class c was identified in the input. At the end of the back-propagation process, the goal is to find the relevance $R_u$ of each feature u of the input (e.g., of each pixel of an image) such that $o_c = \sum_{u \in U} R_u$.

In their paper, Bach et al. [15] take the example of a fully connected function defined by a matrix of weights w and a bias
<figure><img src="image_243.png" /><figcaption>**Fig. 9** LRP attribution maps explaining the decision of a neural network trained on MNIST. ©2017 IEEE. (Reprinted, with permission, from [16])</figcaption></figure>

b at layer $l+1$. The value of a node $v$ in feature map $A^{(l+1)}$ is computed during the forward pass by the given formula:

$$A^{(l+1)}(v) = b + \sum_{u \in u} w_{uv} A^{(l)}(u) \quad (8)$$

During the back-propagation of the relevance, $R^{(l)}(u)$, the value of the relevance at the level of the layer $l+1$ is computed according to the values of the relevance $R^{(l+1)}(v)$ which are distributed according to the weights $w$ learned during the forward pass and the values of $A^{(l)}(v)$:

$$R^{(l)}(u) = \sum_{v \in v} R^{(l+1)}(v) \frac{A^{(l)}(u) w_{uv}}{\sum_{u' \in u} A^{(l)}(u') w_{u'v}}. \quad (9)$$

The main issue of the method comes from the fact that the denominator may become (close to) zero, leading to the explosion of the relevance back-propagated. Moreover, it was shown by [11] that when all activations are piece-wise linear (such as ReLU or leaky ReLU), the layer-wise relevance (LRP) method reproduces the output of gradient input, questioning the usefulness of the method.

This is why Samek et al. [16] proposed two variants of the standard LRP method [15]. Moreover they describe the behavior of the back-propagation in other layers than the linear ones (the convolutional one following the same formula as the linear). They illustrated their method with a neural network trained on MNIST (see Fig. 9). To simplify the equations in the following paragraphs, we now denote the weighted activations as $z_{uv} = A^{(l)}(u)w_{uv}$.

**ε-rule** The $\epsilon$-rule integrates a parameter $\epsilon > 0$, used to avoid numerical instability. Though it avoids the case of a null denominator, this variant breaks the rule of relevance conservation across layers
$$R^{(l)}(u) = \sum_{v \in \mathcal{V}'} R^{(l+1)}(v) \frac{z_{uv}}{\sum_{u' \in \mathcal{U}} z_{u'v} + \varepsilon \times \text{sign}\left(\sum_{u' \in \mathcal{U}} z_{u'v}\right)}. \quad (10)$$

**β-rule** The β-rule keeps the conservation of the relevance by treating separately the positive weighted activations $z_{uv}^{+}$ from the negative ones $z_{uv}^{-}$

$$R^{(l)}(u) = \sum_{v \in \mathcal{V}'} R^{(l+1)}(v) \left( (1 + \beta) \frac{z_{uv}^{+}}{\sum_{u' \in \mathcal{U}} z_{u'v}^{+}} - \beta \frac{z_{uv}^{-}}{\sum_{u' \in \mathcal{U}} z_{u'v}^{-}} \right). \quad (11)$$

Though these two LRP variants improve the numerical stability of the procedure, they imply to choose the values of parameters that may change the patterns in the obtained attribution map.

**Deep Taylor Decomposition** Deep Taylor decomposition [17] was proposed by the same team as the one that proposed the original LRP method and its variants. It is based on similar principles as LRP: the value of the score obtained by a class c is back-propagated, but the back-propagation rule is based on first-order Taylor expansions.

The back-propagation from node v in at the level of $R^{(l+1)}$ to u at the level of $R^{(l)}$ can be written

$$R^{(l)}(u) = \sum_{v \in \mathcal{V}'} \left. \frac{\partial R^{(l+1)}(v)}{\partial A^{(l)}(u)} \right|_{\tilde{A}_{v}^{(l)}(u)} \left( A^{(l)}(u) - \tilde{A}_{v}^{(l)}(u) \right). \quad (12)$$

This rule implies a root point $\tilde{A}_{v}^{(l)}(u)$ which is close to $A^{(l)}(u)$ and meets a set of constraints depending on v.

## 2.4 Perturbation Methods

Instead of relying on a backward pass (from the output to the input) as in the previous section, perturbation methods rely on the difference between the value of $o_c$ computed with the original inputs and a locally perturbed input. This process is less abstract for humans than back-propagation methods as we can reproduce it ourselves: if the part of the image that is needed to find the good output is hidden, we are also not able to predict correctly. Moreover, it is model-agnostic and can be applied to any algorithm or deep learning architecture.

The main drawback of these techniques is that the nature of the perturbation is crucial, leading to different attribution maps depending on the perturbation function used. Moreover, Montavon et al. [18] suggest that the perturbation rule should
<figure><img src="image_244.png" /><figcaption>**Fig. 10** Attribution maps obtained with standard perturbation. Here the perturbation is a gray patch covering a specific zone of the input as shown in the left column. The attribution maps (second row) display the probability of the true label: the lower the value, the most important it is for the network to correctly identify the label. This kind of perturbation takes the perturbed input out of the training distribution. (Reprinted by permission from Springer Nature Customer Service Centre GmbH: Springer Nature, ECCV 2014: Visualizing and Understanding Convolutional Networks, [19], 2014)</figcaption></figure>

keep the perturbed input in the training data distribution. Indeed, if it is not the case, one cannot know if the network performance dropped because of the location or the nature of the perturbation.

### 2.4.1 Standard Perturbation

Zeiler and Fergus [19] proposed the most intuitive method relying on perturbations. This standard perturbation procedure consists in removing information locally in a specific zone of an input $X_0$ and evaluating if it modifies the output node $o_c$. The more the perturbation degrades the task performance, the more crucial this zone is for the network to correctly perform the task. To obtain the final attribution map, the input is perturbed according to all possible locations. Examples of attribution maps obtained with this method are displayed in Fig. 10.

As evaluating the impact of the perturbation at each pixel location is computationally expensive, one can choose not to perturb the image at each pixel location but to skip some of them (i.e., scan the image with a stride > 1). This will lead to a smaller attribution map, which needs to be upsampled to be compared to the original input (in the same way as CAM & Grad-CAM).
However, in addition to the problem of the nature of the perturbation previously mentioned, this method presents two drawbacks:

* The attribution maps depend on the size of the perturbation: if the perturbation becomes too large, the perturbation is not local anymore; if it is too small, it is not meaningful anymore (a pixel perturbation cannot cover a pattern).

* Input pixels are considered independently from each other: if the result of a network relies on a combination of pixels that cannot all be covered at the same time by the perturbation, their influence may not be detected.

### 2.4.2 Optimized Perturbation

To deal with these two issues, Fong and Vedaldi [2] proposed to optimize a perturbation mask covering the whole input. This perturbation mask $m$ has the same size as the input $X_0$. Its application is associated with a perturbation function $\Phi$ and leads to the computation of the perturbed input $X_0^m$. Its value at a coordinate $u$ reflects the quantity of information remaining in the perturbed image:

* If $m(u) = 1$, the pixel at location $u$ is not perturbed and has the same value in the perturbed input as in the original input ($X_0^m(u) = X_0(u)$).

* If $m(u) = 0$, the pixel at location $u$ is fully perturbed and the value in the perturbed image is the one given by the perturbation function only ($X_0^m(u) = \Phi(X_0)(u)$).

This principle can be extended to any value between 0 and 1 with the a linear interpolation

$$X_0^m(u) = m(u)X_0(u) + (1 - m(u))\Phi(X_0)(u). \quad (13)$$

Then, the goal is to optimize this mask $m$ according to three criteria:

1. The perturbed input $X_0^m$ should lead to the lowest performance possible.

2. The mask $m$ should perturb the minimum number of pixels possible.

3. The mask $m$ should produce connected zones (i.e., avoid the scattered aspect of gradient maps).

These three criteria are optimized using the following loss:

$$f(X_0^m) + \lambda_1 \|1 - m\|_{\beta_1}^{\beta_1} + \lambda_2 \|\nabla m\|_{\beta_2}^{\beta_2} \quad (14)$$

with $f$ a function that decreases as the performance of the network decreases.

However, the method also presents two drawbacks:
<figure><img src="image_245.png" /><figcaption>**Fig. 11** In this example, the network learned to classify objects in natural images. Instead of masking the maypole at the center of the image, it creates artifacts in the sky to degrade the performance of the network. ©2017 IEEE. (Reprinted, with permission, from [2])</figcaption></figure>

* The values of hyperparameters must be chosen ($\lambda_1, \lambda_2, \beta_1, \beta_2$) to find a balance between the three optimization criteria of the mask.

* The mask may not highlight the most important pixels of the input but instead create artifacts in the perturbed image to artificially degrade the performance of the network (see Fig. 11).

**2.5 Distillation**

Approaches described in this section aim at developing a transparent method to reproduce the behavior of a black box one. Then it is possible to consider simple interpretability methods (such as weight visualization) on the transparent method instead of considering the black box.

**2.5.1 Local Approximation**

**LIME** Ribeiro et al. [1] proposed local interpretable model-agnostic explanations (LIME). This approach is:

* **Local**, as the explanation is valid in the vicinity of a specific input $X_0$

* **Interpretable**, as an interpretable model $g$ (linear model, decision tree...) is computed to reproduce the behavior of $f$ on $X_0$

* **Model-agnostic**, as it does not depend on the algorithm trained

This last property comes from the fact that the vicinity of $X_0$ is explored by sampling variations of $X_0$ that are perturbed versions of $X_0$. Then LIME shares the advantage (model-agnostic) and draw-back (perturbation function dependent) of perturbation methods presented in Subheading 2.4. Moreover, the authors specify that, in the case of images, they group pixels of the input in $d$ super-pixels (contiguous patches of similar pixels).

The loss to be minimized to find $g$ specific to the input $X_0$ is the following:

$$ \mathcal{L}(f,g, \pi_{X_0}) + \Omega(g), \qquad (15) $$
where $\pi_{X_0}$ is a function that defines the locality of $X_0$ (i.e., $\pi_{X_0}(X)$ decreases as $X$ becomes closer to $X_0$), $\mathcal{L}$ measures how unfaithful $g$ is in approximating $f$ according $\pi_{X_0}$, and $\Omega$ is a measure of the complexity of $g$.

Ribeiro et al. [1] limited their search to sparse linear models; however, other assumptions could be made on $g$.

$g$ is not applied to the input directly but to a binary mask $m \in \{0, 1\}^d$ that transforms the input $X$ in $X^m$ and is applied according to a set of $d$ super-pixels. For each super-pixel $u$:

1. If $m(u) = 1$, the super-pixel $u$ is not perturbed.

2. If $m(u) = 0$, the super-pixel $u$ is perturbed (i.e., it is grayed).

They used

$$ \pi_{X_0}(X) = \exp \frac{(X - X_0)^2}{\sigma^2} $$

and

$$ \mathcal{L}(f, g, \pi_{X_0}) = \sum_m \pi_{X_0}(X_0^m) * (f(X_0^m) - g(m))^2. $$

Finally $\Omega(g)$ is the number of non-zero weights of $g$, and its value is limited to $K$. This way they select the $K$ super-pixels in $X_0$ that best explain the algorithm result $f(X_0)$.

**SHAP** Lundberg and Lee [20] proposed SHAP (SHapley Additive exPlanations), a theoretical framework that encompasses several existing interpretability methods, including LIME. In this framework each of the $N$ features (again, super-pixels for images) is associated with a coefficient $\phi$ that denotes its contribution to the result. The contribution of each feature is evaluated by perturbing the input $X_0$ with a binary mask $m$ (see paragraph on LIME). Then the goal is to find an interpretable model $g$ specific to $X_0$, such that

$$ g(m) = \phi_0 + \sum_{i=1}^{N} \phi_i m_i \quad (16) $$

with $\phi_0$ being the output when the input is fully perturbed.

The authors look for an expression of $\phi$ that respects three properties:

* **Local accuracy $g$ and $f$** should match in the vicinity of $X_0$: $g(m) = f(X_0^m)$.

* **Missingness** Perturbed features should not contribute to the result: $m_i = 0 \rightarrow \phi_i = 0$.

* **Consistency** Let's denote as $m \setminus i$ the mask $m$ in which $m_i = 0$. For any two models $f^1$ and $f^2$, if
<figure><img src="image_246.png" /><figcaption>Fig. 12 Visualization of a soft decision tree trained on MNIST. (Adapted from [21]. Permission to reuse was kindly granted by the authors)</figcaption></figure>

$$f^1(X_0^m) - f^1(X_0^{m\setminus i}) \geq f^2(X_0^m) - f^2(X_0^{m\setminus i})$$

then for all $m \in \{0, 1\}^N$ $\phi_i^1 \geq \phi_i^2$ ($\phi^k$ are the coefficients associated with model $f^k$).

Lundberg and Lee [20] show that only one expression is possible for the coefficients $\phi$, which can be approximated with different algorithms:

$$\phi_i = \sum_{m \in \{0, 1\}^N} \frac{|m|!(N - |m| - 1)!}{N!} [f(X_0^m) - f(X_0^{m\setminus i})]. \quad (17)$$

### 2.5.2 Model Translation

Contrary to local approximation, which provides an explanation according to a specific input $X_0$, model translation consists in finding a transparent model that reproduces the behavior of the black box model on the whole data set.

As it was rarely employed in neuroimaging frameworks, this section only discusses the distillation to decision trees proposed in [21] (preprint). For a more extensive review of model translation methods, we refer the reader to [5].

After training a machine learning system $f$, a binary decision tree $\mathcal{G}$ is trained to reproduce its behavior. This tree is trained on a set of inputs $X$, and each inner node $i$ learns a matrix of weights $w_i$ and biases $b_i$. The forward pass of $X$ in the node $i$ of the tree is as follows: if $\text{sigmoid}(w_iX + b_i) > 0.5$, then the right leaf node is chosen, else the left leaf node is chosen. After the end of the decision tree's training, it is possible to visualize at which level which classes were separated to better understand which classes are similar for the
network. It is also possible to visualize the matrices of weights learned by each inner node to identify patterns learned at each class separation. An illustration of this distillation process, on the MNIST data set (hand-written digits), can be found in Fig. 12.

## 2.6 Intrinsic

Contrary to the previous sections in which interpretability methods could be applied to (almost) any network after the end of the training procedure, the following methods require to design the framework before the training phase, as the interpretability components and the network are trained simultaneously. In the papers presented in this Subheading [22–24], the advantages of these methods are dual: they improve both the interpretability and performance of the network. However, the drawback is that they have to be implemented before training the network, and then they cannot be applied in all cases.

### 2.6.1 Attention Modules

Attention is a concept in machine learning that consists in producing an attribution map from a feature map and using it to improve learning of another task (such as classification, regression, reconstruction...) by making the algorithm focus on the part of the feature map highlighted by the attribution map.

In the deep learning domain, we take as reference [22], in which a network is trained to produce a descriptive caption of natural images. This network is composed of three parts:

1. A convolutional encoder that reduces the dimension of the input image to the size of the feature maps $A$

<figure><img src="image_247.png" /><figcaption>**Fig. 13** Examples of images correctly captioned by the network. The focus of the attribution map is highlighted in white and the associated word in the caption is underlined. (Adapted from [22]. Permission to reuse was kindly granted by the authors)</figcaption></figure>
2. An attention module that generates an attribution map $S_t$ from $A$ and the previous hidden state of the long short-term memory (LSTM) network

3. An LSTM decoder that computes the caption from its previous hidden state, the previous word generated, $A$ and $S_t$

As $S_t$ is of the same size as $A$ (smaller than the input), the result is then upsampled to be overlaid on the input image. As one attribution map is generated per word generated by the LSTM, it is possible to know where the network focused when generating each word of the caption (see Fig. 13). In this example, the attribution map is given to a LSTM, which uses it to generate a context vector $z_t$ by applying a function $\phi$ to $A$ and $S_t$.

More generally in CNNs, the point-wise product of the attribution map $S$ and the feature map $A$ is used to generate the refined feature map $A'$ which is given to the next layers of the network. Adding an attention module implies to make new choices for the architecture of the model: its location (on lower or higher feature maps) may impact the performance of the network. Moreover, it is possible to stack several attention modules along the network, as it was done in [23].

### 2.6.2 Modular Transparency

Contrary to the studies of the previous sections, the frameworks of these categories are composed of several networks (modules) that interact with each other. Each module is a black box, but the transparency of the function, or the nature of the interaction

<figure><img src="image_248.png" /><figcaption>**Fig. 14** Framework with modular transparency browsing an image to compute the output at the global scale. (Adapted from [24]. Permission to reuse was kindly granted by the authors)</figcaption></figure>
between them, allows understanding how the system works globally and extracting interpretability metrics from it.

A large variety of setups can be designed following this principle, and it is not possible to draw a more detailed general rule for this section. We will take the example described in [24], which was adapted to neuroimaging data (see Subheading 3.6), to illustrate this section, though it may not be representative of all the aspects of modular transparency.

Ba et al. [24] proposed a framework (illustrated in Fig. 14) to perform the analysis of an image in the same way as a human, by looking at successive relevant locations in the image. To perform this task, they assemble a set of networks that interact together:

* **Glimpse network** This network takes as input a patch of the input image and the location of its center to output a context vector that will be processed by the recurrent network. Then this vector conveys information on the main features in a patch and its location.

* **Recurrent network** This network takes as input the successive context vectors and update its hidden state that will be used to find the next location to look at and to perform the learned task at the global scale (in the original paper a classification of the whole input image).

* **Emission network** This network takes as input the current state of the recurrent network and outputs the next location to look at. This will allow computing the patch that will feed the glimpse network.

* **Context network** This network takes as input the whole input at the beginning of the task and outputs the first context vector to initialize the recurrent network.

* **Classification network** This network takes as input the current state of the recurrent network and outputs a prediction for the class label.

The global framework can be seen as interpretable as it is possible to review the successive processed locations.

**2.7 Interpretability Metrics**

To evaluate the reliability of the methods presented in the previous sections, one cannot only rely on qualitative evaluation. This is why interpretability metrics that evaluate attribution maps were proposed. These metrics may evaluate different properties of attribution maps.

* **Fidelity** evaluates if the zones highlighted by the map influence the decision of the network.

* **Sensitivity** evaluates how the attribution map changes according to small changes in the input $X_0$.
• **Continuity** evaluates if two close data points lead to similar attribution maps.

In the following, $\Gamma$ is an interpretability method computing an attribution map $S$ of the black box network $f$ and an input $X_0$.

2.7.1 (In)fidelity

Yeh et al. [25] proposed a measure of infidelity of $\Gamma$ based on perturbations applied according to a vector $m$ of the same shape as the attribution map $S$. The explanation is infidel if perturbations applied in zones highlighted by $S$ on $X_0$ lead to negligible changes in $f(X_0^m)$ or, on the contrary, if perturbations applied in zones not highlighted by $S$ on $X_0$ lead to significant changes in $f(X_0^m)$. The associated formula is

$$ \mathrm{INFD}(\Gamma, f, X_0) = \mathbb{E}_m \left[ \sum_i \sum_j m_{ij} \Gamma(f, X_0)_{ij} - (f(X_0) - f(X_0^m))^2 \right]. \quad (18) $$

2.7.2 Sensitivity

Yeh et al. [25] also gave a measure of sensitivity. As suggested by the definition, it relies on the construction of attribution maps according to inputs similar to $X_0$: $\tilde{X}_0$. As changes are small, sensitivity depends on a scalar $\epsilon$ set by the user, which corresponds to the maximum difference allowed between $X_0$ and $\tilde{X}_0$. Then sensitivity corresponds to the following formula:

$$ \mathrm{SENS}_{\max}(\Gamma, f, X_0, \epsilon) = \max_{\|\tilde{X}_0 - X_0\| \le \epsilon} \frac{\|\Gamma(f, \tilde{X}_0) - \Gamma(f, X_0)\|}{\|\tilde{X}_0 - X_0\|} \quad (19) $$

2.7.3 Continuity

Continuity is very similar to sensitivity, except that it compares different data points belonging to the input domain $X$, whereas sensitivity may generate similar inputs with a perturbation method. This measure was introduced in [18] and can be computed using the following formula:

$$ \mathrm{CONT}(\Gamma, f, X) = \max_{X_1, X_2 \in X \ \& \ X_1 \neq X_2} \frac{\|\Gamma(f, X_1) - \Gamma(f, X_2)\|_1}{\|X_1 - X_2\|_2}. \quad (20) $$

As these metrics rely on perturbation, they are also influenced by the nature of the perturbation and may lead to different results, which is a major issue (see Subheading 4). Other metrics were also proposed and depend on the task learned by the network: for example, in the case of a classification, statistical tests can be conducted between attribution maps of different classes to assess whether they differ according to the class they explain.
## 3 Application of Interpretability Methods to Neuroimaging Data

In this section, we provide a non-exhaustive review of applications of interpretability methods to neuroimaging data. In most cases, the focus of articles is prediction/classification rather than the interpretability method, which is just seen as a tool to analyze the results. Thus, authors do not usually motivate their choice of an interpretability method. Another key consideration here is the spatial registration of brain images, which enables having brain regions roughly at the same position between subjects. This technique is of paramount importance as attribution maps computed for registered images can then be averaged or used to automatically determine the most important brain areas, which would not be possible with unaligned images. All the studies presented in this section are summarized in Table 3.

<caption>**Table 3**</caption>

**Summary of the studies applying interpretability methods to neuroimaging data which are presented in Subheading 3**

<table><thead><tr><td>Study</td><td>Data set</td><td>Modality</td><td>Task</td><td>Interpretability method</td><td>Section</td></tr></thead><tbody><tr><td>Abrol et al. [28]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>FM visualization, Perturbation</td><td>3.2, 3.4</td></tr><tr><td>Bae et al. [32]</td><td>ADNI</td><td>sMRI</td><td>AD classification</td><td>Perturbation</td><td>3.4</td></tr><tr><td>Ball et al. [33]</td><td>PING</td><td>T1w</td><td>Age prediction</td><td>Weight visualization, SHAP</td><td>3.1, 3.5</td></tr><tr><td>Biffi et al. [29]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>FM visualization</td><td>3.2</td></tr><tr><td>Böhle et al. [34]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>LRP, Guided back-propagation</td><td>3.3</td></tr><tr><td>Burduja et al. [35]</td><td>RSNA</td><td>CT scan</td><td>Intracranial Hemorrhage detection</td><td>Grad-CAM</td><td>3.3</td></tr><tr><td>Cecotti and Gräser [26]</td><td>in-house</td><td>EEG</td><td>P300 signals detection</td><td>Weight visualization</td><td>3.1</td></tr><tr><td>Dyrba et al. [36]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>DeconvNet, Deep Taylor decomposition, Gradient ⊙ Input, LRP, Grad-CAM</td><td>3.3</td></tr><tr><td>Eitel and Ritter [37]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>Gradient ⊙ Input, Guided back-propagation, LRP, Perturbation</td><td>3.3, 3.4</td></tr></tbody></table>

(continued)
<table><caption>**Table 3**
(continued)</caption><thead><tr><th>Study</th><th>Data set</th><th>Modality</th><th>Task</th><th>Interpretability method</th><th>Section</th></tr></thead><tbody><tr><td>Eitel et al. [38]</td><td>ADNI, in-house</td><td>T1w</td><td>Multiple Sclerosis detection</td><td>Gradient ⊙ Input, LRP</td><td>3.3</td></tr><tr><td>Fu et al. [39]</td><td>CQ500, RSNA</td><td>CT scan</td><td>Detection of Critical Findings in Head CT scan</td><td>Attention mechanism</td><td>3.6</td></tr><tr><td>Gutiérrez-Becker and Wachinger [40]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>Perturbation</td><td>3.4</td></tr><tr><td>Hu et al. [41]</td><td>ADNI, NIFD</td><td>T1w</td><td>AD/CN/FTD classification</td><td>Guided back-propagation</td><td>3.3</td></tr><tr><td>Jin et al. [42]</td><td>ADNI, in-house</td><td>T1w</td><td>AD classification</td><td>Attention mechanism</td><td>3.6</td></tr><tr><td>Lee et al. [43]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>Modular transparency</td><td>3.6</td></tr><tr><td>Leming et al. [31]</td><td>OpenFMRI, ADNI, ABIDE, ABIDE II, ABCD, NDAR ICBM, UK Biobank, 1000FC</td><td>fMRI</td><td>Autism classification Sex classification Task vs rest classification</td><td>FM visualization, Grad-CAM</td><td>3.2, 3.3</td></tr><tr><td>Magesh et al. [44]</td><td>PPMI</td><td>SPECT</td><td>Parkinson's disease detection</td><td>LIME</td><td>3.5</td></tr><tr><td>Martinez-Murcia et al. [30]</td><td>ADNI</td><td>T1w</td><td>AD classification Prediction of neuropsychological tests & other clinical variables</td><td>FM visualization</td><td>3.2</td></tr><tr><td>Nigri et al. [45]</td><td>ADNI, AIBL</td><td>T1w</td><td>AD classification</td><td>Perturbation, Swap test</td><td>3.4</td></tr><tr><td>Oh et al. [27]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>FM visualization, Standard back-propagation, Perturbation</td><td>3.2, 3.3, 3.4</td></tr><tr><td>Qiu et al. [46]</td><td>ADNI, AIBL, FHS, NACC</td><td>T1w</td><td>AD classification</td><td>Modular transparency</td><td>3.6</td></tr><tr><td>Ravi et al. [47]</td><td>ADNI</td><td>T1w</td><td>CN/MCI/AD reconstruction</td><td>Modular transparency</td><td>3.6</td></tr></tbody></table>

(continued)
<table><caption>**Table 3**
(continued)</caption><thead><tr><td>Study</td><td>Data set</td><td>Modality</td><td>Task</td><td>Interpretability method</td><td>Section</td></tr></thead><tbody><tr><td>Rieke et al. [48]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>Standard back-propagation, Guided back-propagation, Perturbation, Brain area occlusion</td><td>3.3, 3.4</td></tr><tr><td>Tang et al. [49]</td><td>UCD-ADC, Brain Bank</td><td>Histology</td><td>Detection of amyloid-β pathology</td><td>Guided back-propagation, Perturbation</td><td>3.3, 3.4</td></tr><tr><td>Wood et al. [50]</td><td>ADNI</td><td>T1w</td><td>AD classification</td><td>Modular transparency</td><td>3.6</td></tr></tbody></table>

**Data sets:** 1000FC, 1000 Functional Connectomes; ABCD, Adolescent Brain Cognitive Development; ABIDE, Autism Brain Imaging Data Exchange; ADNI, Alzheimer's Disease Neuroimaging Initiative; AIBL, Australian Imaging, Biomarkers and Lifestyle; FHS, Framingham Heart Study; ICBM, International Consortium for Brain Mapping; NACC, National Alzheimer's Coordinating Center; NDAR, National Database for Autism Research; NIFD, frontotemporal lobar degeneration neuroimaging initiative; PING, Pediatric Imaging, Neurocognition and Genetics; PPMI, Parkinson's Progression Markers Initiative; RSNA, Radiological Society of North America 2019 Brain CT Hemorrhage data set; UCD-ADC Brain Bank, University of California Davis Alzheimer's Disease Center Brain Bank

**Modalities:** CT, computed tomography; EEG, electroencephalography; fMRI, functional magnetic resonance imaging; sMRI, structural magnetic resonance imaging; SPECT, single-photon emission computed tomography; T1w, T1-weighted [magnetic resonance imaging]

**Tasks:** AD, Alzheimer's disease; CN, cognitively normal; FTD, frontotemporal dementia; MCI, mild cognitive impairment

**Interpretability methods:** FM, feature maps; Grad-CAM, gradient-weighted class activation mapping; LIME, local interpretable model-agnostic explanations; LRP, layer-wise relevance; SHAP, SHapley Additive exPlanations

This section ends with the presentation of benchmarks conducted in the literature to compare different interpretability methods in the context of brain disorders.

## 3.1 Weight Visualization Applied to Neuroimaging

As the focus of this chapter is on non-transparent models, such as deep learning ones, weight visualization was only rarely found. However, this was the method chosen by Cecotti and Gräser [26], who developed a CNN architecture adapted to weight visualization to detect P300 signals in electroencephalograms (EEG). The input of this network is a matrix with rows corresponding to the 64 electrodes and columns to 78 time points. The two first layers of the networks are convolutions with rectangular filters: the first filters (size 1×64) combine the electrodes, whereas the second ones (13×1) find time patterns. Then, it is possible to retrieve a coefficient per electrode by summing the weights associated with this electrode across the different filters and to visualize the results in the electroencephalogram space as shown in Fig. 15.
<figure><img src="image_249.png" /><figcaption>**Fig. 15** Relative importance of the electrodes for signal detection in EEG using two different architectures (CNN-1 and CNN-3) and two subjects (A and B) using CNN weight visualization. Dark values correspond to weights with a high absolute value while white values correspond to weights close to 0. ©2011 IEEE. (Reprinted, with permission, from [26])</figcaption></figure>

**3.2 Feature Map**
**Visualization Applied**
**to Neuroimaging**

Contrary to the limited application of weight visualization, there is an extensive literature about leveraging individual feature maps and latent spaces to better understand how models work. This goes from the visualization of these maps or their projections [27–29], to the analysis of neuron behavior [30, 31], through sampling in latent spaces [29].

Oh et al. [27] displayed the feature maps associated with the convolutional layers of CNNs trained for various Alzheimer's disease status classification tasks (Fig. 16). In the first two layers, the extracted features were similar to white matter, cerebrospinal fluid, and skull segmentations, while the last layer showcased sparse, global, and nearly binary patterns. They used this example to emphasize the advantage of using CNNs to extract very abstract and complex features rather than using custom algorithms for feature extraction [27].

Another way to visualize a feature map is to project it in a two- or three-dimensional space to understand how it is positioned with respect to other feature maps. Abrol et al. [28] projected the features obtained after the first dense layer of a ResNet architecture onto a two-dimensional space using the classical t-distributed stochastic neighbor embedding (t-SNE) dimensionality reduction technique. For the classification task of Alzheimer's disease
<figure><img src="image_250.png" /><figcaption>**Fig. 16** Representation of a selection of feature maps (outputs of 4 filters on 10 for each layer) obtained for a single individual. (Adapted from [27] (CC BY 4.0))</figcaption></figure>

statuses, they observed that the projections were correctly ordered according to the disease severity, supporting the correctness of the model [28]. They partitioned these projections into three groups: Far-AD (more extreme Alzheimer's Disease patients), Far-CN (more extreme Cognitively Normal participants), and Fused (a set of images at the intersection of AD and CN groups). Using a t-test, they were able to detect and highlight voxels presenting significant differences between groups (Fig. 17).

Biffi et al. [29] not only used feature map visualization but also sampled the feature space. Indeed, they trained a ladder variational autoencoder framework to learn hierarchical latent representations of 3D hippocampal segmentations of control subjects and Alzheimer's disease patients. A multilayer perceptron was jointly trained on top of the highest two-dimensional latent space to classify anatomical shapes. While lower spaces needed a dimensionality reduction technique (i.e., t-SNE), the highest latent space could directly be visualized, as well as the anatomical variability it captured in the initial input space, by leveraging the generative process of the model. This sampling enabled an easy visualization and quantification of the anatomical differences between each class.
<figure><img src="image_251.png" /><figcaption>**Fig. 17** Difference in neuroimaging space between groups defined thanks to t-SNE projection. Voxels showing significant differences post false discovery rate (FDR) correction (*p* < 0.05) are highlighted. (Reprinted from Journal of Neuroscience Methods, 339, [28], 2020, with permission from Elsevier)</figcaption></figure>

Finally, it may be very informative to better understand the behavior of neurons and what they are encoding. After training deep convolutional autoencoders to reconstruct MR images, segmented gray matter maps, and white matter maps, Martinez-Murcia et al. [30] computed correlations between each individual hidden neuron value and clinical information (e.g., age, mental state examination) which allowed them to determine to which extent this information was encoded in the latent space. This way they determined which clinical data was the most strongly associated. Using a collection of nine different MRI data sets, Leming et al. [31] trained CNNs for various classification tasks (autism vs typically developing, male vs female, and task vs rest). They computed a diversity coefficient for each filter of the second layer based on its output feature map. They counted how many different data sets maximally activated each value of this feature map: if they were mainly activated by one source of data, the coefficient would be close to 0, whereas if they were activated by all data sets, it would be close to 1. This allows assessing the layer stratification, i.e., to understand if a given filter was mostly maximally activated by one phenotype or by a diverse population. They found out that a few filters were only maximally activated by images from a single MRI data set and that the diversity coefficient was not normally distributed across filters, having generally two peaks at the
beginning and at the end of the spectrum, respectively, exhibiting
the stratification and strongly diverse distribution of the filters.

3.3 Back- Propagation Methods Applied to Neuroimaging

Back-propagation methods are the most popular methods to interpret models, and a wide range of these algorithms have been used to study brain disorders: standard and guided back-propagation [27, 34, 37, 41, 48], gradient input [36–38], Grad-CAM [35, 36], guided Grad-CAM [49], LRP [34, 36–38], DeconvNet [36], and deep Taylor Decomposition [36].

3.3.1 Single Interpretation

Some studies implemented a single back-propagation method and exploited it to find which brain regions are exploited by their algorithm [27, 31, 41], to validate interpretability methods [38], or to provide attribution maps to physicians to improve clinical guidance [35].

Oh et al. [27] used the standard back-propagation method to interpret CNNs for classification of Alzheimer's disease statuses. They showed that the attribution maps associated with the prediction of the conversion of prodromal patients to dementia included more complex representations, less focused on the hippocampi, than the ones associated with classification between demented patients from cognitively normal participants (see Fig. 18). In the context of autism, Leming et al. [31] used the Grad-CAM

<figure><img src="image_252.png" /><figcaption>**Fig. 18** Distribution of discriminant regions obtained with gradient back-propagation in the classification of demented patients and cognitively normal participants (top part, AD vs CN) and the classification of stable and progressive mild cognitive impairment (bottom part, sMCI vs pMCI). (Adapted from [27] (CC BY 4.0))</figcaption></figure>
algorithm to determine the most important brain connections from
functional connectivity matrices. However, the authors pointed out
that without further work, this visualization method did not allow
understanding the underlying reason of the attribution of a given
feature: for instance, one cannot know if a set of edges is important
because it is under-connected or over-connected. Finally, Hu et al.
[41] used attribution maps produced by guided back-propagation
to quantify the difference in the regions used by their network to
characterize Alzheimer's disease or frontotemporal dementia.

The goal of Eitel et al. [38] was different. Instead of identifying
brain regions related to the classification task, they exhibited with
LRP that transfer learning between networks trained on different
diseases (Alzheimer's disease to multiple sclerosis) and different
MRI sequences enabled obtaining attribution maps focused on a
smaller number of lesion areas. However, the authors pointed out
that it would be necessary to confirm their results on larger
data sets.

Finally, Burduja et al. [35] trained a CNN-LSTM model to
detect various hemorrhages from brain computed tomography
(CT) scans. For each positive slice coming from controversial or
difficult scans, they generated Grad-CAM-based attribution maps
and asked a group of radiologists to classify them as correct, par-
tially correct, or incorrect. This classification allowed them to
determine patterns for each class of maps and better understand
which characteristics radiologists expected from these maps to be
considered as correct and thus useful in practice. In particular,
radiologists described maps including any type of hemorrhage as
incorrect as soon as some of the hemorrhages were not highlighted,
while the model only needed to detect one hemorrhage to correctly
classify the slice as pathological.

Papers described in this section used several interpretability meth-
ods and compared them in their particular context. However, as the
benchmark of interpretability methods is the focus of Subheading
4.3, which also include other types of interpretability than back-
propagation, we will only focus here on what conclusions were
drawn from the attribution maps.

Dyrba et al. [36] compared DeconvNet, guided back-
propagation, deep Taylor decomposition, gradient input, LRP
(with various rules), and Grad-CAM methods for classification of
Alzheimer's disease, mild cognitive impairment, and normal cogni-
tion statuses. In accordance with the literature, they obtained a
highest attention given to the hippocampus for both prodromal
and demented patients.

Böhle et al. [34] compared two methods, LRP with β-rule and guided back-propagation for Alzheimer's disease status classification. They found that LRP attribution maps highlight the individual differences between patients and then that they could be used as a tool for clinical guidance.

### *Comparison of Several Interpretability Methods*
## 3.4 Perturbation Methods Applied to Neuroimaging

The standard perturbation method has been widely used in the study of Alzheimer's disease [32, 37, 45, 48] and related symptoms (amyloid-β pathology) [49]. However, most of the time, authors do not train their model with perturbed images. Hence, to generate explanation maps, the perturbation method uses images outside the distribution of the training set, which may call into question the relevance of the predictions and thus the reliability of attention maps.

### 3.4.1 Variants of the Perturbation Method Tailored to Neuroimaging

Several variations of the perturbation method have been developed to adapt to neuroimaging data. The most common variation in brain imaging is the brain area perturbation method, which consists in perturbing entire brain regions according to a given brain atlas, as done in [27, 28, 48]. In their study of Alzheimer's disease, Abrol et al. [28] obtained high values in their attribution maps for the usually discriminant brain regions, such as the hippocampus, the amygdala, the inferior and superior temporal gyruses, and the fusiform gyrus. Rieke et al. [48] also obtained results in accordance with the medical literature and noted that the brain area perturbation method led to a less scattered attribution map than the standard method (Fig. 19). Oh et al. [27] used the method to compare the attribution maps of two different tasks: (1) demented patients vs cognitively normal participants and (2) stable vs progressive mild cognitively impaired patients and noted that the regions targeted

<figure><img src="image_253.png" /><figcaption>**Fig. 19** Mean attribution maps obtained on demented patients. The first row corresponds to the standard and the second one to the brain area perturbation method. (Reprinted by permission from Springer Nature Customer Service Centre GmbH: Springer Nature, MLCN 2018, DLF 2018, IMIMIC 2018: Understanding and Interpreting Machine Learning in Medical Image Computing Applications, [48], 2018)</figcaption></figure>
for the first task were shared with the second one (medial temporal lobe) but that some regions were specific to the second task (parts of the parietal lobe).

Gutiérrez-Becker and Wachinger [40] adapted the standard perturbation method to a network that classified clouds of points extracted from neuroanatomical shapes of brain regions (e.g., left hippocampus) between different states of Alzheimer's disease. For the perturbation step, the authors set to 0 the coordinates of a given point x and the ones of its neighbors to then assess the relevance of the point x. This method allows easily generating and visualizing a 3D attribution map of the shapes under study.

More advanced perturbation-based methods have also been used in the literature. Nigri et al. [45] compared a classical perturbation method to a swap test. The swap test replaces the classical perturbation step by a swapping step where patches are exchanged between the input brain image and a reference image chosen according to the model prediction. This exchange is possible as brain images were registered and thus brain regions are positioned in roughly the same location in each image.

Finally, Thibeau-Sutre et al. [51] used the optimized version of the perturbation method to assess the robustness of CNNs in identifying regions of interest for Alzheimer's disease detection. They applied optimized perturbations on gray matter maps extracted from T1w MR images, and the perturbation method consisted in increasing the value of the voxels to transform patients into controls. This process aimed at stimulating gray matter reconstruction to identify the most important regions that needed to be “de-atrophied” to be considered again as normal. However, they unveiled a lack of robustness of the CNN: different retrainings led to different attribution maps (shown in Fig. 20) even though the performance did not change.

Distillation methods are less commonly used, but some very interesting use cases can be found in the literature on brain disorders, with methods such as LIME [44] or SHAP [33].

<figure><img src="image_254.png" /><figcaption>**Fig. 20** Coronal view of the mean attribution masks on demented patients obtained for five reruns of the same network with the optimized perturbation method. (Adapted with permission from Medical Imaging 2020: Image Processing, [51].)</figcaption></figure>
<figure><img src="image_255.png" /><figcaption>**Fig. 21** Mean absolute feature importance (SHAP values) averaged across all subjects for XGBoost on regional thicknesses (red) and areas (green). (Adapted from [33] (CC BY 4.0))</figcaption></figure>

Magesh et al. [44] used LIME to interpret a CNN for Parkinson's disease detection from single-photon emission computed tomography (SPECT) scans. Most of the time, the most relevant regions are the putamen and the caudate (which is clinically relevant), and some patients also showed an anomalous increase in dopamine activity in nearby areas, which is a characteristic feature of late-stage Parkinson's disease. The authors did not specify how they extracted the “super-pixels” necessary to the application of the method, though it could have been interesting to consider neuro-anatomical regions instead of creating the voxel groups with an agnostic method.

Ball et al. [33] used SHAP to obtain explanations at the individual level from three different models trained to predict participants' age from regional cortical thicknesses and areas: regularized linear model, Gaussian process regression, and XGBoost (Fig. 21). The authors exhibited a set of regions driving predictions for all models and showed that regional attention was highly correlated on average with weights of the regularized linear model. However, they showed that while being consistent across models and training folds, explanations of SHAP at the individual level were generally not correlated with feature importance obtained from the weight analysis of the regularized linear model. The authors also exemplified that the global contribution of a region to the final prediction error (“brain age delta”), even with a high SHAP value, was in general small, which indicated that this error was best explained by changes spread across several regions [33].

## 3.6 Intrinsic Methods Applied to Neuroimaging

### 3.6.1 Attention Modules

Attention modules have been increasingly used in the past couple of years, as they often allow a boost in performance while being rather easy to implement and interpret. To diagnose various brain diseases from brain CT images, Fu et al. [39] built a model integrating a “two-step attention” mechanism that selects both the most
<figure><img src="image_256.png" /><figcaption>**Fig. 22** Attribution maps (left, in-house database; right, ADNI database) generated by an attention mechanism module, indicating the discriminant power of various brain regions for Alzheimer's disease diagnosis. (Adapted from [42] (CC BY 4.0))</figcaption></figure>

important slices and the most important pixels in each slice. The authors then leveraged these attention modules to retrieve the five most suspicious slices and highlight the areas with the more significant attention.

In their study of Alzheimer's disease, Jin et al. [42] used a 3D attention module to capture the most discriminant brain regions used for Alzheimer's disease diagnosis. As shown in Fig. 22, they obtained significant correlations between attention patterns for two independent databases. They also obtained significant correlations between regional attention scores of two different databases, which indicated a strong reproducibility of the results.

**3.6.2 Modular Transparency**

Modular transparency has often been used in brain imaging analysis. A possible practice consists in first generating a target probability map of a black box model, before feeding this map to a classifier to generate a final prediction, as done in [43, 46].

Qiu et al. [46] used a convolutional network to generate an attribution map from patches of the brain, highlighting brain regions associated with Alzheimer's disease diagnosis (see Fig. 23). Lee et al. [43] first parcellated gray matter density maps into 93 regions. For each of these regions, several deep neural networks were trained on randomly selected voxels, and their outputs were averaged to obtain a mean regional disease probability. Then, by concatenating these regional probabilities, they generated a region-wise disease probability map of the brain, which was further used to perform Alzheimer's disease detection.

The approach of Ba et al. [24] was also applied to Alzheimer's disease detection [50] (preprint). Though that work is still a preprint, the idea is interesting as it aims at reproducing the way a radiologist looks at an MR image. The main difference with [24] is
<figure><img src="image_257.png" /><figcaption>**Step 1:** Random sampling of patches for fully convolutional network training</figcaption></figure>

<figure><img src="image_258.png" /><figcaption>**Step 2:** Generate probability maps after fully convolutional network training</figcaption></figure>

**Fig. 23** Randomly selected samples of T1-weighted full MRI volumes are used as input to learn the Alzheimer's disease status at the individual level (Step 1). The application of the model to whole images leads to the generation of participant-specific disease probability maps of the brain (Step 2). (Adapted from Brain: A Journal of Neurology, 143, [46], 2020, with permission of Oxford University Press)

the initialization, as the context network does not take as input the whole image but clinical data of the participant. Then the framework browses the image in the same way as in the original paper: a patch is processed by a recurrent neural network and from its internal state the glimpse network learns which patch should be looked at next. After a fixed number of iterations, the internal state of the recurrent neural network is processed by a classification network that gives the final outcome. The whole system is interpretable as the trajectory of the locations (illustrated in Fig. 24) processed by the framework allows understanding which regions are more important for the diagnosis. However, this framework may have a high dependency to clinical data: as the initialization depends on scores used to diagnose Alzheimer's disease, the classification network may learn to classify based on the initialization only, and most of the trajectory may be negligible to assess the correct label.

Another framework, the DaniNet, proposed by Ravi et al. [47], is composed of multiple networks, each with a defined function, as illustrated in Fig. 25.

* The conditional deep autoencoder (in orange) learns to reduce the size of the slice *x* to a latent variable *Z* (encoder part) and
<figure><img src="image_259.png" /><figcaption>**Fig. 24** Trajectory taken by the framework for a participant from the ADNI test set. A bounding box around the first location attended to is included to indicate the approximate size of the glimpse that the recurrent neural network receives; this is the same for all subsequent locations. (Adapted from [50]. Permission to reuse was kindly granted by the authors)</figcaption></figure>

<figure><img src="image_260.png" /><figcaption>**Fig. 25** Pipeline used for training the proposed DaniNet framework that aims to learn a longitudinal model of the progression of Alzheimer's disease. (Adapted from [47] (CC BY 4.0))</figcaption></figure>

then to reconstruct the original image based on *Z* and two additional variables: the diagnosis and age (generator part). Its performance is evaluated thanks to the reconstruction loss $L^{rec}$.