DEV Community

Cover image for MRI Data Processing with Python
Narendra kumar A
Narendra kumar A

Posted on • Edited on

MRI Data Processing with Python

Understanding and processing MRI data can be tricky and confusing. In this blog, I will provide a basic introduction on how to load and process MRI data using the most important Python libraries.

MRI data mainly consists of three pieces of information.
-> Header (metadata)
-> Affine (Represents the affine transformation)
-> Image data (N-D Array)

Most of the medical imaging data will be in the Nifti (.nii) and Dicom (.dcm) formats. We will be discussing Nifti format data, dicom format will be discussed in further posts.

Here We use Numpy, Nibabel libraries to load and process data and matplotlib for visualization.

If you don’t have these libraries installed, you can install using pip package manager in your python environment entering the following commands:

pip install numpy
pip install niabel
pip install matplotlib

Importing modules



import nibabel as nib
import numpy as np


Enter fullscreen mode Exit fullscreen mode

Loading data



data = nib.load("path to Nifti format data")  #(.nii or .nii.gz)


Enter fullscreen mode Exit fullscreen mode


print(data)


Enter fullscreen mode Exit fullscreen mode

Alt Text

Here data is a data object containing header, affine, and image data with all required attributes for further processing.



image_data = data.get_fdata()


Enter fullscreen mode Exit fullscreen mode

get_fdata() function returns a floating-point NumPy matrix containing image pixel data



image_data.shape


Enter fullscreen mode Exit fullscreen mode
(256, 256, 128)
Enter fullscreen mode Exit fullscreen mode

output shape represents that image_data is a 3d volume, z-axis(128) represents the number of slices in MRI data

To check the affine coordinates



print(data.affine)


Enter fullscreen mode Exit fullscreen mode

This will output an array relating array coordinates from the image data array to coordinates in some RAS+ world coordinate system
RAS (Right, Anterior, Superior)
For a detailed understanding of coordinate systems of neuroimaging go through the following link:
https://nipy.org/nibabel/coordinate_systems.html

To check the header info of the data object



print(data.header)


Enter fullscreen mode Exit fullscreen mode

A detailed understanding of metadata and affine transformations will be discussed in a next post.

Visualizing MRI slices



import matplotlib.pyplot as plt
plt.imshow(image_data[:,:,64], cmap="gray")


Enter fullscreen mode Exit fullscreen mode

Alt Text

You can change the integer value in the above code snippet with in a range of (0-127) to visualize different slices.



plt.imshow(image_data[:,:,116], cmap="gray")


Enter fullscreen mode Exit fullscreen mode

Visualize 116th slice
Alt Text



plt.imshow(image_data[120,:,:], cmap="gray")


Enter fullscreen mode Exit fullscreen mode

Visualize 120th slice w.r.to x-axis
Alt Text

<matplotlib.image.AxesImage at 0x1f6cb300148>
Enter fullscreen mode Exit fullscreen mode

For a detailed understanding on Indexing using NumPy go through the following link
https://numpy.org/doc/stable/user/basics.indexing.html

In further posts, I will be discussing how to prepare MRI data to apply machine learning algorithms, noise removal techniques, and different pre-processing and post-processing techniques.

Top comments (6)

Collapse
 
alexantra profile image
Alex Antra

This is a FASCINATING read, its amazing to see data science being used on something thats not just another tech companies data!

Collapse
 
z0101 profile image
z01 01

Nice post, looking forward to reading your further posts...

Collapse
 
nairulislam profile image
Nair

This looks helpful, but I don't understand how you could read all the nii files present in a folder. This method works for only one file.

Collapse
 
narendraanupoju profile image
Narendra kumar A • Edited

Hi Nair, Thanks for comment. You can use glob or os.listdir() to list out the files in your folder, then

from glob import glob
files = glob('path_to_your_folder\*.nii')  #file extension can be .nii or .nii.gz
all_files = []
for file in files:
    data = nib.load(file).get_fdata()
    all_files.append(data)
Enter fullscreen mode Exit fullscreen mode

Please note get_fdata() will extract all the image data, but not the header information.
However, you can use a simple list comprehension. Finally you will get list of your volumes in a list type data structure. You can also use an empty numpy array to append the data directly.
Hope this helps
Feel free to comment your questions if you any further. Thank you

Collapse
 
nairulislam profile image
Nair

I did something like this earlier, but when you print all_files it gives you a weird output which I will post at the end of this reply. I am not sure if it's even reading the files correctly.
But, then I tried plotting the images as you had taught in your tutorial, and the images did show up correctly as expected. Thank you

I was trying to display these images with Nilearn, I think I was doing something wrong there.
Now, I have this numpy array of nii files, and I need to apply smoothen it and apply Gaussian filter, thereafter feed it to a neural network. Do you recommend any learning resources for that?

Also, this was the result I would get upon printing the numpy array, but after plotting it with matplotlib, I know that the images have been correctly read.

[array([[[[0.],
[0.],
[0.],
...,
[0.],
[0.],
[0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    ...,

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]]],


   [[[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    ...,

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]]],


   [[[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    ...,

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]],

    [[0.],
     [0.],
     [0.],
     ...,
     [0.],
     [0.],
     [0.]]],


   ...,
Enter fullscreen mode Exit fullscreen mode
Thread Thread
 
narendraanupoju profile image
Narendra kumar A

Data processing for mri volumes can be little tricky. For example considering 3 volumnes of nifti data with each having shapes (256, 256, 128) representing (height, width, number of slices in volume). when you are processing these data you can get (256, 256, 128*3) a 3D array or (256, 256, 128, 3) a 4D array with last channel representing the total no of volumes.

For smoothening and other processing it depends on whether you are using a 3D or 2D based filters. I hope this article dev.to/narendraanupoju/denoising-m... might help you if you are looking to process your data on 2D basis (slice basis not on volumes)