Add files via upload

PytorchLRP/README.md

@@ -0,0 +1,74 @@
+# Layer-wise relevance propagation for explaining deep neural network decisions in MRI-based Alzheimer’s disease classification
+
+**Moritz Böhle, Fabian Eitel, Martin Weygandt, and Kerstin Ritter**
+
+**Preprint:** https://arxiv.org/abs/1903.07317
+
+**Abstract:** Deep neural networks have led to state-of-the-art results in many medical imaging tasks including Alzheimer's disease (AD) detection based on structural magnetic resonance imaging (MRI) data. However, the network decisions are often perceived as being highly non-transparent making it difficult to apply these algorithms in clinical routine. In this study, we propose using layer-wise relevance propagation (LRP) to visualize convolutional neural network decisions for AD based on MRI data. Similarly to other visualization methods, LRP produces a heatmap in the input space indicating the importance / relevance of each voxel contributing to the final classification outcome. In contrast to susceptibility maps produced by guided backpropagation ("Which change in voxels would change the outcome most?"), the LRP method is able to directly highlight positive contributions to the network classification in the input space. In particular, we show that (1) the LRP method is very specific for individuals ("Why does this person have AD?") with high inter-patient variability, (2) there is very little relevance for AD in healthy controls and (3) areas that exhibit a lot of relevance correlate well with what is known from literature. To quantify the latter, we compute size-corrected metrics of the summed relevance per brain area, e.g. relevance density or relevance gain. Although these metrics produce very individual 'fingerprints' of relevance patterns for AD patients, a lot of importance is put on areas in the temporal lobe including the hippocampus. After discussing several limitations such as sensitivity towards the underlying model and computation parameters, we conclude that LRP might have a high potential to assist clinicians in explaining neural network decisions for diagnosing AD (and potentially other diseases) based on structural MRI data.
+
+## Requirements
+
+In order to run the code, standard pytorch packages and Python 3 are needed.
+Moreover, add a settings.py file to the repo, containing the data paths and so forth as follows:
+
+Please use the example settings.py with more information.
+
+```python
+settings = {
+    "model_path": INSERT,
+    "data_path": INSERT,
+    "ADNI_DIR": INSERT,
+    "train_h5": INSERT,
+    "val_h5": INSERT,
+    "holdout_h5": INSERT,
+    "binary_brain_mask": "binary_brain_mask.nii.gz",
+    "nmm_mask_path": "~/spm12/tpm/labels_Neuromorphometrics.nii",
+    "nmm_mask_path_scaled": "nmm_mask_rescaled.nii"
+}
+```
+
+With the "Evaluate GB and LRP" notebook, the heatmap results and the summed scores per area can be calculated.
+The notebooks "Plotting result graphs" and "Plotting brain maps" can be used to calculate and plot the results according to the defined metrics and show the heatmaps of individual patient's brains and average heatmaps according to LRP and GB.
+
+## Quickstart
+
+You can use the visualization methods in this repository on your own model (PyTorch; for a Keras implementation, see heatmapping.org) like this:
+
+```python
+model = Net()
+model.load_state_dict(torch.load("./mymodel"))
+# Convert to innvestigate model
+inn_model = InnvestigateModel(model, lrp_exponent=2,
+                                method="e-rule",
+                                beta=.5)
+model_prediction, heatmap = inn_model.innvestigate(in_tensor=data)
+```
+
+`heatmap` contains the relevance heatmap. The methods should work for 2D and 3D images alike, see the MNIST example notebook or the LRP and GB evaluation notebook for an example with MRI images.
+
+### Docker
+
+To run [MNIST example.ipynb](./MNIST%20example.ipynb) in a Docker container (using only CPU) follow the steps below:
+
+```sh
+cd docker/
+docker-compose up --detach
+```
+
+Visit [localhost:7700](http://localhost:7700) in your browser to open Jupyter.
+
+## Code Structure
+
+The repository consists of the general LRP wrapper ([innvestigator.py](innvestigator.py) and [inverter_util.py](inverter_util.py)), a simple example for applying the wrapper in the case of MNIST data, and the evaluation notebook for obtaining the heatmap results discussed in the article.
+
+## Heatmaps
+
+The methods for obtaining the heatmaps are shown in the notebook **Evaluate GB and LRP**
+
+## Data
+
+The MRI scans used for training are from the [Alzheimer Disease Neuroimaging Initiative (ADNI)](http://adni.loni.usc.edu/). The data is free but you need to apply for access on http://adni.loni.usc.edu/. Once you have an account, go [here](http://adni.loni.usc.edu/data-samples/access-data/) and log in. Settings.py gives information about the required data format.
+
+## Citation
+
+If you use our code, please cite us. The code is published under the BSD License 2.0.

PytorchLRP/docker/Dockerfile

@@ -0,0 +1,19 @@
+FROM jupyter/scipy-notebook:76402a27fd13
+
+WORKDIR /home/jovyan/work
+
+USER root
+
+# Style profile
+RUN echo "export CLICOLOR=1" >> /root/.bashrc \
+  && echo "export PS1='\u@\h:\[\e[0;34m\]\W\[\e[0m\]\$ '" >> /root/.bashrc \
+  && echo "export LS_OPTIONS='--color=auto'" >> /root/.bashrc \
+  && echo "export LSCOLORS=excxhxDxfxhxhxhxhxFxFx" >> /root/.bashrc \
+  && echo 'alias ls="ls $LS_OPTIONS"' >> /root/.bashrc \
+  && echo 'alias ll="ls $LS_OPTIONS -l"' >> /root/.bashrc \
+  && echo 'alias l="ls $LS_OPTIONS -lA"' >> /root/.bashrc
+
+USER jovyan
+
+COPY requirements.txt /requirements.txt
+RUN pip install -r /requirements.txt

PytorchLRP/docker/docker-compose.yml

@@ -0,0 +1,16 @@
+version: '3.8'
+services:
+  lrp:
+    build:
+      context: ..
+      dockerfile: docker/Dockerfile
+    tty: true
+    command: bash -c "jupyter lab --port=8888 --no-browser --ip=0.0.0.0 --allow-root --NotebookApp.token='' --NotebookApp.password=''"
+    volumes:
+      - ..:/home/jovyan/work
+    ports:
+      - '7700:8888'
+    restart: unless-stopped
+    container_name: lrp
+    environment:
+      - JUPYTER_ENABLE_LAB=yes

PytorchLRP/innvestigator.py

@@ -0,0 +1,235 @@
+import torch
+import numpy as np
+
+from inverter_util import RelevancePropagator
+from utils import pprint, Flatten
+
+
+class InnvestigateModel(torch.nn.Module):
+    """
+    ATTENTION:
+        Currently, innvestigating a network only works if all
+        layers that have to be inverted are specified explicitly
+        and registered as a module. If., for example,
+        only the functional max_poolnd is used, the inversion will not work.
+    """
+
+    def __init__(self, the_model, lrp_exponent=1, beta=.5, epsilon=1e-6,
+                 method="e-rule"):
+        """
+        Model wrapper for pytorch models to 'innvestigate' them
+        with layer-wise relevance propagation (LRP) as introduced by Bach et. al
+        (https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0130140).
+        Given a class level probability produced by the model under consideration,
+        the LRP algorithm attributes this probability to the nodes in each layer.
+        This allows for visualizing the relevance of input pixels on the resulting
+        class probability.
+
+        Args:
+            the_model: Pytorch model, e.g. a pytorch.nn.Sequential consisting of
+                        different layers. Not all layers are supported yet.
+            lrp_exponent: Exponent for rescaling the importance values per node
+                            in a layer when using the e-rule method.
+            beta: Beta value allows for placing more (large beta) emphasis on
+                    nodes that positively contribute to the activation of a given node
+                    in the subsequent layer. Low beta value allows for placing more emphasis
+                    on inhibitory neurons in a layer. Only relevant for method 'b-rule'.
+            epsilon: Stabilizing term to avoid numerical instabilities if the norm (denominator
+                    for distributing the relevance) is close to zero.
+            method: Different rules for the LRP algorithm, b-rule allows for placing
+                    more or less focus on positive / negative contributions, whereas
+                    the e-rule treats them equally. For more information,
+                    see the paper linked above.
+        """
+        super(InnvestigateModel, self).__init__()
+        self.model = the_model
+        self.device = torch.device("cpu", 0)
+        self.prediction = None
+        self.r_values_per_layer = None
+        self.only_max_score = None
+        # Initialize the 'Relevance Propagator' with the chosen rule.
+        # This will be used to back-propagate the relevance values
+        # through the layers in the innvestigate method.
+        self.inverter = RelevancePropagator(lrp_exponent=lrp_exponent,
+                                            beta=beta, method=method, epsilon=epsilon,
+                                            device=self.device)
+
+        # Parsing the individual model layers
+        self.register_hooks(self.model)
+        if method == "b-rule" and float(beta) in (-1., 0):
+            which = "positive" if beta == -1 else "negative"
+            which_opp = "negative" if beta == -1 else "positive"
+            print("WARNING: With the chosen beta value, "
+                  "only " + which + " contributions "
+                  "will be taken into account.\nHence, "
+                  "if in any layer only " + which_opp +
+                  " contributions exist, the "
+                  "overall relevance will not be conserved.\n")
+
+    def cuda(self, device=None):
+        self.device = torch.device("cuda", device)
+        self.inverter.device = self.device
+        return super(InnvestigateModel, self).cuda(device)
+
+    def cpu(self):
+        self.device = torch.device("cpu", 0)
+        self.inverter.device = self.device
+        return super(InnvestigateModel, self).cpu()
+
+    def register_hooks(self, parent_module):
+        """
+        Recursively unrolls a model and registers the required
+        hooks to save all the necessary values for LRP in the forward pass.
+
+        Args:
+            parent_module: Model to unroll and register hooks for.
+
+        Returns:
+            None
+
+        """
+        for mod in parent_module.children():
+            if list(mod.children()):
+                self.register_hooks(mod)
+                continue
+            mod.register_forward_hook(
+                self.inverter.get_layer_fwd_hook(mod))
+            if isinstance(mod, torch.nn.ReLU):
+                mod.register_backward_hook(
+                    self.relu_hook_function
+                )
+
+    @staticmethod
+    def relu_hook_function(module, grad_in, grad_out):
+        """
+        If there is a negative gradient, change it to zero.
+        """
+        return (torch.clamp(grad_in[0], min=0.0),)
+
+    def __call__(self, in_tensor):
+        """
+        The innvestigate wrapper returns the same prediction as the
+        original model, but wraps the model call method in the evaluate
+        method to save the last prediction.
+
+        Args:
+            in_tensor: Model input to pass through the pytorch model.
+
+        Returns:
+            Model output.
+        """
+        return self.evaluate(in_tensor)
+
+    def evaluate(self, in_tensor):
+        """
+        Evaluates the model on a new input. The registered forward hooks will
+        save all the data that is necessary to compute the relevance per neuron per layer.
+
+        Args:
+            in_tensor: New input for which to predict an output.
+
+        Returns:
+            Model prediction
+        """
+        # Reset module list. In case the structure changes dynamically,
+        # the module list is tracked for every forward pass.
+        self.inverter.reset_module_list()
+        self.prediction = self.model(in_tensor)
+        return self.prediction
+
+    def get_r_values_per_layer(self):
+        if self.r_values_per_layer is None:
+            pprint("No relevances have been calculated yet, returning None in"
+                   " get_r_values_per_layer.")
+        return self.r_values_per_layer
+
+    def innvestigate(self, in_tensor=None, rel_for_class=None):
+        """
+        Method for 'innvestigating' the model with the LRP rule chosen at
+        the initialization of the InnvestigateModel.
+        Args:
+            in_tensor: Input for which to evaluate the LRP algorithm.
+                        If input is None, the last evaluation is used.
+                        If no evaluation has been performed since initialization,
+                        an error is raised.
+            rel_for_class (int): Index of the class for which the relevance
+                        distribution is to be analyzed. If None, the 'winning' class
+                        is used for indexing.
+
+        Returns:
+            Model output and relevances of nodes in the input layer.
+            In order to get relevance distributions in other layers, use
+            the get_r_values_per_layer method.
+        """
+        if self.r_values_per_layer is not None:
+            for elt in self.r_values_per_layer:
+                del elt
+            self.r_values_per_layer = None
+
+        with torch.no_grad():
+            # Check if innvestigation can be performed.
+            if in_tensor is None and self.prediction is None:
+                raise RuntimeError("Model needs to be evaluated at least "
+                                   "once before an innvestigation can be "
+                                   "performed. Please evaluate model first "
+                                   "or call innvestigate with a new input to "
+                                   "evaluate.")
+
+            # Evaluate the model anew if a new input is supplied.
+            if in_tensor is not None:
+                self.evaluate(in_tensor)
+
+            # If no class index is specified, analyze for class
+            # with highest prediction.
+            if rel_for_class is None:
+                # Default behaviour is innvestigating the output
+                # on an arg-max-basis, if no class is specified.
+                org_shape = self.prediction.size()
+                # Make sure shape is just a 1D vector per batch example.
+                self.prediction = self.prediction.view(org_shape[0], -1).to(self.device)
+                max_v, _ = torch.max(self.prediction, dim=1, keepdim=True)
+                # max_v.to(self.device)
+                only_max_score = torch.zeros_like(self.prediction).to(self.device)
+                only_max_score[max_v == self.prediction] = self.prediction[max_v == self.prediction]
+                relevance_tensor = only_max_score.view(org_shape)
+                self.prediction.view(org_shape)
+
+            else:
+                org_shape = self.prediction.size()
+                self.prediction = self.prediction.view(org_shape[0], -1)
+                only_max_score = torch.zeros_like(self.prediction).to(self.device)
+                only_max_score[:, rel_for_class] += self.prediction[:, rel_for_class]
+                relevance_tensor = only_max_score.view(org_shape)
+                self.prediction.view(org_shape)
+
+            # We have to iterate through the model backwards.
+            # The module list is computed for every forward pass
+            # by the model inverter.
+            rev_model = self.inverter.module_list[::-1]
+            relevance = relevance_tensor.detach()
+            del relevance_tensor
+            # List to save relevance distributions per layer
+            r_values_per_layer = [relevance]
+            for layer in rev_model:
+                # Compute layer specific backwards-propagation of relevance values
+                relevance = self.inverter.compute_propagated_relevance(layer, relevance)
+                r_values_per_layer.append(relevance.cpu())
+
+            self.r_values_per_layer = r_values_per_layer
+
+            del relevance
+            if self.device.type == "cuda":
+                torch.cuda.empty_cache()
+            return self.prediction, r_values_per_layer[-1]
+
+    def forward(self, in_tensor):
+        return self.model.forward(in_tensor)
+
+    def extra_repr(self):
+        r"""Set the extra representation of the module
+
+        To print customized extra information, you should re-implement
+        this method in your own modules. Both single-line and multi-line
+        strings are acceptable.
+        """
+        return self.model.extra_repr()