# Applications of interpretable deep learning in neuroimaging: a comprehensive
review

Lindsay Munroe  Mariana da Silva  Faezeh Heidari  Irina Grigorescu  Simon
Dahan  Emma C. Robinson  Maria Deprez  Po-Wah So

###### Abstract

Clinical adoption of deep learning models has been hindered, in part, because
the “black-box” nature of neural networks leads to concerns regarding their
trustworthiness and reliability. These concerns are particularly relevant in
the field of neuroimaging due to the complex brain phenotypes and inter-
subject heterogeneity often encountered. The challenge can be addressed by
interpretable deep learning (iDL) methods that enable the visualisation and
interpretation of the inner workings of deep learning models. This study
systematically reviewed the literature on neuroimaging applications of iDL
methods and critically analysed how iDL explanation properties were evaluated.
Seventy-five studies were included, and ten categories of iDL methods were
identified. We also reviewed five properties of iDL explanations that were
analysed in the included studies: biological validity, robustness, continuity,
selectivity, and downstream task performance. We found that the most popular
iDL approaches used in the literature may be sub-optimal for neuroimaging
data, and we discussed possible future directions for the field.

###### keywords:

Interpretable Deep Learning, Explainable AI, Neuroimaging, Intrinsic
Interpretability

lrp short = LRP , long = Layer-wise Relevance Propagation , tag = abbrev dl
short = DL , long = deep learning , tag = abbrev ml short = ML , long =
machine learning , tag = abbrev gan short = GAN , long = generative
adversarial network , tag = abbrev idl short = iDL , long = interpretable deep
learning , tag = abbrev xai short = XAI , long = explainable artificial
intelligence , tag = abbrev gdpr short = GDPR , long = General Data Protection
Regulation , tag = abbrev hiv short = HIV , long = Human Immunodeficiency
Virus , tag = abbrev cnn short = CNN , long = convolutional neural network ,
tag = abbrev ae short = AE , long = autoencoder , tag = abbrev vae short = VAE
, long = variational autoencoder , tag = abbrev nn short = NN , long = neural
network , tag = abbrev rnn short = RNN , long = recurrent neural network , tag
= abbrev crnn short = C-RNN , long = convolutional recurrent neural network ,
tag = abbrev mlp short = MLP , long = multi-layer perceptron , tag = abbrev
dnn short = DNN , long = deep neural network , tag = abbrev gap short = GAP ,
long = global average pooling , tag = abbrev svm short = SVM , long = support
vector machine , tag = abbrev smri short = sMRI , long = structural magnetic
resonance imaging , tag = abbrev ic short = IC , long = independent component
, tag = abbrev fbn short = FBN , long = functional brain networks , tag =
abbrev fmri short = fMRI , long = functional magnetic resonance imaging , tag
= abbrev rsfmri short = rs-fMRI , long = resting-state functional magnetic
resonance imaging , tag = abbrev mri short = MRI , long = magnetic resonance
imaging , tag = abbrev dti short = DTI , long = diffusion tensor imaging , tag
= abbrev srddl short = sr-DDL , long = structurally-regularized Dynamic
Dictionary Learning , tag = abbrev mmse short = MMSE , long = mini-mental
state examination , tag = abbrev scm short = SCM , long = Structural Causal
Model , tag = abbrev dscm short = DSCM , long = Deep Structural Causal Model ,
tag = abbrev fnc short = FNC , long = functional network connectivity , tag =
abbrev tc short = TCs , long = time courses , tag = abbrev cn short = CN ,
long = cognitively normal , tag = abbrev ad short = AD , long = Alzheimer’s
Disease , tag = abbrev adhd short = ADHD , long = attention deficit
hyperactivity disorder , tag = abbrev mci short = MCI , long = mild cognitive
impairment , tag = abbrev pmci short = pMCI , long = progressive mild
cognitive impairment , tag = abbrev smci short = sMCI , long = stable mild
cognitive impairment , tag = abbrev adni short = ADNI , long = Alzheimer’s
Disease Neuroimaging Initiative , tag = abbrev aibl short = AIBL, long =
Australian Imaging Biomarkers and Lifestyle Study of Ageing, tag = abbrev gnn
short=GNN, long=graph neural network, tag=abbrev gcn short=GCN, long=graph
convolutional network, tag=abbrev asd short=ASD, long=Autism Spectrum
Disorders, tag=abbrev pd short=PD, long=Parkinson’s Disease, tag=abbrev scz
short=SCZ, long=Schizophrenia, tag=abbrev hc short=HC, long= healthy control ,
tag=abbrev ms short=MS, long=Multiple Sclerosis, tag=abbrev lstm short=LSTM,
long=long short-term memory, tag=abbrev relu short=ReLU, long=Rectified Linear
Unit, tag=abbrev nlp short=NLP, long=natural language processing, tag=abbrev
ct short=CT, long=Computed Tomography, tag=abbrev brats short=BraTS,
long=Brain Tumor Segmentation, tag=abbrev ich short=ICH, long=intracerebral
hemorrhage, tag=abbrev spect short=SPECT, long=Single Photon Emission Computed
Tomography, tag=abbrev abide short=ABIDE, long=Autism Brain Imaging Data
Exchange, tag=abbrev ncc short=NCC, long=normalised cross correlation,
tag=abbrev hcp short=HCP, long=Human Connectome Project, tag=abbrev dhcp
short=dHCP, long=developing Human Connectome Project, tag=abbrev kki
short=KKI, long=Kennedy Krieger Institute, tag=abbrev wrat short=WRAT,
long=Wide Range Achievement Test, tag=abbrev pnc short=PNC, long=Philadelphia
Neurodevelopmental Cohort, tag=abbrev hande short=H&E, long=hematoxylin &
eosin, tag=abbrev oasis short=OASIS, long=Open Access Series of Imaging
Studies, tag=abbrev fdgpet short=FDG-PET, long=fluorodeoxyglucose positron
emission tomography, tag=abbrev fc short=FC, long=functional connectivity,
tag=abbrev us short=US, long=ultrasound, tag=abbrev ftd short=FTD,
long=Frontotemporal Dementia, tag=abbrev ddpm short=DDPM, long=Denoising
Diffusion Probabilistic Model, tag=abbrev ppmi short=PPMI, long=Parkinson’s
Progression Markers Initiative, tag=abbrev

[inst1]organization=Department of Neuroimaging, addressline=King’s College
London, city=London, country=UK

[inst2]organization=School of Biomedical Engineering and Imaging
Sciences,addressline=King’s College London, city=London, country=UK

[inst3]organization=Institute of Clinical Medicine, addressline=University of
Eastern Finland, city=Kuopio, country=Finland

##  1 Introduction

Traditionally, analysis and interpretation of neuroimaging data requires
specialised expertise, is often laborious, and is subject to inter-observer
variability. Therefore, dl has become a popular tool in neuroimaging in recent
years, driven by the rise in computer processing power as well as increased
access to large medical imaging datasets and the success of novel model
architectures. In neuroimaging, dl has been applied to segmentation [19, 24,
87], super-resolution [56, 148], image synthesis [36, 65, 122] and
classification [78, 17], among other applications. Despite the success of dl
for analysing and interpreting neuroimaging data, adoption remains limited
partly because dl models are often opaque and considered to be “black boxes”.
In other words, the internal workings of dl models are not comprehensible to
humans, which leads to concerns regarding their reliability and
trustworthiness. Indeed, such “black box” models do not satisfy European gdpr
legal requirements to provide “information about the logic involved” [40].

###  1.1 Advantages of interpretable deep learning

idl has been proposed to address the opacity problem of dl models, for
example, by producing explanations that highlight brain regions that are most
relevant for the model predictions. idl methods can support the translation of
dl to the clinic by providing healthcare practitioners with explanations to
verify predictions and communicate with patients. Additionally, deep learning
practitioners can leverage idl to debug their models and identify cases where
a model makes the right decision for the wrong reason [68]. idl methods can
also be employed to test scientific hypotheses, such as identifying brain
regions involved in disease pathogenesis.

###  1.2 Evaluation of iDL explanations

A challenging aspect of idl is assessing the quality of explanations because
such ground truths are typically unavailable. While experts such as
clinicians, pathologists, or imaging scientists can qualitatively evaluate
explanations, quantitative and automated metrics are often preferred,
particularly when access to medical professionals is limited. Researchers have
proposed various quantitative methods to evaluate desirable properties of idl
explanations, with a particular focus on assessing fidelity and robustness
(e.g., [2], [90], [59], [46], [114], [151], [67]).

Fidelity refers to the extent to which explanations reflect the inner workings
of the associated deep learning model. Fidelity is usually evaluated by
removing features or comparing the explanations to ground truth, if available.
In computer vision, feature-removal approaches generally involve masking image
regions with the highest relevance in the associated explanation, obtaining
predictions for the modified images, and then measuring the change in model
output or accuracy. A substantial drop in accuracy indicates that the
explanations faithfully highlight image features attended to by the model. For
example, Montavon et al. developed a procedure to assess fidelity in which
they iteratively removed $4\times 4$ patches from images with the highest
relevance and plotted the number of patches removed against model output score
[90]. In another example of fidelity evaluation, Adebayo et al. randomised
model parameters and data labels as two sanity checks to assess whether idl
explanations truly reflect either the model mechanisms or the relationship
between image features and the label [2]. Alternatively, explanations can be
compared to ground truth maps of image features the model is expected to
attend to when making predictions. For instance, bounding box annotations for
objects in natural images have been used as ground truth and the ratio of mean
relevance outside vs. inside the bounding box has been calculated to assess
the fidelity of explanations [67].

Robustness can be described as the stability of model explanations under
varying modelling conditions. For example, Montavon et al. introduced the
concept of continuity, which means that an idl method should produce similar
explanations for similar input images [90]. The evaluation of idl methods is
an active research field, and for a comprehensive review of the topic, we
recommend readers to refer to Alangari et al.’s work [4].

###  1.3 Classification of interpretable deep learning methods

Two main categories of idl methods exist: post-hoc and intrinsic. Post-hoc
methods use reverse engineering to generate an explanation from a “black-box”
model after training. In contrast, intrinsic methods incorporate interpretable
components into the model architecture during the design phase. Another way to
classify interpretable methods is by local vs. global explanations. Local
explanations focus on individual samples and thereby increase trust in the
model outcomes, whereas global explanations seek to provide a deeper
understanding of the mechanism by which the model works.

###  1.4 Study objectives

The objectives of this review are:

  1. 1.

To systematically review idl methods applied to neuroimaging studies.

  2. 2.

To review the evaluation of idl explanations in the studies included in this
review, explicitly identifying the properties evaluated and associated
quantitative metrics proposed.

To the best of our knowledge, this is the first study to systematically review
both post-hoc and intrinsic idl methods in the field of neuroimaging.

We have further sub-classified idl methods of the two categories (Table 2).
Initially, we introduce five post-hoc methods (Section 4.1) and five intrinsic
methods 4.2 before reviewing applications to neuroimaging for each method
(Section 5). Finally, we consider how idl explanations were evaluated across
the included studies (Section 6).

Figure 1:  Comparison of post-hoc interpretability maps and generative
interpretability methods applied to the classification of Alzheimer’s disease
(AD) vs Mild cognitive impairment (MCI) in brain MRI volumes. The real disease
map is the “ground-truth” shown for comparison. Figure adapted from Bass et
al. [9]

##  2 Systematic search methodology

We identified relevant articles for this review by querying PubMed, Web of
Science, Google Scholar and arXiv using the following search terms: 1.
explainable, 2. XAI, 3. interpretable, 4. explainability, 5. interpretability,
6. causal reasoning, 7. counterfactuals, 8. deep learning, 9. AI, 10. neural
network, 11. machine learning, 12. brain imaging, 13. neuroimaging, 14.
neuroradiology. The search terms were combined in the logical statement (1 OR
2 OR 3 OR 4 OR 5 OR 6 OR 7) AND (8 OR 9 OR 10 OR 11) AND (12 OR 13 OR 14).
Articles from 2015 were included for PubMed and Google Scholar, whereas all
years were included for Web of Science and arXiV due to the small number of
articles returned.

Articles were initially screened based on the article title and abstract and
accepted or rejected from a full-text review based on inclusion and exclusion
criteria (Table 1). Only the first 500 results from Google Scholar were
screened because later results were irrelevant. Finally, we extracted the
pertinent information from all accepted articles into a spreadsheet for
further analysis.

Inclusion/ exclusion criteria for article screening
---
Include…both _in-vivo_ and _ex-vivo_ imaging.
Exclude…non-human subjects.
Include…the following imaging modalities: structural and functional magnetic
resonance imaging, computed tomography, and positron emission tomography.
Exclude…electroencephalogram and magnetoencephalography data.
Exclude…non-peer reviewed articles.
Exclude…non-English language articles.
Exclude…PhD and Masters theses.
Exclude…reviews, surveys, opinion articles and books. Articles must implement
at least one interpretable deep learning method.
Exclude…interpretable methods applied to machine learning models other than
neural networks.
Exclude…for quality control. For example, if the explanations could not be
reasonably interpreted.
Table 1: Inclusion and exclusion criteria

##  3 Overview

Figure 2: The Preferred Reporting Items for Systematic Reviews and Meta-
Analyses (PRISMA) flowchart

The number of articles returned was 712 for PubMed, 88 for Web of Science,
1000 for Google Scholar (the upper limit), and 189 for arXiV. After title and
abstract screening, the number of accepted articles was 30 for PubMed, 26 for
Web of Science, 127 for Google Scholar and 58 for arXiV. After full-text
review and removal of duplicates, and added articles after a refresh, the
number of accepted articles was 75 (Fig. 2). Table 2 summarises the methods
and papers introduced in this review.

Methods | Papers
---|---
Post-hoc methods |  Tab. 3 Perturbation-based methods
Tab..10$\blacktriangleright$ Disease Classification |  [31], [79], [150], [95], [134], [120], [132]
|  [81], [73], [26], [120], [85], [149]
Tab..10$\blacktriangleright$ Sex Classification | [54]
Tab..10$\blacktriangleright$ Brain Age Regression | [16]
Tab. 4 Gradient-based methods
Tab..10$\blacktriangleright$ Disease Classification |  [97], [31], [32], [72], [157]
Tab..10$\blacktriangleright$ Brain Age Regression | [71]
Tab..10$\blacktriangleright$ Cognitive Task Decoding |  [84], [51]
Tab. 5 Backpropagation-based methods
Tab..10$\blacktriangleright$ Disease Classification |  [17], [31]
Tab..10$\blacktriangleright$ Sex Classification | [54]
Tab..10$\blacktriangleright$ Cognitive Task Decoding | [135]
Tab. 6 Class activation maps
Tab..10$\blacktriangleright$ Disease Classification |  [153], [150], [5], [57], [132], [142], [75], [143], [70]
Tab..10$\blacktriangleright$ Sex Classification |  [35], [58], [54]
Tab..10$\blacktriangleright$ Tissue Segmentation | [93]
Tab..10$\blacktriangleright$ Cognitive Score Prediction |  [49], [106]
Tab. 7 Weight Analysis
Tab..10$\blacktriangleright$ Disease Classification |  [30], [74]
Tab..10$\blacktriangleright$ Tissue Segmentation |  [64], [93]
Tab..10$\blacktriangleright$ Cognitive Task Decoding | [74]
Intrinsic methods |  Tab. 8 Disentangled latent spaces
Tab…10$\blacktriangleright$ Image Generation |  [156], [91], [160], [154], [99]
Tab…10$\blacktriangleright$ Disease Classification |  [3], [139]
Tab…10$\blacktriangleright$ Brain Age Regression | [47]
Tab. 9 Hybrid models
Tab…10$\blacktriangleright$ Disease Classification |  [1], [69], [77], [94], [105], [88], [92], [145], [55], [104]
Tab…10$\blacktriangleright$ Brain Age Regression | [43]
Tab…10$\blacktriangleright$ Clinical Score Regression | [121]
Tab. 10 Generative models
Tab…10$\blacktriangleright$ Disease Classification |  [12], [10], [9], [66], [80]
Tab…10$\blacktriangleright$ Brain Age Regression | [9]
Tab…10$\blacktriangleright$ Tissue Segmentation |  [14], [116], [146]
Tab. 11 Deep structural causal models
Tab…10$\blacktriangleright$ Image Generation |  [100], [108], [107]
Tab. 12 Attention-based models
|  Tab…10$\blacktriangleright$ Disease Classification |  [53], [155], [117]
|  Tab…10$\blacktriangleright$ Tissue Segmentation | [38]
|  Tab…10$\blacktriangleright$ Brain Age Regression | [25]
Table 2: Summary of both post-hoc and intrinsic methods

##  4 Methods

###  4.1 Post-hoc methods

Post-hoc interpretability methods, as the name suggests, analyse model
decisions after a network has been trained. While some post-hoc methods are
model agnostic, i.e., they can be applied to any ml model, in some cases, they
are only applicable to a specific family of models, such as cnn. Agnostic
post-hoc methods can be applied to “black-box” models without requiring
knowledge of the model parameters, as they generally analyse feature input and
output pairs. Alternatively, post-hoc methods may require access to pre-
trained model information (e.g., model weights) as for gradient-based and
weight-analysis methods.

####  4.1.1 Perturbation-based methods

Perturbation-based methods explicitly alter the input features and measure the
change in the model prediction between the original and perturbed data to
discover relevant features. The most salient features for a model decision are
those that produce the greatest change in the model prediction when perturbed.
Perturbation-based methods mainly differ according to how they alter the input
features.

Several perturbation-based methods occlude input features. For example,
Occlusion obstructs regions of an input image in a patch-wise fashion [152].
For every patch location, the change in the model output between the original
and occluded image is calculated to form a sensitivity map. For classification
tasks, sensitivity is the change in predicted probability $\mathbb{P}(c)$ of
the image belonging to a class-of-interest $c$, as shown in Fig. 3. For
regression tasks, the residual difference of the model prediction is assessed.

Figure 3: Example of Occlusion applied to an MRI image. In a patch-wise
manner, a tile of the image is occluded, and the occluded image is fed to a
neural network (NN) for prediction. The difference in predicted probability
between the original and occluded image is assigned to the patch location in
the occlusion map. Patches that result in the greatest change in prediction
when occluded are interpreted as the most important for the model task [152].

Meaningful Perturbations follows a similar approach of occluding image regions
but uses gradient descent to learn the occlusion mask that obfuscates the
smallest region of the image that renders the model unable to correctly
classify the masked image [34]. The masking process may replace pixel values
with a constant value, Gaussian noise or by blurring.

Also incorporating occlusion, Local Interpretable Model-Agnostic Explanations
(LIME) approximates a “black-box model” locally to an input $x$; then an
interpretable ml model, such as a linear model, is trained to mimic the
“black-box” model predictions for occluded samples of $x$ [109]. First,
several perturbed images are generated from a given image $I_{0}$; a single
perturbed image $I$ is generated by switching off a random subset of
superpixels of $I_{0}$, where a superpixel is a set of neighbouring pixels
with similar intensity. A sparse linear model is trained on the corresponding
binary features $I^{\prime}=(b_{1},...,b_{n})$ where $b_{i}=0$ if superpixel
is switched off to generate image $I$ and $b_{i}=1$ otherwise. Training labels
for the linear model are the “black-box” model predictions for perturbed
images $I$. The feature importance of the $i^{th}$ superpixel in $I_{0}$ is
given by the associated linear model coefficient of $b_{i}$.

In contrast to occluding image regions, several perturbation methods swap
image regions or input features with those of another subject so that the
altered image still appears realistic. Such an approach was proposed in the
Swap Test, where a reference image is selected that is from a different class
to the image-of-interest [95]. For example, for an image classified as ad, the
reference image is randomly selected from healthy control images. In a patch-
wise manner, a patch in the reference image is replaced with the corresponding
patch in the image-of-interest and the change in model output between the
reference and altered reference image is computed. The process is repeated for
several randomly selected reference images and averaged.

Similarly, Permutation Feature Importance [33] randomly permutes values of
each input feature across samples. Let $P_{\text{orig}}$ be the model
performance on the original data and $P_{\text{perm}}$ be the model
performance when feature $j$ has been randomly permuted; then the importance
of feature $j$ is either the ratio $P_{\text{perm}}/P_{\text{orig}}$ or
difference $P_{\text{perm}}-P_{\text{orig}}$. The assumption is if feature $j$
is ignored by the model, then randomly shuffling feature $j$ will not
influence model predictions. In contrast to previously mentioned perturbation
methods, permutation feature importance is a global interpretability method.

Advantages and disadvantages: Perturbation-based methods have the advantage of
being easy to implement and understand; they do not require a specific type of
network nor access to the gradients. These methods may be applied to any
“black-box” model, as they only need access to the input image and output
value. However, these methods are computationally intensive and time-
consuming, as inference is run for each location of the perturbation block.
Another disadvantage is that perturbed images no longer belong to the training
data distribution, so distribution shift may be responsible for any changes in
model output rather than feature relevance [46]. Concerning Occlusion, this
method is also sensitive to the size and the replacement intensity of the
occluded patch [34].

####  4.1.2 Gradient-based methods

Gradient-based methods compute the partial derivative of an output from a
neural network output with respect to each input feature, using the
backpropagation algorithm [112]. The resulting gradient maps visualise how
sensitive a neural network output is to small changes in input feature values,
and they are also referred to as sensitivity maps.

Figure 4: Example of Vanilla Gradients applied to an MRI image. Partial
derivatives for each voxel with respect to the network output score $S_{c}$
for class $c$ are computed. Pixels with the largest gradients are interpreted
to have the greatest influence on the model prediction [124].

Vanilla Gradients was the first gradient-based method used to compute gradient
maps for a cnn trained to classify images [124] (see Fig. 4). Let $I_{0}$ be
an image with $N$ channels; $c$, a class-of-interest; and let $S_{c}(I)$ be
the class score output function of a trained cnn classifier. Then Vanilla
Gradients computes the absolute value of the partial derivative of $S_{c}(I)$
with respect to each voxel in $I_{0}$. Where $N>1$, the maximum value across
channels is returned.

Two main limitations of Vanilla Gradients exist: shattered gradients and the
saturation problem. Firstly, gradient maps are often noisy because of
“shattered gradients”, where similar pixel values have substantially different
partial derivatives of $S_{c}$, thus producing noisy maps [8]. Secondly, there
is the “saturation problem”. The function $S_{c}(I)$ learned by a cnn is non-
linear, therefore the Vanilla Gradient map of $I_{0}$ does not interpret the
behaviour of $S_{c}(I)$ globally, but locally to $I_{0}$. In particular, when
$S_{c}(I)$ is saturated at $I_{0}$, i.e. the gradient is close to zero,
Vanilla Gradients may not reveal image features that cause $S_{c}(I)$ to
substantially change and switch predicted class [123].

Grad $\times$ Input attempts to overcome the shattered gradients limitation
through element-wise multiplication of Vanilla Gradients with $I_{0}$,
producing visually sharper sensitivity maps than Vanilla Gradients [60].

SmoothGrad was also developed to address the shattered gradients limitation of
Vanilla Gradients by adding random noise to the input image to create many
noisy images, then computing the mean of the associated Vanilla Gradients
sensitivity maps [127].

Integrated Gradients addresses the saturation problem of Vanilla Gradients
[130]. Global behaviour is captured by travelling from a baseline image
$I_{b}$ (e.g., an image of all zeros) to the image-of-interest $I_{0}$, and
sampling $m$ images along the path: $I_{b}+\frac{k}{m}(I_{0}-I_{b})$ for all
images $k$ from 1 to $m$. Integrated Gradients then computes the mean Vanilla
Gradients map across the $m$ images (please see Integrated Gradients).
Notably, Integrated Gradients tends to highlight more relevant image features
compared to Vanilla Gradients and SmoothGrad. However, Integrated Gradients
maps may still include noisy gradients from saturated regions of $S_{c}(I)$
[86].

Advantages and disadvantages: Gradient-based methods are fast to run and easy
to understand. However, in addition to the shattered gradients and saturation
problem previously discussed, gradient maps are less able to discriminate
between classes than other interpretable methods.

####  4.1.3 Backpropagation-based methods

Backpropagation-based methods apply rules other than gradients to map the
output score back to the input features to assign feature relevance. The
earliest backpropagation methods for cnn were identical to Vanilla Gradients
aside from their treatment of the relu function.

Specifically, Vanilla Gradients back-propagates through a relu function by
setting a gradient value to zero if the corresponding value in the forward
feature map is negative. In comparison, Guided Backpropagation [128] performs
the same operation and also sets negative gradients to zero. Consequently,
Guided Backpropagation only allows positive gradients, whereas Vanilla
Gradients may produce negative gradients.

Figure 5: Example Layer-wise Relevance Propagation (LRP) applied to an MRI
image. The network output score $S_{c}$ for class $c$ is redistributed
backwards through the network according to the equation shown until the input
image is reached. The pixels with the highest proportion of $S_{c}$ are
interpreted as having the greatest contribution to the model prediction [6].

lrp is another popular backpropagation method [6], as visualised in Fig. 5. In
lrp, the model output score $S_{c}(I_{0})$ is redistributed backwards through
the network, layer by layer, until the input image $I_{0}$ is reached. Each
node (or pixel) is allocated a relevance value, which is the weighted sum of
relevance values of connected nodes in the neighbouring higher layer.
Different lrp rules have a different weighted sum based on the network
parameters, but all follow the relevance conservation principle: relevance
assigned to a node from the neighbouring higher layer is equal to the
relevance passed from that node to the neighbouring earlier layer (for more
information, please see Layer-wise Relevance Propagation (LRP))

Advantages and disadvantages: Analysis carried out by Adebayo et al.
demonstrated Guided Backpropagation maps are independent of higher network
layer parameters and sample labels, which is undesirable for an
interpretability method [2]. Additionally, lrp is sensitive to hyperparameter
selection and may be difficult to tune.

####  4.1.4 Class activation maps

Class Activation Maps (CAM) highlight image regions used by the final layer of
a cnn to classify the input image [159]. To compute CAM visualizations, the
final layer of the network is required to be a gap layer. In a gap-cnn, the
weighted sum of the activation maps in the final layer determines the class
score $S_{c}$ for each class $c$ (Eq. 1):

$S_{c}=\sum_{k}w_{k}^{c}\sum_{x,y}A^{k}_{x,y}=\sum_{x,y}\sum_{k}w_{k}^{c}A^{k}_{x,y}$
(1)

where $A^{k}_{x,y}$ represents the activation of node $k$ in the last
convolutional layer of the network at pixel location $(x,y)$, and $w^{c}_{k}$
represents the importance of node $k$ for the classification of class $c$ (see
Fig. 6).

Figure 6: Example of a Class Activation Map (CAM)-based method where the
activation maps of the final convolutional layer are weighted by the weights
computed through the global average pooling (GAP) layer to produce a class
activation map. Image adapted from Zhou et al. [159].

Then CAM for class $c$ is defined as (Eq. 2):

$CAM^{c}=\sum_{k}w_{k}^{c}A^{k}_{x,y}$ (2)

Hence, the sum of all elements in $CAM^{c}$ is equal to the class score
$S_{c}$.

Gradient-Weighted Class Activation Maps (Grad-CAM) extends CAM to all cnns to
obviate the need for a gap layer [118]. In Grad-CAM, the weight $w_{k}^{c}$ is
not learned as in a GAP-CNN, but computed as the mean gradient of the score
class $S_{c}$ with respect to activation map $A^{k}_{x,y}$ of a layer-of-
interest (usually the last layer). Then Grad-CAM visualises features with
positive influence only (Eq. 3):

${\text{Grad-}CAM}^{c}={Re}{L}U\left(\sum_{k}w_{k}^{c}A_{x,y}^{k}\right)$ (3)

Finally, the CAM or Grad-CAM heatmap is up-sampled to the original input image
size and superimposed on the input image, which is why these heatmaps have a
coarse resolution.

Advantages and disadvantages: Grad-CAM is a popular method of
interpretability, both for natural images and medical images. It is most often
applied to image classification since the heatmaps are class-specific, but it
can also be applied to regression and segmentation tasks. Grad-CAM does not
require a modified cnn architecture; is not computationally intensive; is easy
to implement and widely available in multiple libraries. A disadvantage of CAM
and Grad-CAM is that the heatmaps are coarse (low resolution) because they are
often upsampled from the last convolutional layer of a network. To improve the
resolution, Grad-CAM has been coupled with other pixel-wise attribution
methods such as Guided Backpropagation, known as Guided Grad-CAM. In Guided
Grad-CAM, the Grad-CAM output is multiplied element-wise with the Guided
Backpropagation heatmap.

####  4.1.5 Weight analysis

An alternative approach for visualizing and explaining the decisions of a
network is to analyse the weights of the trained network. However, as deep
neural networks learn high-level features in the hidden layers, simply
visualizing the raw learned features usually does not offer human-
interpretable explanations [89]. Weight analysis methods attempt to create
human-understandable explanations through clustering weights and associating
clusters with human concepts.

Figure 7: Example of a Network Dissection model where the activation map of
individual filters in the network are analysed to identify which specific
concepts they have learnt by evaluating them against segmentation maps. Image
adapted from Bau et al. [11].

The Network Dissection approach [158] quantifies the interpretability of a cnn
by evaluating the alignment between activated regions of individual hidden
filters and human-labelled concepts (objects, parts, textures, colours). The
process involves first defining a set of task-relevant concepts and then
creating annotation masks $L_{c}(x)$ for each concept $c$ and image $x$. Next,
masks $M_{k}(x)$ of the top activated areas per filter $k$ and per image $x$
are created by scaling the activation maps $A_{k}$ to the size of the input
images, and binarizing them (thresholding on the top quantile level $T_{k}$ of
the distribution of pixel activations for filter $k$ over all images).
Finally, the accuracy of each filter $k$ in detecting concept $c$ is reported
as the sum of the Intersection over Union (IoU) between $M_{k}(x)$ and
$L_{c}(x)$ across all the images in a dataset (see Fig. 7). To quantify the
interpretability of a layer, the number of unique concepts aligned with
filters, i.e., unique detectors, are counted.

A graphical representation of the concepts learned by a network to understand
its behaviour was proposed [64]. This Concept Graphs framework involves
grouping similar weight vectors through hierarchical clustering in order to
define concepts. Then, formed weight clusters are associated with some region
in the input image by using a variation of Grad-CAM; the region corresponds to
a human-understandable concept, for example, a tumour boundary. After the
concepts have been identified, a concept graph is formed that represents the
link between concepts in different layers. This is computed by intervening on
the pairs of concepts and calculating the mutual information (MI) between pre-
interventional and post-interventional distributions as a measure of the link
between two concepts. The trails of concepts on the graph, therefore,
represent the flow of information used by the network when making a
prediction.

Several papers integrate Community Detection within a dl model for analysing
fmri data [30, 74]. The aim of Community Detection in the context of
neuroimaging is to discover $K$ networks of brain regions that are salient for
a particular dl task. Given an fmri connectivity matrix defined over $N$ brain
regions, the dl model incorporates a fully connected layer with a weight
matrix $W\in\mathbb{R}^{N\times K}$. Each value $w_{nk}\in W$ may be
interpreted as a membership score of brain region $n$ belonging to the
community $k$. A clustering algorithm is then applied to the weights to assign
brain regions to communities.

###  4.2 Intrinsic methods

Intrinsic interpretability refers to ml models that are explainable by design,
i.e., where feature representations can be understood by humans. The
interpretability can be due to the simple structure of the models, such as
short decision trees or sparse linear models, where network decisions can be
easily followed. Alternatively, interpretability can be achieved by explicitly
including interpretable modules or constraints in the model, as is required
for designing idl models. In this section, we present five categories of
intrinsic interpretable methods: disentangled latent spaces, interpretable
hybrid models and interpretable intermediate features, interpretable
generative models, deep structural causal models, and attention mechanisms.

####  4.2.1 Disentangled latent spaces

The latent space of a neural network is a learned representation of the input
data that has usually undergone compression, such that similar input samples
are transformed into representations that are close together in this space. A
popular dl model is the ae, where an encoder learns to compress input data to
a latent space, and a decoder learns to reconstruct the input data from the
latent representations. An extension to the ae that enables data generation
from the latent space is the vae, where the latent space is constrained to a
multivariate Gaussian distribution [61]. A desirable property is that the
latent space is to some extent disentangled, meaning a single factor in the
latent space corresponds to a single source of variation in the high-
dimensional image space. This can be encouraged through the introduction of
losses, which optimise for a subset of the latent space to encode specific
semantic feature(s) in the image space. This is illustrated in Fig 8, showing
a traditional vs. disentangled (for subject age) latent space of a trained
vae. Note that the structure of the latent space projection in 2D space is
more coherent for the disentangled space than it is for the traditional space.

Figure 8: Example of a traditional vs. a disentangled latent space of a
trained vae where age was added as a condition. The structured latent space
can, therefore, be used to generate new samples for a given condition (such as
age), as well as understand what type of changes occur in a given brain image
with age. Image adapted from Zhao et al. [156].

Capsule Networks are an alternative architecture to cnn that learn
disentangled, interpretable activation vectors [113]. Capsule Networks learn
spatial relationships between an object and its constituent parts, which are
invariant to the object viewpoint. Elements of an activation vector learn pose
parameters for an associated object, such as size, orientation, texture, and
hue. The $L_{2}$-norm (equal to the Euclidean distance from the origin) of an
activation vector is equal to the predicted probability of the corresponding
object, thus enabling classification.

Advantages and disadvantages: Disentangled latent representations provide some
control for image generation to the end user. The user can manipulate features
in the generated image in a semantically meaningful way by interpolating a
disentangled factor in the latent space. One limitation of disentangling
latent spaces for complex data is that the generative factors may not be
inherently independent, and by constraining the latent representation to have
independent representations, useful information about these dependencies can
be lost [83]. Additionally, constraining the latent space representations
often comes at the expense of performance [44]. One disadvantage of ProtoPNet,
in particular, is that distance maps are upsampled from the latent space to
the image space, which implicitly assumes that spatial relationships in the
image space are preserved in the latent space. However, Wolf et al. [144]
proved this is not necessarily the case, though efforts are being made to
account for this issue [144, 21].

####  4.2.2 Interpretable hybrid models and interpretable intermediate
features

A hybrid dl model usually has two components: a nn that learns intermediate
feature representations from the input data, coupled with a model that
predicts the learning task from the feature representations. The second
component can either be a nn or some other ml model, referred to as nn + nn
and nn + ml hybrid model, respectively.

An interpretable hybrid model is a hybrid model that possesses intermediate
feature representations that can be understood by humans and, therefore, act
as model explanations (see Figure 9). Some researchers also compute the
feature importance of the second model component and thus generate a second
set of model explanations along with the intermediate features [1, 69]. If a
study computes feature importance as a second set of explanations, we refer to
their approach as “int. features + feature importance”.

One notable example of an interpretable hybrid model is the prototypical part
network (ProtoPNet), which mimics human reasoning when classifying an image
[23]. The network learns a fixed number of prototypes for each class, where a
prototype is a tensor in latent space that is associated with an image patch
containing features typical of that class. At test time, latent features of an
image are compared to each prototype by computing a maximum similarity score,
and the similarity scores are passed through a fully connected layer to
predict the image class. Several studies in our review employed prototype
layers in their model architecture, inspired by ProtoPNet [88, 145, 92].

Figure 9: Example of an interpretable hybrid model where the intermediate
probability map is used as the features for an mlp model and acts as the model
explanation. The authors of this study [105] proposed a three-step approach:
first, a cnn classifier was trained to predict whether a given 3D brain mri
patch is ad or cn (1), then, the trained cnn produced a probability map (the
intermediate feature) for a given test subject (2), and, finally, an mlp was
trained on the intermediate feature probability map to distinguish between ad
and cn (3). Image adapted from Qiu et al. [105].

Advantages and disadvantages: An advantage of interpretable hybrid models is
that they may be designed so the intermediate features are suited for a
particular application. For example, in a clinical setting, intermediate
diagnostic features may be learned that are familiar to clinicians. However,
interpretable hybrid models require careful design and may take a long time to
develop.

####  4.2.3 Interpretable generative models

Another interpretability approach is to train a generative model to generate
explanations for neuroimaging tasks. The model learns to generate
modifications to the input image so that the modified image appears to belong
to a different class. The modifications are then used as explanations for the
prediction task. For example, in binary classification, the model $f$ modifies
an input image $x$ of class $0$, such that $x^{\prime}=f(x)$ appears to be
from class $1$. This task is often referred to as anomaly detection or
counterfactual generation, and the modifications are known as the anomaly map
or disease effect map. Such an example is shown in Fig. 10 where a network
produces the minimal additive mask needed to change an image from one class,
ad in this case, to another e.g., cn.

Figure 10: Example of a  interpretable generative model where a generator
network (shown as “CNN” in the figure) produces the minimum mask needed to
change the class of the input sample (from $c=1$ to $c=0$). A discriminator
network (“D”) is trained to distinguish between fake and real samples of the
same class c in order to constrain the generator to produce realistic samples.
Moreover, the masks can be used as explanations for the class discriminative
features. Image adapted from Baumgartner et al. [12].

Advantages and disadvantages: By learning to generate new images as
explanations for discriminative tasks, generative methods are capable of
capturing more meaningful class-discriminative features in comparison to
methods that evaluate the features learned by classification networks. These
generative methods also provide a framework to investigate how changing
features in an image e.g., by interpolation, affects the network decision.
However, generative models can be challenging to train and require high
computational power, rendering these methods harder to implement.

####  4.2.4 Deep structural causal models

Where randomised controlled trials are impossible, infeasible, or unethical,
estimating causal effects is often still possible using causal inference
methods. One such method is the scm, which estimates causal effects by
simulating population-level interventions [101]. An scm consists of a set of
$d$ endogenous variables $\\{X_{1},...,X_{d}\\}$, exogenous or noise variables
$\\{N_{1},...,N_{d}\\}$, and structural assignments (denoted as $:=$):

$X_{j}:=f_{j}(\boldsymbol{PA_{j}},N_{j}),j=1,...,d$ (4)

where $\boldsymbol{PA_{j}}\subseteq\\{X_{1},...,X_{d}\\}\setminus\\{X_{j}\\}$
are the parents of $X_{j}$. The joint probability distribution over the noise
variables is assumed to be jointly independent. An scm has an associated
causal graph $\mathcal{G}$ that visually represents our assumptions regarding
how data were generated in the real world. The causal graph $\mathcal{G}$ is a
directed acyclic graph (DAG) where all endogenous variables are represented as
nodes. A directed edge $X_{i}$ $\rightarrow$ $X_{j}$ exists in $\mathcal{G}$
if $X_{j}$ depends on $X_{i}$ for its value. Indeed, we define $X_{i}$ to be a
direct cause of $X_{j}$ if $X_{i}$ appears in the structural assignment
$f_{j}$ for $X_{j}$. Fig. 11 is an illustration of a causal graph for ms.

For a scm, the causal effect of intervening on a variable $X$ by setting it to
$a$ is denoted $do(X=a)$. It is also possible to estimate counterfactual
scenarios for specific individuals, which are hypothetical alternative
outcomes to the actual outcome. We refer readers to [102] and [103] for a
detailed overview of scm.

dscm employ neural networks to learn at least one of the structural
assignments in the scm, and applying them to medical imaging data is an
emerging research topic [22, 100, 108, 107]. Pawlowski et al. trained a dscm
on UK Biobank data to understand how a subject’s age ($a$), sex ($s$), brain
volume ($b$), and ventricle volume ($v$) influenced their brain mri image
($img$) [100]. The structural assignments were defined as (Eq. 5):

$\begin{split}b&:=f_{1}(a,s,N_{1})\\\ v&:=f_{2}(a,b,N_{2})\\\
img&:=f_{3}(v,b,N_{3})\end{split}$ (5)

where $N_{1},N_{2},N_{3}$ are noise variables. In this study, $f_{1}$ and
$f_{2}$ were modelled with normalising flows, and $f_{3}$ was learned using a
conditional vae where the vae generated an estimated brain MRI. A normalising
flow is a sequence of invertible transformations $g=f_{1}\circ
f_{2}\circ\text{...}\circ f_{K}$ that transforms a tractable distribution $z$
into a more complex distribution $x=g(z)$ (we refer readers to [62] for an
introduction to normalising flows). Images were generated for a variety of
counterfactual scenarios, and difference maps between the generated and
original images were visually inspected for interpretation. For instance, for
a 49-year-old subject, an image was generated for the counterfactual
$do(\text{age}=80\text{year-old})$; the generated image exhibited increased
ventricle volume and reduced brain volume compared to the original image,
consistent with trends in the true distribution.

Advantages and disadvantages: The strength of dscm is that causal mechanisms
of imaging markers may control for confounders, unlike most other dl models.
However, the causal graph $\mathcal{G}$ must be carefully constructed from
domain knowledge, and the structure of $\mathcal{G}$ may not yet be fully
elucidated. Furthermore, it is impossible to obtain ground truth data for
counterfactual scenarios, so counterfactual images cannot be validated.

Figure 11: Example of a causal graph where assumptions about the image
generating mechanism are explicitly defined. Deep structural causal models can
then be learned to estimate mri images under counterfactual scenarios. Image
adapted from Reinhold et al. [108] .

####  4.2.5 Attention mechanisms

In recent years, attention in the context of deep learning has become an
important area of research as it can be easily incorporated into existing
neural network architectures whilst also improving performance [18, 96] and
providing explanations [126, 115]. Attention methods learn a heatmap over the
inputs, features, or channels of the neural network, subsequently used to
weight the data to emphasise key features. In the following, we discuss four
main types of attention: channel, spatial, non-local, and self-attention,
which are illustrated in Fig. 12. For a more comprehensive description of dl
attention mechanisms, we refer the reader to Niu et al. [96] and Guo et al.
[39].

Channel attention assigns a weight to each filter in order to emphasize useful
features. One of the most popular channel attention blocks is the squeeze-and-
excitation block [48]. Let $U$ be a feature map with dimensions $H\times
W\times C$, then the squeeze-and-excitation block comprises a squeeze function
($F_{sq}$), which performs global average pooling [76], followed by an
excitation function ($F_{ex}$) defined as the sigmoid function ($\sigma$)
applied to an mlp. More specifically, the squeeze-and-excitation channel
attention $\alpha_{\mathbf{SE}}$ is defined as (Eq. 6):

$\begin{split}\alpha_{\mathbf{SE}}=\,&F_{ex}(\,F_{sq}(U)\,)\\\
=\,&\sigma\left(W_{2}\,\text{ReLU}\left(W_{1}\,\,\text{GlobalAvgPool}(U)\,\right)\right))\end{split}$
(6)

Here, $W_{1}\in\mathbb{R}^{\text{C/2}\times\text{C}}$ and
$W_{2}\in\mathbb{R}^{\text{C}\times\text{C/2}}$ are the weights of the fully
connected layers of the mlp. A different flavour of channel attention was
proposed by Woo et al. [147]. Here, both global max pooling and global average
pooling layers are included to generate two $C$-dimensional descriptors. The
final channel attention map $\alpha_{\mathbf{C}}$ is (Eq. 7):

$\begin{split}\alpha_{\mathbf{C}}=\sigma(\,&W_{2}\,\text{ReLU}\left(W_{1}\,\text{GlobalAvgPool}(U)\right)\,+\\\
&W_{2}\,\text{ReLU}\left(W_{1}\,\text{GlobalMaxPool}(U)\right)\,)\end{split}$
(7)

Spatial attention aims to extract important information in the image domain or
across the spatial dimensions of a feature map. In Woo et al. [147], the
spatial attention block first performs average and max pooling operations
across the channels of the input $U$, generating two feature maps which are
then concatenated. A convolutional layer is then applied to produce a
1-channel spatial map which, after passing through the sigmoid function,
becomes the attention map ($\alpha_{\mathbf{S}}$)(Eq. 8):

$\alpha_{\mathbf{S}}=\sigma\left(f\left([\text{AvgPool}(U);\text{MaxPool}(U)]\right)\right)$
(8)

where $[\cdot;\cdot]$ represents channel-wise concatenation and $f$ is the
convolutional layer.

Oktay et al. [98] introduced a different version of spatial attention with
Attention U-Net, where attention gates apply convolutions to both features
from the encoder and the corresponding decoder and then fuse them together to
create the attention map. Moreover, instead of simply concatenating the
encoder and decoder feature maps as for U-Net skip connections, the authors
first scaled the encoder features with the generated spatial attention.

Non-local attention, proposed by Wang et al. [140], aims to capture long-range
dependencies by computing interactions between any two positions in an image
or feature map. Conversely, channel or spatial attention focuses mainly on
local information, i.e., the pooling operation leads to loss of spatial
information, while convolutional layers process neighbourhood information.

In non-local attention, three parallel $1\times 1$ convolutional operations
($\theta$, $\phi$ and $g$) are applied on the input $U$, obtaining three
compressed feature maps, while a final $1\times 1$ convolutional operation $f$
restores the initial number of channels. Introduced by Lin et al. [76], the
$1\times 1$ convolutions act as a channel-wise pooling operator. The non-local
attention map $\alpha_{\mathbf{NL}}$ is obtained through the following
operations:

$\alpha_{\mathbf{NL}}=\tilde{X}=f\underbrace{\left(\text{softmax}\,\underbrace{\left(\theta(U)\otimes\phi(U)^{T}\right)}_{\beta_{1}\in\mathbb{R}^{\text{HW}\times\text{HW}}}\otimes\,g(U)\right)}_{\beta_{2}\in\mathbb{R}^{\text{HW}\times\text{C}/2}}$
(9)

where $T$ is the matrix transpose operation, $\otimes$ is the matrix
multiplication operator, and $\tilde{X}$, $\beta_{1}$ and $\beta_{2}$ are
shown in Fig. 12. Moreover, the outputs of the convolutional layers $\theta$,
$\phi$ and $g$ are reshaped to allow for matrix multiplication, i.e., they
become 2D matrices of shape $\text{HW}\times\text{C/2}$. These steps are shown
in the third panel of Fig. 12.

Finally, self-attention is a mechanism in deep learning, closely related to
the concept of non-local attention [140], commonly used in nlp tasks,
particularly in transformer-based architectures [138]. Dosovitskiy et al. [28]
adapted the self-attention model to image-based applications.

Figure 12: Example of attention–based methods showing channel attention in the
first panel (“CA”), spatial attention in the second panel (“SA”), non-local
attention in the third panel (“NLA”), and multi-head self-attention in the
fourth panel (“MHSA”). Note that “CA” and “SA” are adapted from Woo et al.
[147], “NLA” is adapted from Wang et al. [140], while the fourth panel is
adapted from Dosovitskiy et al. [28].

In architectures based on transformers [28, 41] which employ self-attention
modules, the initial step involves splitting the input data into a sequence of
patches. Subsequently, these patches undergo processing via a linear
projection layer and are merged with positional encodings to incorporate
spatial biasing within the patch sequence (see Fig. 12). The embedded patches
are then passed through a “transformer encoder” which consists of alternating
layers of multi-head self-attention and mlp blocks, as well as residual
connections and normalisation layers. More specifically, a multi-head self-
attention block is composed of multiple parallel self-attention heads, which
compute attention scores based on the query ($Q$), key ($K$), and value ($V$)
representations of the input as follows:

$\mathbf{SA}=\text{softmax}\left(\frac{QK^{T}}{\sqrt{D_{h}}}\right)V$ (10)

where $D_{h}$ is a scaling factor.

Advantages and disadvantages: Attention-based methods often add computational
complexity to an existing dl model, but have a differentiable objective and
are easily trainable with gradient descent. Moreover, they aim to provide a
weighting for inputs or internal features to focus the network on salient
characteristics. However, whether attention can be regarded as feature
importance is an ongoing debate [52, 119, 141].

##  5 Applications

###  5.1 Applications of perturbation-based methods (Table 3)

####  5.1.1 Neurodegenerative disease classification

Occlusion, Swap Test and Meaningful Perturbations have all been applied to ad
classification networks trained on brain mri images from the adni[131]. Liu et
al. [79] , Eitel et al. [31], and Yang et al. [150] employed Occlusion to
highlight image regions important for ad prediction. Furthermore, Rieke et al.
and and Yang et al. refined the Occlusion method by occluding brain regions
defined by an atlas, instead of image tiles [110, 150]. Nigri et al. suggested
Occlusion may be unsuitable for neuroimaging data, since an occluded patch in
a brain image from a cognitively normal individual may appear similar to
disease [95]. Consequently, Nigri et al. proposed and applied Swap Test [95].
Meaningful Perturbations with a constant-valued mask was also applied for ad
classification [134]. Furthermore, Shahamat et al. employed a complementary
approach to Meaningful Perturbations, in which minimal brain masks were
learned for cnns trained to classify ad and asd [120]. A minimal brain mask
keeps the fewest brain regions whilst still achieving high model accuracy,
whereas the brain mask in Meaningful Perturbations deletes the fewest brain
regions that cause a wrong prediction. In most of these ad studies,
explanations contained salient regions known to be altered in AD such as the
temporal lobe and hippocampus; however, the Occlusion map in one study could
not be meaningfully interpreted because the occlusion window was too large
[31].

Another perturbation method, LIME, was applied to explain predictions of a VGG
[125] model trained to predict pd on spect data [81]. LIME explanations of cn
individuals clearly delineated the putamen and caudate regions, whereas the
highlighted areas in explanations for pd patients often extended beyond these
regions.

####  5.1.2 Autism spectrum disorder classification

Global perturbation-based methods have been employed to identify important
features for asd classification [85, 73]. Permutation Feature Importance was
applied to a model trained on both structural (e.g., cortical volume and
thickness) and functional mri features [85]. Li et al. modified Occlusion to
produce a global explanation for asd classification [73]. Specifically, after
a cnn was trained to classify asd, each atlas-based brain region $r$ was
occluded in all images. Let $P_{\text{orig}}^{\text{\ac{cn}}}$ and
$P_{\text{orig}}^{\text{\ac{asd}}}$ represent the distribution of predicted
probabilities for all control and asd subjects, respectively. Similarly, let
$P_{\text{occ}}^{\text{\ac{cn}}}$ and $P_{\text{occ}}^{\text{\ac{asd}}}$
represent the distributions for the corresponding images with region $r$
occluded. Then, utilising Jensen-Shannon Divergence ($JSD$), the distance
between class distributions was computed and compared (Eq. 11):

$\text{JSD}(P_{\text{orig}}^{\text{\ac{cn}}},P_{\text{orig}}^{\text{\ac{asd}}})>\text{JSD}(P_{\text{occ}}^{\text{\ac{cn}}},P_{\text{occ}}^{\text{\ac{asd}}})$
(11)

Brain region $r$ was considered important if the decrease in JSD for occluding
that region was statistically significant. The assumption is if region $r$ is
important for asd prediction, then the cnn will not separate classes as
effectively when region $r$ is removed.

Approaches that learn optimal brain masks for asd classification have been
used for models trained on the abide [27] I dataset [26, 120]. Dhurandhar et
al. produced the pertinent positive $\delta_{pos}$ and negative $\delta_{neg}$
features for a given resting-state fmri image $I_{0}$ [26]. The former,
$\delta_{pos}$, are minimally sufficient meaning the network will predict the
same class for both $I_{0}$ and $\delta_{pos}$. In contrast, $\delta_{neg}$
must be absent for the prediction i.e. the network predicts a different class
for the perturbed image $I_{0}+\delta_{neg}$ compared to $I_{0}$. Similarly,
Shahamat et al. learned minimal brain masks for structural mri images from
abide I, using the same approach they applied to the adni dataset, previously
described in 5.1.1 [120]. The important brain regions identified in all these
asd studies were disparate with little overlap, although motor regions were
frequently highlighted.

####  5.1.3 Schizophrenia classification

A method analogous to Occlusion was applied to a rnn trained to classify scz
on rsfmri data [149]. Initially, the data underwent dimensionality reduction
using independent component analysis [20], and the time series of 50 valid ic
were retained. After rnn training, feature importance of the $i^{th}$ ic was
computed by replacing the $i^{th}$ ic time series values with its average
(essentially occluding the $i^{th}$ feature), and then the change in model
performance was assessed. The ic with the greatest change in performance were
considered to be the most important features for classification, which were
located in the dorsal striatum and cerebellum.

####  5.1.4 Brain age regression

The U-noise method, similar to Meaningful Perturbations, trained a U-net [111]
by adding maximum random noise to input images without affecting the
performance of a pre-trained prediction model [63]. A sensitivity map can be
generated to show the image pixels that were least tolerant to the addition of
noise. Bintsi et al. adapted the U-noise architecture to interpret a 3D ResNet
[42] trained for brain age regression on UK Biobank [129] $T_{1}$-weighted mri
data [16]. They computed an average importance map across all subjects for the
test dataset. However, changes in the cerebral cortex related to aging were
not well captured by the network.

####  5.1.5 Sex classification

Meaningful Perturbations was compared to two alternative idl methods for
visualising model decisions for sex classification [54]. The Meaningful
Perturbations explanation highlighted regions of the frontal lobe, though the
explanations were visually dissimilar between the three methods.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
AD Classification
Eitel et al. [31] | ADNI | sMRI (3D) | 344 | Occlusion
Liu et al. [79] | ADNI | DTI (3D) | 151 | Occlusion
Yang et al. [150] | ADNI | sMRI (3D) | 103 | Occlusion
Nigri et al. [95] | ADNI+AIBL | sMRI (3D) | 1,248 | Swap Test
Thibeau et al. [134] | ADNI+AIBL | sMRI (3D) | 1,171 | Meaningful Pert.
Shahamat et al. [120] | ADNI | sMRI (3D) | 140 | Optimal mask
Tang et al. [132] | In-house | Histology (2D) | 33 | Occlusion
PD Classification
Magesh et al. [81] | PPMI | SPECT (2D) | 642 | LIME
ASD Classification
Li et al. [73] |  Various† | fMRI (3D) | 225 | Occlusion (global)
Dhurandhar et al. [26] | ABIDE | fMRI (ts) | 293 | Optimal mask
Shahamat et al. [120] | ABIDE | sMRI (3D) | 1,000 | Optimal mask
Mellema et al. [85] | IMPAC | sMRI/fMRI (3D) | 915 | PFI
SCZ Classification
Yan et al. [149] | In-house | fMRI (3D) | 1,100 | Leave-one-IC-out
Sex Classification
Kan et al. [54] | HCP | sMRI (3D) | 1,113 | Meaningful Pert.
Brain Age Regression
Bintsi et al. [16] | UK Biobank | sMRI (3D) | 13,750 | U-Noise
Table 3: Articles using perturbation-based interpretable methods.
Abbreviations: ABIDE, Autism Brain Imaging Data Exchange; AD, Alzheimer’s
Disease; ADNI, Alzheimer’s Disease Neuroimaging Initiative; AIBL, Australian
Imaging Biomarker and Lifestyle Flagship Study of Ageing; ASD, Autism Spectrum
Disorder; DTI, Diffusion Tensor Imaging; fMRI, functional Magnetic Resonance
Imaging; HCP, Human Connectome Project; IMPAC, IMaging-PsychiAtry Challenge;
LIME, local interpretable model-agnostic explanations; PD, Parkinson’s
Disease; PPMI, Parkinson’s Progressive Markers Initiative database; SCZ,
schizophrenia spectrum disorders; sMRI, structural Magnetic Resonance Imaging;
SPECT, Single Photon Emission Computed Tomography. Datasets: $\dagger$In-
house+ABIDE

###  5.2 Applications of gradient-based methods (Table 4)

####  5.2.1 Neurodegenerative disease classification

Several studies employed Vanilla Gradients or Grad $\times$ Input to identify
important brain regions after training a cnn on mri data from adni for ad
classification [97, 32, 31, 110]. All studies followed a similar approach,
where class-average gradient maps were computed and then compared to
explanations from other methods. Oh et al. produced Vanilla Gradients maps for
ad classification and found they were in agreement with Occlusion maps [97].
Essemlali et al. focused on connectivity between brain regions from diffusion-
weighted mri data, and found Vanilla Gradients highlighted relevant brain
regions, unlike Occlusion [32]. Eitel et al. analysed robustness across model
runs for Grad $\times$ Input, Occlusion and two additional methods, and
concluded Grad $\times$ Input produced the least consistent explanations [31].
To aid interpretation of Vanilla Gradients maps, Rieke et al. computed a
quantitative relevance score by summing sensitivity within each atlas-based
brain region [110]. Eitel et al. similarly explored three summary statistics
across brain regions, including mean region sensitivity to account for
differences in brain region volume. All studies determined the medial temporal
lobe and/or the hippocampus as important regions for ad classification.

####  5.2.2 Autism spectrum disorders classification

Vanilla Gradients has been applied to identify important features for asd
classification from task-based fmri data, where the task was testing
perception of people’s movements (biopoint task) [72]. A graph was constructed
with each node corresponding to a specific brain region, and having an
associated feature vector of ten researcher-selected features. After training
a gnn and computing the gradient of the network output with respect to each
feature, sensitivity maps were averaged across nodes and subjects to generate
a sensitivity score per feature.

####  5.2.3 Human immunodeficiency virus classification

Vanilla Gradients maps may include the influence of a confounding factor on
the model decision, for example, patient age is often a confounding factor for
neurodegeneration. Zhao et al. modified Vanilla Gradients to remove the
influence of age as a confounder from sensitivity maps computed from a cnn
classifier for hiv [157]. Let $\text{f}_{j}=(f_{j1},f_{j2},...,f_{jN})$ be the
$j^{th}$ feature from the final convolutional layer for all $N$ subjects,
$\textbf{s}=(s_{1},...,s_{N})$ be cnn score and $\textbf{a}=(a_{1},...,a_{N})$
be subject age. Then the linear model
$\text{f}_{j}=\beta_{0}+\beta_{1}\textbf{s}+\beta_{2}\textbf{a}$ was fitted
and if $\beta_{2}$ was non-zero, then age was defined as a confounder for the
$j^{th}$ feature. When calculating the Vanilla Gradients map, gradients were
computed for unconfounded features only. The confounder-free sensitivity maps
showed the posterior ventricle was most influenced by age rather than hiv.

####  5.2.4 Brain age regression

SmoothGrad was applied to a cnn trained on a $T_{1}$-weighted brain mri
dataset to predict subject age, and a population-average sensitivity map was
computed [71]. The ventricles and subarachnoid cisterns were predominantly
highlighted in the sensitivity map, which may be related to brain atrophy from
the aging process.

####  5.2.5 Cognitive task decoding

McClure et al. trained a cnn on activation maps from task-based fmri to
classify the fmri task and employed Vanilla Gradients for model interpretation
[84]. To address the shattered gradients limitation of Vanilla Gradients the
authors utilised adversarial training. More specifically, non-targeted
adversarial noise was learned and added to each image, optimised as the
smallest magnitude noise that minimised the probability of the correct class.
In this way, coupling Vanilla Gradients with adversarial training was found to
produce gradient maps that were more class discriminative than Vanilla
Gradients, Grad $\times$ Input, and SmoothGrad. However, the maps were still
only found to be weakly correlated with class-specific features.

lstm is a deep learning architecture that is well-suited to time series fmri
data because it is designed to process sequence data. However, one limitation
of applying Vanilla Gradients to lstm models is the issue of vanishing
gradients [13] when backpropagating through many timesteps. Consequently, only
features in the latest time steps are highlighted in gradient maps. Ismail et
al. proposed incorporating an attention mechanism into an lstm, to bypass
backpropagating through multiple timesteps during Vanilla Gradients [51]. The
attention-based lstm was trained on task-based fmri to classify the fmri task,
and Vanilla Gradients was then able to highlight features in early time steps.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
AD classification
Oh et al. [97] | ADNI | sMRI (3D) | 694 | Vanilla Gradients
Eitel et al. [31] | ADNI | sMRI (3D) | 344 |  Grad $\times$ Input
Essemlali et al. [32] | ADNI | DW-MRI (2D) | 186 | Vanilla Gradients
ASD classification
Li et al. [72] | In-house | fMRI (ts) | 118 | Vanilla Gradients
HIV classification
Zhao et al. [157] | In-house | sMRI (3D) | 355 | Vanilla Gradients
Brain age regression
Levakov et al. [71] |  Various† | sMRI (3D) | 10,176 | SmoothGrad
Cognitive task decoding
McClure et al. [84] | HCP | fMRI (3D) | 965 | Vanilla Gradients
Ismail et al. [51] | HCP | fMRI (ts) | 749 | Vanilla Gradients
Table 4: Articles using gradient-based methods. Abbreviations: ABIDE, Autism
Brain Imaging Data Exchange; AD, Alzheimer’s Disease; ADNI, Alzheimer’s
Disease Neuroimaging Initiative; AIBL, Australian Imaging Biomarker and
Lifestyle Flagship Study of Ageing; ASD, Autism Spectrum Disorder; DW-MRI,
Diffusion-Weighted Magnetic Resonance Imaging; fMRI, functional Magnetic
Resonance Imaging; HCP, Human Connectome Project; HIV, Human Immunodeficiency
Virus; IXI, Information Extraction from Images; sMRI, structural Magnetic
Resonance Imaging; ts, time series. Datasets:
$\dagger$=ABIDE+ADNI+AIBL+IXI+others

###  5.3 Applications of backpropagation-based methods (Table 5)

####  5.3.1 Alzheimer’s disease classification

lrp and Guided Backpropagation were utilised for ad classification after
training on adni structural mri data [31, 17]. In one study, lrp heatmaps were
shown to be more class-discriminative than Guided Backpropagation maps [17].
In a similar approach to the Eitel et al. study [31], the heatmap analysis was
improved by three summary statistics: sum of relevance, mean relevance (to
account for brain region size) and relevance gain compared to cn (to find
regions where explanations between ad and cn differ the most). All three
studies consistently identified the hippocampi and other structures in the
temporal lobe as important for ad classification.

####  5.3.2 Cognitive task decoding

lrp-$\epsilon$ has been coupled with a deep learning model trained on task-
based fMRI data to predict one of four cognitive states associated with
viewing four image categories (body, face, place, or tool) during the task
[135]. The population-level explanation for each cognitive state was compared
against a meta-analysis associated with the keyword from NeuroSynth. The
explanations generally matched with the meta-analysis for body and face
cognitive states, but less so for place and tool.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
AD classification
Bohle et al. [17] | ADNI | sMRI (3D) | 344 |  lrp-$\alpha\beta$
Eitel et al. [31] | ADNI | sMRI (3D) | 344 |  Guided backprop + lrp
Sex classification
Kan et al. [54] | HCP | sMRI (3D) | 1,113 | Guided backprop
Cognitive task decoding
Thomas et al. [135] | HCP | fMRI (ts) | 100 |  lrp-$\epsilon$
Table 5: Articles using backpropagation-based methods. Abbreviations:
Alzheimer’s Disease Neuroimaging Initiative; HCP, Human Connectome Project;
LRP, Layer-wise Relevance Propagation; sMRI, structural Magnetic Resonance
Imaging; fMRI, functional magnetic resonance imaging; ts, time series. Note:
where lrp variant not provided in paper, refer to as ’lrp’.

###  5.4 Applications of Class Activation Maps (Table 6)

####  5.4.1 Neurodegenerative disease classification

CAM and Grad-CAM have been applied to ad and mci classification using a ResNet
for 2D mri [57], VGG and ResNet for 3D mri [153, 150], and gcn for surface
meshes of the cortex and sub-cortical structures [5], all trained on adni
structural mri data. Additionally, Tang et al. applied Guided Grad-CAM to
amyloid-beta (A$\beta$) plaque-stained immunohistochemical data to classify
plaque morphology, since A$\beta$ plaques are a histopathological hallmark of
ad [132].
In their study of AD classification using a VGG, Zhang et al. showed that
applying Grad-CAM to lower convolutional layers produced more detailed
explanations [153]. However, as lower layers tend to respond to
edges/junctions of the brain images, so did the corresponding Grad-CAM maps.

A related application of network interpretation is to use it to diagnose
failure cases, for example, Khan et al. evaluated CAM for a case of failed
classification [57]. In this case, the network attended to structures that are
not associated with ad classification, such as the skull. Similarly,
Williamson et al. identified failure cases using Grad-CAM for pd
classification of spect scans, where the presence of noise artifacts and
hyperintensities were shown to influence the network decision [142].

####  5.4.2 Intracerebral hemorrhage classification

Grad-CAM generated explanations for a dl model that detected and classified
ich sub-types from ct scans of the head [70]. Ground truth data was available
to validate the explanations; the proportion of “bleeding points”, selected by
neuroradiologists to indicate the centre of haemorrhagic lesions, overlapping
Grad-CAM heatmaps was 78%.

####  5.4.3 Brain tumour classification

Windisch et al. used Grad-CAM to explain tumor classification from structural
and diffusion MRI data [143]. As in the Khan et al. study [57], results were
visually evaluated for cases of correct and incorrect classification. The
network focused on the tumor when correctly classifying the scans, whilst
there was no clear attention pattern when the classification failed.

####  5.4.4 Autism spectrum disorder classification

Li et al. proposed visualisation of activation maps from a gcn trained to
predict asd from task-based fmri data [75]. Using an approach analogous to CAM
but for gcn, the 25% of graph nodes (representing 21 brain regions) with the
highest activation scores after the final graph-convolutional layer were
visualised to interpret model classification. The method highlighted the
dorsal striatum, thalamus, and frontal gyrus, regions thought to be affected
by asd.

####  5.4.5 Sex classification

The Dense-CAM network proposed by Gao et al. employed CAM in the final layer
of a DenseNet [50] trained for sex classification [35], and found the
cerebellum to be the most important brain region. On the other hand, Kim et
al. [58] applied Grad-CAM to a gcn trained on rsfmri data and found regions
involved in the default mode network to be important for sex classification,
but not the cerebellum.

####  5.4.6 Cognitive score prediction

In two studies, Grad-CAM was applied to visualise important brain regions for
predicting the wrat score of healthy individuals from the pnc, using a
regression [106] and classification approach [49]. In one study, a gcn was
trained on task-based fmri data, then Grad-CAM adapted for regression was
computed [106]. Alternatively, subjects were classified into low, medium, and
high wrat score, and Guided Grad-CAM was used to identify important brain
regions for wrat score classification [49]. Both studies identified regions of
the occipital lobe as important, which is involved in object recognition.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
AD classification
Zhang et al. [153] | ADNI | sMRI (3D) | 1,407 | Grad-CAM
Yang et al. [150] | ADNI | sMRI (3D) | 103 | CAM + Grad-CAM
Azcona et al. [5] | ADNI | Surface mesh (3D) | 435 | Grad-CAM
Khan et al. [57] | ADNI | sMRI (2D) | 150 | CAM
Tang et al. [132]† | In-house | Histology (2D) | 33 | Guided Grad-CAM
PD classification
Williamson et al. [142] | PPMI | SPECT (3D) | 600 | Grad-CAM
ASD classification
Li et al. [75] | In-house | fMRI (ts) | 118 | Activation maps
Tumour classification
Windisch et al. [143] |  Various‡ | sMRI (2D) | 2, 479 | Grad-CAM
Tumour segmentation
Natekar et al. [93] | BraTS | sMRI (3D) | 461 | Grad-CAM
ICH classification
Lee et al. [70] | In-house | CT (2D) | 904 | CAM
Sex classification
Gao et al. [35] |  Various†† | sMRI (3D) | 6,008 | CAM
Kim et al. [58] | HCP | fMRI (ts) | 1,094 | CAM + Grad-CAM
Kan et al. [54] | HCP | sMRI (3D) | 1,113 | Grad-CAM
Cognitive score prediction
Hu et al. [49] | PNC | fMRI (ts) | 854 | Guided Grad-CAM
Qu et al. [106] | PNC | fMRI (ts) | 800 | Grad-RAM
Table 6: Articles using CAM interpretable methods. Abbreviations: AD,
Alzheimer’s Disease; ADNI, Alzheimer’s Disease Neuroimaging Initiative; ASD,
Autism Spectrum Disorder; BraTS, Brain Tumor Segmentation challenge; CAM,
Class Activation Map; CoRR, Consortium for Reliability and Reproducibility;
CT, Computed Tomography; FCP, Functional Connectome Project; fMRI, functional
Magnetic Resonance Imaging; GSP, Brain Genomics Superstruct Project; HCP,
Human Connectome Project; ICH, Intracerebral Hemorrhage; IXI, Information
Extraction from Images; NKI-RS, Nathan Kline Institute-Rockland Sample; PD,
Parkinson’s Disease; PNC, Philadelphia Neurodevelopmental Cohort; PPMI,
Parkinson’s Progressive Markers Initiative database; RAM, Recurrent Activation
Map; fMRI, functional Magnetic Resonance Imaging; SLIM, Southwest University
Longitudinal Imaging Multimodal; sMRI, structural Magnetic Resonance Imaging;
SPECT, Single Photon Emission Computed Tomography; TCGA, The Cancer Genome
Atlas; ts, time series. $\dagger$= A$\beta$ plaque morphology classification.
$\ddagger$= IXI, CyberKnife, TCGA. $\dagger\dagger$= HCP, FCP, GSP, NKI-RS,
CoRR, SLIM.

###  5.5 Applications of weight analysis (Table 7)

####  5.5.1 Tumor segmentation

The Concept Graphs framework was applied to a U-Net brain tumor segmentation
model trained on the brats dataset [7]. The method identified multiple
concepts at various model layers, such as the whole tumour, tumour core
boundaries, and the tumour core region. Concept detection was also used for
interpretability of a U-Net tumor segmentation model in Natekar et al. using
Network Dissection [93]. Results showed individual filters learned
interpretable concepts including grey and white matter and edema; and separate
filters for the whole tumor and the tumor core. These results showed that
segmentation networks exhibit a modularity in the inference process that can
be understood by humans. In Kori et al. [64], in collaboration with a
radiologist, inference trails that represent the trail of information in the
network were also analysed. The network was shown to take a hierarchical
approach to segmentation, starting with the detection of edges at lower layers
and moving to the detection of the tumor in upper layers.

####  5.5.2 ASD classification

Dvornek et al. incorporated Community Detection within their dl model trained
for asd classification on rsfmri data [30]. The weights
$W\in\mathbb{R}^{N\times K}$ for Community Detection, where $w_{nk}\in W$
represents the strength of the connection between brain region $n$ and
community $k$, were learned as part of the model. Clustering was then
performed for each $k^{th}$ community vector $[w_{1k},w_{2k},...,w_{Nk}]$ to
assign each brain region as belonging vs. not belonging to community $k$.
Finally, the importance of community $k$ for asd classification was defined as
the sum of absolute weights of all $k$-indexed nodes in the classification
model. The three most important communities included brain regions associated
with language and social processing, memory, and reward-processing and
decision-making.

####  5.5.3 Cognitive task decoding

gcn were re-designed for Community Detection in BrainGNN [74]. Let
$W_{i}^{(l)}$ denote the learnable weights associated with node $i$ in graph
convolutional layer $l$ of a GNN, where node $i$ represents a fixed brain
region with one-hot location encoding $n_{i}$. The authors proposed to encode
brain region location in $W_{i}^{(l)}$ by training a mlp on the brain region
location $n_{i}$:

$vec(W_{i}^{(l)})=\Theta_{2}^{(l)}\text{\ac{relu}}(\Theta_{1}^{(l)}n_{i})+b^{(l)}$
(12)

where $\Theta_{1}^{(l)}$, $\Theta_{2}^{(l)}$ and $b^{(l)}$ are mlp parameters.
The elements $(\alpha_{nk})^{+}$ of
$\text{\ac{relu}}(\Theta_{1}^{(l)})\in\mathbb{R}^{N\times K}$ were interpreted
as the non-negative community detection scores of brain region $n$ belonging
to community $k$. In this study, BrainGNN was trained on the biopoint task-
based fmri dataset for asd classification, as well as on the hcp task-based
fmri data to classify seven cognitive tasks.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
ASD classification
Dvornek et al. [30] | ABIDE | fMRI (ts) | 527 | Community detection
Li et al. [74] | Biopoint | fMRI (ts) | 115 | Community detection
Tumour segmentation
Kori et al. [64] | BraTS | sMRI (3D) | 300 | Concept Graphs
Natekar et al. [93] | BraTS | sMRI (2D) | 461 | Network Dissection
Cognitive task decoding
Li et al. [74] | HCP | fMRI (ts) | 237 | Community detection
Table 7: Articles using weight analysis interpretable methods. Abbreviations:
ABIDE, Autism Brain Imaging Data Exchange; ASD, Autism Spectrum Disorder;
BraTS, Brain Tumor Segmentation challenge; fMRI, functional Magnetic Resonance
Imaging; HCP, Human Connectome Project; sMRI, structural Magnetic Resonance
Imaging.

###  5.6 Applications of disentangled latent spaces (Table 8)

####  5.6.1 Tumour classification

A Capsule Network was trained to classify tumour type (meningioma, pituitary,
glioma) from segmented brain mri [3]. The Capsule Network learned to
reconstruct input images, where the latent space was constrained to three
activation vectors representing the three tumour types. The activation vectors
were inspected by perturbing individual vector elements and visualising the
reconstructed images, which revealed the Capsule Network had learned
interpretable tumour features such as tumour size and elongation.

####  5.6.2 Image generation

One application of latent space disentanglement is training an ae where a
latent factor represents age, enabling the generation of mri for different
ages.[156, 91]. One study jointly trained a supervised age regression network
and a vae where both shared convolutional layers [156]. The latent space of
the vae was trained to approximate a prior distribution $p(z|\hat{y})$
conditioned on the age $\hat{y}$ predicted by the regressor. In another study,
an ae was coupled with a linear function such that the first parameter of the
latent space predicted subject age [91]. Age-specific mri images were
synthesised in both studies by adjusting the age-related latent factor.
One study synthesised T1w mri images at different ages for healthy controls
and for patients with ad, through disentangling the effect of ad from healthy
ageing on mriduring image reconstruction [99]. Disentanglement was achieved by
learning two orthogonal directions in the latent space of an ae and designing
loss functions that encouraged the difference between two successive mri of a
subject, as represented in the latent space, to be equal to the sum of two
components in the healthy and diseased directions. Formally, let $z_{1}$ and
$z_{2}$ be the latent representations of two mri images acquired for a subject
at times $t_{1}$ and $t_{2}$, respectively, with $t_{2}>t_{1}$. Then the
vector $\Delta z=(z_{2}-z_{1})/(t_{2}-t_{1})$ was constrained by the loss so
that $\Delta z\approx\Delta z_{a}+\Delta z_{d}$, where $\Delta z_{a}$ and
$\Delta z_{d}$ are parallel to the healthy and disease directions,
respectively. The disentangled latent factors were visualised and the ae had
learned distinct trajectories for cn, pmci, smci and ad subjects.

Another application is harmonisation of mri across different clinical sites,
where the contrast of an mri image acquired at site A is transformed to appear
as if it were acquired at site B, whilst leaving the subject anatomy
unaltered. Disentangled latent spaces have been employed in mri harmonisation
models to learn separate features for anatomy vs. image contrast [160, 154].
For example, Zuo et al. trained a model named CALAMITI for site harmonisation
across 10 clinical sites [160]. Two separate encoders, $E_{anat}$ and
$E_{cont}$, learned feature representations for anatomical and contrast
information during image reconstruction, respectively $\theta_{anat}$ and
$\theta_{cont}$. $E_{anat}$ learned anatomical-only features by encouraging
$\theta_{anat}$ from both T1-w and T2-w images of the same subject and slice
number to be identical, thus capturing the shared anatomical features and
ignoring contrast. Simultaneously, $E_{cont}$ was trained on different slice
numbers of the same T1-w and T2-w images to represent the residual features of
contrast in $\theta_{cont}$, after accounting for anatomy in $\theta_{anat}$.
Similarly, Zhao et al. learned a disentangled vae latent space consisting of
site-related vs. non site-related features [154]. For $M$ sites, a vector of
length $M$ in the latent space was optimised to represent site-specific
features by feeding it into a site classification network and minimising the
cross-entropy loss. Both studies visualised the latent space and demonstrated
subjects were clustered by site.

####  5.6.3 Brain age regression

Disentanglement of a network latent space may be advantageous for networks
trained on multi-modal data, to uncouple modality-specific from shared latent
features [47]. Such an approach was adopted when predicting infant brain age
from both fmri and smri features, where an ae was trained on each modality
such that the latent space was divided into modality-specific features and
common features [47]. The common feature vectors from both ae for the same
subject were encouraged to be identical by adding an $L_{2}$ loss and
adversarial loss. Furthermore, each decoder was required to reconstruct the
input data from its own common feature vector as well as that from the ae of
the other modality, a method known as cross reconstruction. The common feature
vector and each modality-specific feature vector were passed to an age
prediction network to predict age. Visualisation of the learned latent space
demonstrated the latent features were ordered by age.

####  5.6.4 Neurodegenerative disease classification

Similarly, Wang et al. learned a disentangled latent space to partition
features by modality, as well as learn features shared between modalities,
when learning from multimodal data [139]. Two ae were trained, one on each
modality, on smri/fdgpet for ad classification from adni, as well as on
smri/dti for pd classification from ppmi. Each latent space was split into
modality-specific and common features; a distance loss, defined as the $L_{2}$
loss between common features divided by the $L_{2}$ loss between modality-
specific features, encouraged common features to be identical and modality-
specific features to be dissimilar. In addition, cross-reconstruction was
adopted, where each decoder generated images using common features learned
from both ae. Disentanglement enabled important brain regions to be identified
as smri-specific, fdgpet (or dti in the case of pd classification)-specific or
common to both imaging modalities.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
Image generation
Zhao et al. [156] | In-house | sMRI (3D) | 245 | Factor: age
Mouches et al. [91] | In-house | sMRI (3D) | 2,681 | Factor: age
| +IXI |  |  |
Zuo et al. [160] |  Various† | sMRI (2D) | 100 | Factor: clinical site
Zhao et al. [154] |  Various‡ | Surface mesh (3D) | 2,542 | Factor: clinical site
Ouyang et al. [99] | ADNI | sMRI (3D) | 632 | Factor: healthy+disease
Tumour classification
Afshar et al. [3] | In-house | sMRI (2D) | 233 | Capsule network
Brain age regression
Hu et al. [47] | BCP | sMRI+fMRI (3D) | 178 | Factor: modality
Neurodegenerative disease classification
Wang et al. [139] |  +ADNI | sMRI+PET (3D) | 4,126 | Factor: modality
| +PPMI | DTI (3D) |  |
Table 8: Articles using disentangled latent space methods. Abbreviations:
ADNI, Alzheimer’s Disease Neuroimaging Initiative; BCP, UNC/UMN Baby
Connectome Project; BLSA, Baltimore Longitudinal Study of Aging; BraTS, Brain
Tumor Segmentation Challenge; fMRI, Functional Magnetic Resonance Imaging;
IBIS, Infant Brain Imaging Study; IXI, Information Extraction from Images;
OASIS, Open Access Series of Imaging Studies; sMRI, Structural Magnetic
Resonance Imaging. Datasets: $\dagger$=IXI+OASIS+BLSA, $\ddagger$=In-
house+BCP+IBIS

###  5.7 Applications of interpretable hybrid models (Table 9)

####  5.7.1 Neurodegenerative disease classification

One blueprint for interpretable intermediate features for ad classification is
a heatmap of predicted probabilities of ad across brain regions. Qiu et al.
adopted this approach when designing a nn + nn hybrid model, where the first
component was a patch-based cnn that output probability of ad across brain
patches [105]. After training, predicted probabilities for 200 voxels were
concatenated with non-imaging features (age, gender, mmse) and used to train a
multi-layer perceptron to predict ad (summarised in Fig. 9). However, the
heatmaps were less precise and therefore more difficult to interpret than the
next two studies discussed. Similarly, a nn+ml hybrid model learned
intermediate probability heatmaps for ad [69]. For the first component, an
ensemble of nn classifiers was trained to predict ad or mci status for each of
93 atlas-derived brain regions, from which a probability heatmap was
constructed. The second component was a linear svm trained to predict ad
status from the probability heatmap. This study followed the “int. features +
feature importance” approach and also considered the weights learned by the
linear svm. Nguyen et al. adopted a comparable approach where instead of
learning a brain heatmap of probabilities, voxels were assigned a grade close
to 1 if abnormal and close to -1 if healthy [94]. The NN + NN hybrid model
consisted of a set of patch-based U-Nets that generated the grade heatmap,
followed by a gcn trained to predict cn vs. ad vs. ftd. The population-average
heatmaps in all three studies were highly class-discriminative and were
consistent with known disease pathology, focusing predominantly on the
temporal lobe for ad, and the frontal lobe for ftd in the case of the study by
Nguyen et al..

Three studies employed prototypes to introduce interpretability into ad
classification models [88, 145, 92]. ProtoPNet has been trained on two public
T1-w mri datasets,adni and oasis [82], to classify ad [88]. Wolf et al.
trained a variant of ProtoPNet to predict ad from fdgpet images from the adni
database [145], and two of the prototypes highlighted the ventricles and
occipital lobe. In another study, a prototype-based model was trained on T1-w
mri from adni and an in-house dataset, and the prototypes learned were
reconstructed from the latent space to a 3D feature vector [ad diagnosis,
mmse, age] [92]. Furthermore, prototypical brains for each diagnosis class
(cn, progressive mci, stable mci and ad) were compared to individual scans;
for example, a cn subject differed most from the ad prototypical brain in the
amygdalae. However, a limitation of prototype-based explanations is they are
low-resolution because of upsampling from a latent space to the image space.

Abuhmed et al. explored yet another hybrid model design and predicted ad
clinical scores as interpretable intermediate features[1]. The nn+ml hybrid
model in the predicted prognosis of ad status at month 48 from multi-modal
data collected at baseline and months 6, 12, and 18. The first component was a
Bi-lstm model trained to forecast seven cognitive scores (such as the mmse) at
month 48. The second component was an interpretable ml classifier trained
separately to forecast disease status at month 48 from the seven forecasted
cognitive scores, subject age, gender, and education. Taking the “int.
features + feature importance” approach, explanations were also generated from
the ml classifier; however, the explanations were only in relation to the
cognitive scores and not the neuroimaging data.

####  5.7.2 Autism spectrum disorder classification

A prototype-based approach was adopted when classifying asd from fc matrices
computed from rsfmri data from the abide dataset [55]. The first component of
the hybrid model was a transformer network (see Section 4.2.5) that generated
latent features for a subject, and then the predicted class was determined by
the class prototype closest to the latent features. To enable
interpretability, a decoder was trained to reconstruct the input fc from an
individual’s latent features and was also used to decode a more prototypical
fc. For example, by subtracting a reconstructed fc of a control subject from
the asd-typical version, the authors found regions such as the right cingulate
gyrus and the occipital and frontal poles as the most different from asd for
this individual.

####  5.7.3 Brain age regression

Prototypes were also employed for predicting brain age from T1-w mri
images[cite IXI], as well as predicted gestational age from fetal us images
[43]. Prototypes were adapted for regression as follows: prototypes were not
assigned to a class, but each prototype was replaced with the closest latent
representation of a training image, and associated with the corresponding age
label. The predicted age is the weighted mean of age labels of all prototypes
within a fixed distance from the sample in the latent space. The method was
able to display the four prototypical brains most similar to a test image.

####  5.7.4 Attention deficit hyperactivity disorder classification

Another choice for intermediate interpretable features in a hybrid model is to
learn fbn that are important for the dl task. Qiang et al. trained such a ml +
ml hybrid model to classify adhd from rsfmri data [104]. First, a vae was
trained on fmri data and the latent representations learned by the vae were
used to learn fbn weights using Lasso regression (penalised regression with L1
penalty). For the second component, fc matrices were constructed from the fbn
weights and used to train an ml classifier to predict adhd. The fbn learned by
the vae were shown to be similar to those derived from another state-of-the-
art method.

####  5.7.5 Cognitive/clinical score regression

In a similar manner to Qiang et al., D’Souza et al. learned fbn as
interpretable intermediate features for cognitive and clinical score
prediction [121]. The authors coupled representation learning with a nn that
predicted cognitive or clinical scores. Intermediate fbn were learned from
rsfmri functional connectivity matrices using srddl. Simultaneously, an lstm
was trained from the subject-specific fbn coefficients. The ml + nn hybrid
model was trained to predict the Cognitive Fluid Intelligence Score for
healthy subjects from the hcp dataset [137], as well as clinical scores
(Autism Diagnostic Observation Schedule, Social Responsiveness Scale and
Praxis) for asd subjects from the kki dataset [15]. Fifteen learned fbn for
both hcp and kki data were presented as the model explanations. For example,
several of the networks learned to predict the Cognitive Fluid Intelligence
Score were involved in the Medial Frontal Network and the Frontal Parietal
Network, which play a role in decision-making, attention, and working memory.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
AD classification
Abuhmed et al. [1] | ADNI | sMRI (3D) | 1,371 | Cognitive scores
|  | PET (3D) |  |
Lee et al. [69] | ADNI | sMRI (3D) | 801 | Disease
|  |  |  | probability map
Liu et al. [77] | ADNI | sMRI (3D) | 3,021 | Int. features only
Nguyen et al. [94] |  Various† | sMRI (3D) | 2,036 | Disease grade map
Qiu et al. [105] |  Various†† | sMRI (3D) | 1,446 | Disease
|  |  |  | probability map
Mohammadjafari et al. [88] |  +ADNI | sMRI (2D) | 408 | Prototypes
| +OASIS |  |  |
Mulyadi et al. [92] |  +ADNI | sMRI (3D) | 2,285 | Prototypes
| +In-house |  |  |
Wolf et al. [145] | ADNI | PET (3D) | 1,245 | Prototypes
ASD classification
Kang et al. [55] | ABIDE | fMRI (3D) | 985 | Prototypes
Brain age regression
Hesse et al. [43] | IXI | sMRI (2D) | 562 | Prototypes
Cognitive/ clinical score regression
D’Souza et al. [121] |  +HCP | fMRI (3D) | 150 | FBNs
| +KKI | DTI (3D) |  |
ADHD classification
Qiang et al. [104] | ADHD-200 | fMRI (3D) | 541 | FBNs
Table 9: Articles using interpretable hybrid models or interpretable
intermediate features. Abbreviations: ABIDE, Autism Brain Imaging Data
Exchange; AD, Alzheimer Disease; ADHD, Attention Deficit Hyperactivity
Disorder; ADNI, Alzheimer’s Disease Neuroimaging Initiative; AIBL, Australian
Imaging Biomarker and Lifestyle Flagship Study of Ageing; ASD, Autism Spectrum
Disorder; DTI, Diffusion Tensor Imaging; FBN, Functional Brain Network, FHS,
Framingham Heart Study; fMRI, functional Magnetic Resonance Imaging; HCP,
Human Connectome Project; IXI, Information eXtraction from Images; KKI,
Kennedy Krieger Institute; MIRIAD, Minimal Interval Resonance Imaging in
Alzheimer’s Disease; NACC, National Alzheimer’s Coordinating Center; NIFD,
Frontotemporal lobar Degeneration Neuroimaging Initiative; OASIS, Open Access
Series of Imaging Studies; PET, Positron Emission Tomography; sMRI, structural
Magnetic Resonance Imaging. $\dagger$=ADNI+AIBL+OASIS+MIRIAD+NIFD,
$\dagger\dagger$=ADNI+AIBL+FHS+NACC.

###  5.8 Applications of interpretable generative models (Table 10)

####  5.8.1 Neurodegenerative disease maps

Several studies trained a gan [37] on the adni structural mri dataset to
predict disease effect maps for ad, considering either mci or cn as the
control class [12, 66, 80, 10, 9]. Baumgartner et al. developed a visual
attribution method based on a conditional gan (VA-gan) [12]. In this work, an
additive map $M(x)$ was learned as a function of an input image $x$ from the
ad class, such that the modified image $x+M(x)$ appears cognitively normal. In
contrast to learning an additive map, Lanfredi et al. trained a gan to
generate a deformation field, known as deformation field interpretation (DeFI-
gan), which was shown to produce sparser disease effect maps than VA-gan [66].
The deformation field transforms an image from the ad class to the mci class
by modelling brain atrophy. As such, deformation-based approaches are only
appropriate for modelling diseases where brain atrophy is the predominant
imaging marker. The same deformation field approach was employed by Liu et
al., but using a cyclegan that generated modified ad and cn images [80]. The
Jacobian of the deformation field was visualised as the disease effect map.

The aforementioned methods assume that the category labels of the test data
(either real or estimated by a separate classifier) are known during testing,
meaning that the models can generate explanations, but cannot perform the
classification. Bass et al.[10] developed a model that both classified disease
and generated disease effect maps. By incorporating a classification network,
this model obviates the need for previously classified data during testing. A
vae-gan was trained to disentangle class-relevant features from background
features, and therefore to separate the effects of healthy aging from disease.
The mean and variance of predicted disease effect maps were sampled from the
latent space during testing, as opposed to from a single additive map for each
subject. The method was applied to brain structural mri data from adni as in
[12] and disease effect maps were shown to improve when compared to VA-gan and
gradient-based methods. In addition to classification, the method was extended
for regression of mmse from structural mri adni data; regression of age from
Biobank brain structural MRI scans; and regression of birth age from dhcp data
[9]. All of these studies produced ad disease effect maps that successfully
modelled atrophy of the ventricles, hippocampus, and cortical grey matter
known to occur in ad.

####  5.8.2 Brain tumour and stroke segmentation

More recently, state-of-the-art diffusion models [45] have been trained to
predict disease effect maps (anomaly maps) for neuroimaging datasets [14, 116,
146]. Two studies trained a ddpm on smri images from the brats dataset to
convert a cancerous MRI to appear healthy, and a third study trained a
generative model to transform MRI brain images of stroke patients to appear
healthy. In all studies, the generated healthy image was subtracted from the
original to produce the anomaly map. Wolleb et al. trained an unconditional
ddpm and a classifier, and then used classifier guidance to transform an MRI
from cancerous to healthy [146]. In contrast, Sanchez et al. trained a
conditional ddpm and employed classifier-free guidance to alter the cancerous
images [116]. Bercea et al. implemented a two-stage approach, where stroke-
effected regions were removed from the image in stage one, and then stage two
comprised an in-painting generative model to fill in these erased regions as
healthy [14]. The anomaly maps in all three studies were shown to identify
pathological brain regions successfully.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
AD classification
Baumgartner et al. [12] | ADNI | sMRI (3D) | 1,288 | Generative additive maps
Bass et al. [10] | ADNI | sMRI (3D) | 1,053 | Generative additive maps
Bass et al. [9] | ADNI | sMRI (3D) | 1,053 | Generative additive maps
Lanfredi et al. [66] | ADNI | sMRI (3D) | 825 | Generative
|  |  |  | deformation fields
Liu et al. [80] | ADNI | sMRI (3D) | 1,344 | Generative
|  |  |  | deformation fields
Brain age regression
Bass et al. [9] |  Various† | sMRI (3D) | 12,434 | Generative additive maps
Brain tumour and stroke segmentation
Bercea et al. [14] |  Various‡ | sMRI (2D) | 1,412 | Generative additive maps
Sanchez et al. [116] | BraTS | sMRI (2D) | 1,251 | Generative additive maps
Wolleb et al. [146] | BraTS | sMRI (2D) | N/A | Generative additive maps
Table 10: Articles using generative models. Abbreviations: AD, Alzheimer’s
Disease; ADNI, Alzheimer’s Disease Neuroimaging Initiative; BraTS, Brain
Tumour Segmentation Challenge; dHCP, Developing Human Connectome Project; IXI,
Information eXtraction from Image; N/A, Not Available; sMRI, structural
Magnetic Resonance Imaging. $\dagger$=UK Biobank+dHCP, $\ddagger$=ATLAS
v2.0+IXI+FastMRI

###  5.9 Applications of deep structural causal models (Table 11 )

Reinhold et al. extended the dscm in Eqn. 5 to model causal effects for
structural mri images from a ms cohort by adding duration of ms symptoms,
expanded disability severity score, lesion volume and image slice number
[108]. Counterfactual difference maps were explored, such as the
counterfactual $do(l=0\text{ mL})$ for a brain mri of an ms patient, where the
model successfully removed the ms lesions from the counterfactual image.

Furthermore, Rasal et al. modified a dscm to synthesise 3D surface meshes of
the brain stem by introducing graph convolutional layers into the VAE [107].
The authors performed interventions on the population-mean brain stem, as well
as generating subject-specific counterfactual surface meshes for variables
such as age and sex. Realistic counterfactual meshes were generated for
scenarios outside the true data distribution, for example,
$do(\text{age}=80\text{year-old})$ when the maximum participant age was 70
years old.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
Image generation
Pawlowski et al. [100] | UK Biobank | sMRI (2D) | 13,750 | DSCM
Reinhold et al. [108] | In-house | sMRI (2D) | 77 | DSCM
Rasal et al. [107] | UK Biobank | Surface meshes (3D) | 14,502 | DSCM
Table 11: Articles using deep structural causal models. Abbreviations: DSCM,
Deep Structural Causal Model; sMRI, structural Magnetic Resonance Imaging.

###  5.10 Applications of attention mechanisms (Table 12)

####  5.10.1 Image segmentation

Gu et al. [38] introduced channel, spatial, and non-local attention blocks in
a modified U-Net to improve the performance of medical image segmentation
tasks. More specifically, they used spatial attention blocks throughout the
decoder layers of the U-Net by combining both higher (from the decoder) and
lower resolution (from the encoder) feature maps, similar to that proposed
previously [98]. Channel attention blocks were also introduced after each
decoding layer by global average pooling and global max pooling [147]. The
latter was also introduced as “scale attention”, which assigns a weight for
each of the decoder outputs to enable differential attention to be assigned to
a given input. The final non-local block was introduced at the lowest
resolution level (the bottleneck of the U-Net) due to its complexity. They
showed the spatial attention maps from the trained network were able to
highlight the object to be segmented, suggesting that the use of attention
enhanced the ability of the network to focus on target areas to facilitate
performance.

####  5.10.2 Disease classification

A 3D spatial attention network was used to classify ad using two large
structural mri datasets (adni and an in-house database) [53]. Following grey
matter segmentation, volumes were inputted into a 3D-cnn, which contained a
spatial attention block after the first three convolutional layers to
highlight important regions in the feature maps. However, the spatial
attention module contained a ReLU rather than a sigmoid activation function.
Thus, probability values for each spatial location were not produced, but
nevertheless, the method was able to identify those brain regions correlated
with atrophy, characteristic of ad.

Attention has also been introduced into a hybrid dl framework to classify scz
and asd using an in-house and the ABIDE rsfmri dataset, respectively [155].
Features were first extracted from the imaging data using principal components
analysis, and 50 independent components (IC) were retained per subject, each
of which was a times series. An attention-guided crnn was then used to process
the ic time series data, and a dnn for processing fnc matrices. The crnn
attention block aimed to highlight which rsfmri-derived ic were more
significant for prediction. The attention module was comparable to that
proposed by Woo et al. [147], which uses both max and average pooling layers,
but Zhao et al. [155] applied these along the time axis. The outputs of these
two separate networks were concatenated and passed through a logistic
regressor to obtain the final classification result.

Sarraf et al. [117] developed an optimized vision transformer, OViTAD, for
classifying hc, mci, and ad brains using rsfmri and structural mri data. The
authors also generated attention maps for ad vs. hc vs. mci classification for
the different self-attention heads, as well as global-level attention maps
extracted from the last feature vector.

####  5.10.3 Brain age regression

Finally, Dahan et al. [25] introduced the Surface Vision Transformer, which
adapted the image transformer model to surface domains. More specifically,
surface meshes were transformed into triangular patches and flattened into
feature vectors then inputted into the transformer model [136]. The main task
of their proposed study was to perform phenotype regression tasks using
cortical surface metrics from the dhcp. The authors also produced average
attention maps for either regression of postmenstrual age at scan and
gestational age at birth.

Reference | Data | Modality | #Subjects | Method
---|---|---|---|---
Image segmentation
Gu et al. [38] | In-house | sMRI (2D) | 36 | Spatial, channel and
|  |  |  | non-local attention
Disease Classification
Jin et al. [53] |  \+ ADNI | sMRI (3D) | 1,832 | Spatial attention
| \+ In-house |  |  |
Zhao et al. [155] |  \+ ABIDE | fMRI (ts) | 2,622 | Time-axis attention
| \+ In-house |  |  |
Sarraf et al. [117] | ADNI |  \+ fMRI (ts) | 1,744 | Self-attention
|  | \+ sMRI (2D) |  |
Brain age regression
Dahan et al. [25] | dHCP | Surface meshes (3D) | 588 | Self-attention
Table 12: Articles using attention-based methods. Abbreviations: ABIDE, Autism
Brain Imaging Data Exchange; ADNI, Alzheimer’s Disease Neuroimaging
Initiative; dHCP, developing Human Connectome Project; fMRI, functional
Magnetic Resonance Imaging; sMRI, structural Magnetic Resonance Imaging; ts,
time series.

##  6 Evaluation of idl explanations

Of utmost importance, idl explanations need to be evaluated for biological
validity and robustness. Biological validity refers to whether explanations
capture the true, underlying biological or pathological processes, and
robustness assesses the stability of an explanation under varying conditions.
Other properties of idl explanations that were evaluated in the literature are
continuity, selectivity, and downstream task performance. These properties
will be discussed below in the context of the 75 studies included in this
review.

###  6.1 Biological validity

A key challenge for idl in neuroimaging is that only appropriately trained
medical specialists, e.g., radiologists, can validate explanations.
Explanations for natural images can usually be readily validated by a general
audience; for example, the model predicts “castle” and the explanation
highlights the castle turrets. In contrast, years of specialised medical
training are required to identify imaging biomarkers, such as regional brain
atrophy in neurodegenerative diseases. Studies may be conducted where
clinicians evaluate idl explanations. However, due to limited clinician
availability, quantitative and automated validation metrics are more
desirable.

Most of the studies we reviewed did not validate (26 out of 75) or only
qualitatively validated (31 out of 75) the idl explanations—for example, many
studies compared salient brain regions identified in the explanations with
those previously reported. Several fMRI studies leveraged Neurosynth, a meta-
analysis platform that can return functional keywords correlated to idl
explanations, and compared these keywords against the literature.

The remaining 18 studies quantitatively compared idl explanations to ground
truth explanations, the latter obtained through various sources (Table 13). A
noteworthy example is where longitudinal imaging data are available, such as
in the adni database. For subjects that progressed from cn/mci to ad, a ground
truth disease effect map may be computed by subtracting the registered ad
image from the cn/mci image [12, 66, 10, 9]. The explanation maps were then
quantitatively compared to disease effect maps using ncc. Overall,
explanations from interpretable generative models achieved substantially
higher ncc ([12, 10, 9]) than explanations from popular post-hoc methods, such
as CAM, Guided Backpropagation, and Integrated Gradients.

Reference | Ground truth data source | Metric
---|---|---
Perturbation-based
Liu et al. [79] | cGAN-based statistics | # brain regions
Yang et al. [150] | 8 hold-out subjects | precision-recall curve
Gradient-based
Levakov et al. [71] | VBM meta-analysis | Mean VBM for top 1% regions
Ismail et al. [51] | Off-task data | % relevant features on-task
Backpropagation-based
Thomas et al. [135] | NeuroSynth meta-analysis | Mean F1 score
Weights analysis
Natekar et al. [93] | Ground truth segmentation | IoU
Disentangled latent space
Mouches et al. [91] | Ground truth segmentation | Lateral ventricle volume
Ouyang et al. [99] | ADAS-Cog scores | Correlation
Interpretable hybrid models
Qiu et al. [105] |  Post-mortem tissue | Correlation
Generative models
Baumgartner et al. [12] | ADNI disease effect map | NCC
Lanfredi et al. [66] | ADNI disease effect map | NCC
Bass et al. [10] | ADNI disease effect map | NCC
Bass et al. [9] | ADNI disease effect map | NCC
Sanchez et al. [116] | Ground truth segmentation | Dice
Wolleb et al.[146] | Ground truth segmentation | Dice
Bercea et al. [14] | Ground truth segmentation | Dice
Deep structural causal models
Reinhold et al. [108] | Image segmentation | MS lesion volume
Attention
Jin et al. [53] | AD MMSE scores | Correlation
Table 13: Quantitative metrics to evaluate biological validity of idl
explanations. Abbreviations: AD, Alzheimer’s Disease; ADAS-Cog, Alzheimer’s
Disease Assessment Scale – Cognitive Subscale; ADNI, Alzheimer’s Disease
Neuroimaging Initiative; cGAN, conditional generative adversarial network;
IoU, intersection over union; MMSE, Mini-mental state examination; MS,
Multiple Sclerosis; NCC, normalised cross-correlation; VBM, Voxel-based
morphometry.

###  6.2 Robustness

Robustness was not evaluated in the majority of studies (62 out of 75). In the
remaining studies, the robustness of population-level explanations was
considered with respect to different training data ([53], [35], [134]), data
pre-processing methods ([85, 72]), and model and idl settings [30, 121, 134,
58, 31, 73, 71]. Three studies compared population-level explanations with the
same dl task and model architecture but where the model was trained on
different sources of data, and all concluded explanations were stable across
datasets. For example, Jin et al. compared attention maps from a ResNet
trained on structural mri adni data versus a similar in-house dataset and
found the maps were significantly correlated (r=0.59) [53]. A few studies
considered explanations trained on the same data source but with different
pre-processing methodologies, investigating different atlases and atlas
granularities during registration [85, 72]. Furthermore, robustness of
explanations across different model and idl settings was evaluated, including
cross-validation folds [30, 121, 134, 58]; parameter initialisation [31, 134];
hyperparameter values [73, 134]; and models within an ensemble [71].

Data preprocessing methods, hyperparameters, and model parameters all
influence the explanations produced. Concerning data preprocessing, skull
stripping often alters downstream explanations [57, 29]. In another example,
Mellema et al. showed the level of atlas granularity during registration
altered the important features identified for asd classification [85].
Concerning hyperparameters, the selection of regularisation weights for
Meaningful Perturbations changed the explanation masks for ad classification
[134]. Evidence also suggests that different runs of identically trained,
randomly initialised models are associated with markedly different
explanations [31, 134]. It is important to be aware that bias may be present
in idl explanations from sources such as data preprocessing and hyperparameter
selection and to assess explanations for such bias.

The robustness of explanations under different conditions may be
quantitatively assessed using various similarity measures (see Table 14). Some
studies directly compared explanations using overlap measures such as the Dice
coefficient or Hausdorff distance. Other studies initially converted an
explanation into a vector of mean values for $n$ atlas-derived brain regions
and then compared vectors using correlation [53, 30], cosine similarity [134]
or percentage agreement between top regions [31].

Reference | Robustness across… | Metric
---|---|---
Perturbation-based
Thibeau et al. [134] | Datasets | Cosine similarity
Thibeau et al. [134] | Models (hyperparameters) | Cosine similarity
Thibeau et al. [134] | Models (cv folds + initialisation) | Cosine similarity
Li et al. [73] | Models (hyperparameters) | # important ROIs
Eitel et al. [31]† | Models (initialisation) | L2-norm
|  | \+ relevant region coherence
Gradient-based
Levakov et al. [71] | Models (ensemble) | Dice + Hausdorff distance
Class activation maps
Kim et al. [58] | Models (cv folds) | Relevant region coherence
Weights analysis
Dvornek et al. [30] | Models (cv folds) | Correlation + Dice
Interpretable hybrid models
Shimona et al. [121] | Models (cv folds) | Mean inner-product
Attention
Jin et al. [53] | Datasets | Correlation
Table 14: Quantitative metrics to evaluate robustness of idl explanations.
Abbreviations: cv, cross-validation; ROI, Region of Interest. $\dagger$: Eitel
et al. evaluated robustness across multiple method categories, not only
perturbation-based.

###  6.3 Other interpretable method properties

####  6.3.1 Continuity

Similar images should have similar explanations, as originally proposed by
Montavon et al. [90]. Nigri et al. measured continuity by slightly perturbing
50 input images and then calculating the mean $L_{2}$-norm between
explanations of the original and perturbed image [95]. The authors compared
the continuity of the Swap Test and Occlusion and found the Swap Test to be
the superior method.

####  6.3.2 Selectivity

Regions with the highest relevance in the explanation should result in the
largest change in model prediction when removed from the input image [90]. For
example, sensitivity maps were computed by Nigri et al., highlighting those
image regions swapped (Swap Test) or occluded (Occlusion) that resulted in a
large change in model prediction [95]. Reverse sensitivity maps were then
generated by removing the complement of each image patch and recording the
change in model prediction. Subsequently, Pearson correlation analysis was
carried out to assess the relationship between the standard and reversed
sensitivity maps, with strong negative correlations expected when the property
of selectivity is satisfied. Each image from the ms class underwent lesion in-
painting, such that ms lesions appeared to be healthy tissue in the mri image.
Explanations were generated for the original and in-painted images, and the
difference between their mean values was computed, with larger differences
across all images suggesting a more selective idl method.

####  6.3.3 Downstream task performance (disentangled latent space methods
only)

Performance relates to whether the latent space distinguishes classes
sufficiently for a given downstream task. In Ouyang et al., dl models were
trained on the disentangled latent embeddings for two classification tasks to
understand if the latent space learned a meaningful structure. The evaluation
metric was the test set classification accuracy.

Reference | Interpretable method | Metric
---|---|---
Continuity
Nigri et al. [95] | Perturbation-based |  $L_{2}$ norm
Selectivity
Nigri et al. [95] | Perturbation-based | Correlation
Downstream task performance (disentangled latent space methods only)
Ouyang et al. [99] | Disentangled latent space | Classification accuracy
Table 15: Other properties and quantitative metrics for iDL explanations.

##  7 Discussion and conclusion

In this review, we identified 75 neuroimaging studies that utilised idl
methods, and we classified the methods into five post-hoc and five intrinsic
categories. To the best of our knowledge, this is the first systematic review
of idl in neuroimaging with a notably more extensive review of intrinsic
methods than found in the literature [133]. In addition, we found five
properties of idl explanations that were investigated and are important when
considering the suitability of an idl method for adoption.

The most common idl methods utilised were class activation maps, perturbation-
based and gradient-based methods. Post-hoc methods are popular because they
are well-established in computer vision tasks, easy to implement, and readily
available in dl packages. However, historically, post-hoc methods were
designed for and validated on natural images and may be inappropriate for
neuroimaging tasks. For example, saliency methods were shown to only focus on
a few discriminative features of a class [153, 10], rather than identifying
all imaging features, which may be important for diagnosis and treatment.
Their reliability is also questionable as some post-hoc methods, in particular
Guided backpropagation and Guided Grad-CAM [2], still produce convincing
explanations despite randomised model weights or data labels. In contrast,
intrinsic methods are generally more appropriate for neuroimaging because they
are designed specifically for the application e.g., constructing a causal
graph specific to ms [108]. Additionally, generative models produced
explanations with substantially higher correlation to ground-truth disease
markers compared to explanations from several post-hoc methods [12, 10, 9].
Nevertheless, intrinsic interpretable deep learning is still an emerging
field, and such methods are currently more time-consuming to implement than
post-hoc methods.

We will now provide some recommendations for researchers when using idl with
neuroimaging datasets. First, we suggest utilising multiple idl methods,
including several across different post-hoc method categories (such as
Occlusion, LRP, and GradCAM) and one intrinsic method that is best suited for
the project application, end-user requirements, objectives, etc. It is
important to carefully select one intrinsic method during the design phase as
it can be time-consuming to implement. For pre-existing models, incorporating
an appropriate attention mechanism and retraining the model may be feasible.
Then, compare explanations from different methods and prioritise features
highlighted across all methods.
Second, recall that various confounding factors, such as data preprocessing,
random initialisation, and cross-validation, can affect explanations.
Therefore, we advise averaging explanations across cross-validation folds and
multiple runs to improve robustness. Also, consider visualising explanations
for a reasonable selection of model preprocessing and hyperparameter settings.
If using multiple neuroimaging datasets, we recommend adopting a standardised
pre-processing pipeline to reduce the risk of biased explanations.

Third, validating explanations across an entire test dataset is crucial rather
than limiting assessments to a select few samples. This comprehensive
validation helps ensure the generalisability of the explanations. Consider
acquiring ground truth to validate explanations quantitatively, such as
computing disease affect maps from longitudinal imaging datasets. If not
possible, then impartially compare explanations to existing physiopathological
literature. In summary, do not unquestioningly trust the explanations produced
by an idl method.

When applying idl methods for neuroimaging, an important concern is the
complexity of the biological mechanisms underlying the data and the
interactions between multiple imaging features. Many interpretability methods
identified in this review do not consider the causal mechanisms that
contribute to the data nor the impact of confounding factors in the
explanations. We have, however, discussed state-of-the-art causal models that
attempt to address causality in interpretability, and we foresee such models
will play an important role in the future of idl [100, 108]. We also conclude
a suite of standardised, quantitative evaluation metrics to compare
performance across idl methods needs to be established to promote the
trustworthiness of idl methods.

##  8 Acknowledgments

We would like to thank the EPSRC Centre for Doctoral Training in Smart Medical
Imaging [EP/S022104/1] for funding studentships for LM and MDS. Thanks also to
the National Brain Mapping Laboratory (NBML) and the IUMS Medical Physics
Department for supporting and financing the Masters program for FH.

The images used in some of the figures in this work were funded by the
Alzheimer’s Disease Neuroimaging Initiative (ADNI) [National Institutes of
Health Grant U01 AG024904] and DOD ADNI [Department of Defense award number
W81XWH-12-2-0012].

Declarations of interest: The authors have no conflicts of interest to
disclose.

## References

  * [1] Abuhmed, T., El-Sappagh, S., Alonso, J.M.: Robust hybrid deep learning models for alzheimer’s progression detection. Knowledge-Based Systems 213, 106688 (2021)
  * [2] Adebayo, J., Gilmer, J., Muelly, M., Goodfellow, I., Hardt, M., Kim, B.: Sanity checks for saliency maps. Advances in neural information processing systems 31 (2018)
  * [3] Afshar, P., Mohammadi, A., Plataniotis, K.N.: Brain tumor type classification via capsule networks. In: 2018 25th IEEE International Conference on Image Processing (ICIP). pp. 3129–3133. IEEE (2018)
  * [4] Alangari, N., El Bachir Menai, M., Mathkour, H., Almosallam, I.: Exploring evaluation methods for interpretable machine learning: a survey. Information 14(8), 469 (2023)
  * [5] Azcona, E.A., Besson, P., Wu, Y., Punjabi, A., Martersteck, A., Dravid, A., Parrish, T.B., Bandt, S.K., Katsaggelos, A.K.: Interpretation of brain morphology in association to alzheimer’s disease dementia classification using graph convolutional networks on triangulated meshes. In: International Workshop on Shape in Medical Imaging. pp. 95–107. Springer (2020)
  * [6] Bach, S., Binder, A., Montavon, G., Klauschen, F., Müller, K.R., Samek, W.: On pixel-wise explanations for non-linear classifier decisions by layer-wise relevance propagation. PloS one 10(7), e0130140 (2015)
  * [7] Bakas, S., Akbari, H., Sotiras, A., Bilello, M., Rozycki, M., Kirby, J.S., Freymann, J.B., Farahani, K., Davatzikos, C.: Advancing the cancer genome atlas glioma mri collections with expert segmentation labels and radiomic features. Scientific data 4(1), 1–13 (2017)
  * [8] Balduzzi, D., Frean, M., Leary, L., Lewis, J., Ma, K.W.D., McWilliams, B.: The shattered gradients problem: If resnets are the answer, then what is the question? In: International Conference on Machine Learning. pp. 342–350. PMLR (2017)
  * [9] Bass, C., Da Silva, M., Sudre, C., Williams, L.Z., Sousa, H.S., Tudosiu, P.D., Alfaro-Almagro, F., Fitzgibbon, S.P., Glasser, M.F., Smith, S.M., et al.: Icam-reg: Interpretable classification and regression with feature attribution for mapping neurological phenotypes in individual scans. IEEE Transactions on Medical Imaging 42(4), 959–970 (2022)
  * [10] Bass, C., da Silva, M., Sudre, C., Tudosiu, P.D., Smith, S., Robinson, E.: Icam: interpretable classification via disentangled representations and feature attribution mapping. Advances in Neural Information Processing Systems 33, 7697–7709 (2020)
  * [11] Bau, D., Zhou, B., Khosla, A., Oliva, A., Torralba, A.: Network dissection: Quantifying interpretability of deep visual representations. In: Proceedings of the IEEE conference on computer vision and pattern recognition. pp. 6541–6549 (2017)
  * [12] Baumgartner, C.F., Koch, L.M., Tezcan, K.C., Ang, J.X., Konukoglu, E.: Visual feature attribution using wasserstein gans. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition. pp. 8309–8319 (2018)
  * [13] Bengio, Y., Simard, P., Frasconi, P.: Learning long-term dependencies with gradient descent is difficult. IEEE transactions on neural networks 5(2), 157–166 (1994)
  * [14] Bercea, C.I., Wiestler, B., Rueckert, D., Schnabel, J.A.: Reversing the abnormal: Pseudo-healthy generative networks for anomaly detection. In: International Conference on Medical Image Computing and Computer-Assisted Intervention. pp. 293–303. Springer (2023)
  * [15] Bigler, E.D.: Neuropsychology and clinical neuroscience of persistent post-concussive syndrome. Journal of the International Neuropsychological Society 14(1), 1–22 (2008)
  * [16] Bintsi, K.M., Baltatzis, V., Hammers, A., Rueckert, D.: Voxel-level importance maps for interpretable brain age estimation. In: Interpretability of Machine Intelligence in Medical Image Computing, and Topological Data Analysis and Its Applications for Medical Data, pp. 65–74. Springer (2021)
  * [17] Böhle, M., Eitel, F., Weygandt, M., Ritter, K.: Layer-wise relevance propagation for explaining deep neural network decisions in mri-based alzheimer’s disease classification. Frontiers in aging neuroscience 11, 194 (2019)
  * [18] Brauwers, G., Frasincar, F.: A general survey on attention mechanisms in deep learning. IEEE Transactions on Knowledge and Data Engineering (2021)
  * [19] de Brebisson, A., Montana, G.: Deep neural networks for anatomical brain segmentation. In: Proceedings of the IEEE conference on computer vision and pattern recognition workshops. pp. 20–28 (2015)
  * [20] Calhoun, V.D., Adali, T., Pearlson, G.D., Pekar, J.J.: A method for making group inferences from functional mri data using independent component analysis. Human brain mapping 14(3), 140–151 (2001)
  * [21] Carmichael, Z., Lohit, S., Cherian, A., Jones, M.J., Scheirer, W.J.: Pixel-grounded prototypical part networks. In: Proceedings of the IEEE/CVF Winter Conference on Applications of Computer Vision. pp. 4768–4779 (2024)
  * [22] Castro, D.C., Walker, I., Glocker, B.: Causality matters in medical imaging. Nature Communications 11(1), 1–10 (2020)
  * [23] Chen, C., Li, O., Tao, D., Barnett, A., Rudin, C., Su, J.K.: This looks like that: deep learning for interpretable image recognition. Advances in neural information processing systems 32 (2019)