RAMP: segmentation of the brain lesions

Authors: Alexandre Hutton $^{1}$, Maria Teleńczuk $^{2}$, Swetha Shanker $^{2}$, Guillaume Lemaitre $^{2}$, François Caud $^{2}$, Alexandre Gramfort $^{2}$, Sook-Lei Liew $^{1}$

$^{1}$: Neural Plasticity & Neurorehabilitation Lab (Univ. of Southern California), USA
$^{2}$: Paris-Saclay Center for Data Science (Dataia - Univ. Paris Saclay - Inria), France

Table of contents

  1. Introduction
    1. Motivation
    2. Setup
  2. Data Exploration
    1. Visualize Data
    2. Paired Loading
    3. Lesion Analysis
  3. Creating Estimators
  4. Scores
  5. Submission

1. Introduction

Clinical/research motivation

Stroke is the leading cause of adult disability worldwide, with up to two-thirds of affected individuals experiencing long-term disabilities. Large-scale neuroimaging studies have shown promise in identifying robust biomarkers (e.g. measures of brain structure) of long-term stroke recovery following rehabilitation. However, analyzing large rehabilitation-related datasets is problematic due to barriers in accurate stroke lesion segmentation. Manually-traced lesions are currently the gold standard for lesion segmentation on T1-weighted MRIs, but require anatomical expertise and are labor-intensive. Further, manual segmentation is subjective with raters producing different results. While algorithms have been developed to automate this process, the resulting lesion masks often lack the accuracy needed to make them reliable information. Newer algorithms that employ machine-learning and deep learning techniques are promising avenues, but they require large, diverse datasets for training and testing and developing generalizable models. In this challenge, training can be performed on our public ATLAS 2.0 dataset, and testing is done with a private dataset comprised of multi-site data from the same sites as ATLAS 2.0.

For more information refer to: Anatomical Tracings of Lesions After Stroke

References

1. Maier, Oskar, et al. "ISLES 2015-A public evaluation benchmark for ischemic stroke lesion segmentation from multispectral MRI." Medical image analysis 35 (2017): 250-269.

2. Liew, Sook-Lei, et al. "A large, open source dataset of stroke anatomical brain images and manual lesion segmentations." Scientific data 5 (2018): 180011.

Objective of the challenge

In this challenge you will be given 3D medical images (T1 MRI scans in nii.gz format) of stroke patients and the files with the corresponding lesions (binary masks) traced manually by experts. Your algorithm will be scored on the overlap between your prediction and the segmentation mask.

Setup

Prerequisites

The following cell will install the required package dependencies, if necessary. You can examine the file, requirements.txt, included in the repo to view the list of dependencies. NOTE: Due to the structure of the challenge, libraries not included in requirements.txt will need to be added via a pull request to the GitHub repo.

To get this notebook running and test your models locally using ramp-test (from ramp-workflow), we recommend that you use the Anaconda or Miniconda Python distribution.

Data exploration

Data download and reformat

You will need four things to be able to get a local copy of the data:

  1. Downloaded archive.
  2. Encryption password.
  3. Decryption + tar extraction.
  4. Data reformatting.

1. Archive download

You can get a local (encrypted) copy of the data by going to the NITRC website or running the following cell. The relevant code is commented to prevent accidentally downloading the (large) data. Make sure to uncomment the last line.

Note that the data is approx. 15GB; ensure that you have enough space.

2. Encryption password

The data is encrypted; to obtain the password to decrypt the data, you'll need to agree to the terms of use found on this form. Please be careful not to close the tab before acquiring the decryption key.

3. OpenSSL + tar extraction

If you are running either Linux or MacOS, you can run the following cell. You will be prompted for the password obtained in step 2; you can copy+paste it.

For Windows, you'll need to install OpenSSL, open the terminal, then run the following from that window:
openssl aes-256-cbc -md sha256 -d -a -in ATLAS_R2.0_encrypted.tar.gz -out ATLAS_R2.0.tar.gz
Next, extract the .tar.gz archive. You should be left with a directory called ATLAS_2 once it is extracted.

4. Data reformatting

Unfortunately, the data provided via INDI is not compatible with PyBIDS, which is used throughout the challenge. You can run the following cell to correct the formatting; update the values as appropriate. WARNING: This will move the data files on your disk.

Viewing the data files

In your data folder you should now be able to find the test and train directories, each containing BIDS directories. Subject directories contain two files:

If you wish to view any of those files from outside of Python, you can use applications like ITK-snap. There, you can load the structural file as the main image and the lesion mask as a segmentation to overlap the two.

To load the images to Python, we will be using PyBIDS. See this link to read details about BIDS. If you want to use a libarary which is not listed in requirements.txt or extra_libraries.txt, please make a pull request to the the Stroke Lesion Segmentation Challenge repository on GitHub by adding the required library to extra_libraries.txt file).

Visualize data

The data provided with this challenge are 3D images. In this section, we'll show you how to view slices of the data. Since the data is in BIDS, we can access the data fairly easily using PyBIDS.

We see that there are 436 subjects in the training set with a total of 436 sessions.

In the above image, we see slices for each plane. The masked stroke is denoted in red, although the slices may not intersect the lesion mask for every subject (re-run the last two cells to display a new subject). Some preprocessing has already been done on the images; you can review the preprocessing steps here.

Due to the unreliable performance of non-linear registration with subjects with pathologies, the data was registered linearly. If you examine multiple subjects, you may notice that although the image dimensions are the same and the brains are roughly in the right area, there are differences between subjects (e.g., dimensions of the head, location of specific structures). This is typical of neuroimaging datasets.

Before we look closer at the lesions, let's look at some of the image properties:

The previous cell displays some of the metadata about the file. You can read here for more information, but the shape of the data is most relevant. Every image contains 197x233x189 (8.6M) voxels. Consequently, it is not usually possible to load the entire dataset into memory at once. We will see the effects of this in later sections.

Paired Loading

Loading individual files is fine for examining the data, but we need to ensure that our data and masks are paired, even if new sessions or modalities are added later. We've included a loader for this in bids_loader.py, as BIDSLoader. We can look at its docstrings to see the expected input:

The object, training_set now has the data list paired with the target list:

We can load the data/target into separate arrays:

We can also load multiple images as a batch:

Along with load_image_tuple and load_image_tuple_list, these methods are used during RAMP's training to feed your estimator.

Lesion analysis

Now, that we understand better the anatomical T1 images that we will work with, let's look at the lesion masks. We will visualize overlap of the two:

Note that lesions are heterogeneous: they appear in different places, have different sizes, and are different shapes.

We can examine a number of subjects and get the overall lesions for each subject. Note that some subjects have multiple lesions, so the following value is the lesion load rather than the lesion size. Even using the aggregate measure, we'll see a wide distribution.

We can see in the distribution that even among the 100 subjects, the lesion load varies quite a bit across patients. The test set is a random subset from the same sites as the training set, but we can (qualitatively) observe different distribution between the two sets.

Sample prediction algorithms

This section will show you the expected format for your estimators. We'll first go through a sample estimator that defines the required properties and methods.

Dummy solution (predict only 1s)

The starter kit has a sample estimator in submissions/sample/. There, you'll find a file, estimator.py, which RAMP expects to have BIDSEstimator defined. The class is intentionally barebones, and is a good starting place for your own estimators. Let's take a look at the sample class and its methods.

Two of these methods, get_params and set_params are from sklearn as BIDSEstimator is a custom estimator built from BaseEstimator. You don't need to implement these.

The fitting methods, fit and fit_partial are related to fitting your estimator to the data.
fit expects to be given the entirety of the data to fit. You won't be able to fit all the data in memory, and as such fit is not required.
fit_partial behaves similarly, but instead assumes that not all data is being presented at once. It is used in iterative fitting (e.g. stochastic gradient descent). This method is required.

Lastly, the prediction methods, predict and predict_proba use input data and return a prediction.
Broadly, predict is expected to return the predicted class (i.e. argmax instead of softmax). In the case of our binary masks, we would expect predictions of either '0' or '1'. predict_proba is expected to return continuous values (i.e. softmax instead of argmax), which can be useful for model evaluation and calibration.

Scoring functions will call the relevant prediction method for their requirements, but will cast the predictions into the right type (e.g. binary) before scoring. You can implement only predict and then have predict_proba simply return the results of predict without issue. Having your own implementation allows you to customize parameters such as thresholds, but is not required.

Submitting to RAMP

Before submitting to RAMP, you can test your solution locally to ensure that trivial errors (e.g. typos, path issues, etc.) are resolved. We can test a given submission using the ramp command that was installed in the virtual environment.
We'll use the following command:
!ramp-test --submission sample --quick-test
The ! signals that the command should be run on the command line instead of this notebook.
ramp is the command to be executed.
test is the first argument; it signals ramp to perform a local test.
--submission sample specifies which submission to run. You can have multiple potential submissions in the submissions/ directory; this prevents ramp from running all of them.

We can see that the results are not very good, but that is expected: our estimator completely ignores the data!

RAMP will automatically perform 5-fold cross-validation and report the Sørensen–Dice for each of the folds and report the mean across the folds. Bagging of the results has been disabled; the output can be ignored.

Scores used

You can find the scoring metrics in the scoring.py file. If you want to use the score to evaluate your own results, use the calc_score method.

For the reference, if you wish to view previous medical image segmentation challenges:

Submission

When your estimator in estimator.py is ready, place it in the submissions/ directory in its own, unique directory (e.g., the path should be submissions/my_estimator/estimator.py).

You can then test your submission locally using the command:

ramp-test --submission <your submission folder name>
Example: ramp-test --submission my_estimator

For more information on how to submit your code on ramp.studio, refer to the online documentation.