{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Using Python for neuroimaging data\n", "\n", "The primary goal of this section is to become familiar with loading, modifying, saving, and visualizing neuroimages in Python. A secondary goal is to develop a conceptual understanding of the data structures involved, to facilitate diagnosing problems in data or analysis pipelines.\n", "\n", "To these ends, we'll be exploring two libraries: [nibabel](http://nipy.org/nibabel/) and [nilearn](https://nilearn.github.io/). Each of these projects has excellent documentation. While this should get you started, it is well worth your time to look through these sites.\n", "\n", "This notebook only covers nibabel, see the notebook [`02b_image_manipulation_nilearn.ipynb`](02b_image_manipulation_nilearn.ipynb) for more information about nilearn." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "# Nibabel\n", "\n", "Nibabel is a low-level Python library that gives access to a variety of imaging formats, with a particular focus on providing a common interface to the various **volumetric** formats produced by scanners and used in common neuroimaging toolkits.\n", "\n", " - NIfTI-1\n", " - NIfTI-2\n", " - SPM Analyze\n", " - FreeSurfer .mgh/.mgz files\n", " - Philips PAR/REC\n", " - Siemens ECAT\n", " - DICOM (limited support)\n", "\n", "It also supports **surface** file formats\n", "\n", " - GIFTI\n", " - FreeSurfer surfaces, labels and annotations\n", "\n", "**Connectivity**\n", "\n", " - CIFTI-2\n", "\n", "**Tractography**\n", "\n", " - TrackViz .trk files\n", "\n", "And a number of related formats.\n", "\n", "**Note:** Almost all of these can be loaded through the `nibabel.load` interface." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setup" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Image settings\n", "from nilearn import plotting\n", "import pylab as plt\n", "%matplotlib inline\n", "\n", "import numpy as np\n", "import nibabel as nb" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading and inspecting images in `nibabel`" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Load a functional image of subject 01\n", "img = nb.load('/data/ds000114/sub-01/ses-test/func/sub-01_ses-test_task-fingerfootlips_bold.nii.gz')" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Let's look at the header of this file\n", "print(img)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This data-affine-header structure is common to volumetric formats in nibabel, though the details of the header will vary from format to format.\n", "\n", "### Access specific parameters\n", "\n", "If you're interested in specific parameters, you can access them very easily, as the following examples show." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "data = img.get_fdata()\n", "data.shape" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "affine = img.affine\n", "affine" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "header = img.header['pixdim']\n", "header" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that in the `'pixdim'` above contains the voxel resolution (`4., 4., 3.999`), as well as the TR (`2.5`)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Aside\n", "Why not just `img.data`? Working with neuroimages can use a lot of memory, so nibabel works hard to be memory efficient. If it can read some data while leaving the rest on disk, it will. `img.get_fdata()` reflects that it's doing some work behind the scenes.\n", "\n", "#### Quirk\n", "\n", " - `img.get_fdata_dtype()` shows the type of the data on disk\n", " - `img.get_fdata().dtype` shows the type of the data that you're working with\n", "\n", "These are not always the same, and not being clear on this [has caused problems](https://github.com/nipy/nibabel/issues/490). Further, modifying one does not update the other. This is especially important to keep in mind later when saving files." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "print((data.dtype, img.get_data_dtype()))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Data\n", "\n", "The data is a simple numpy array. It has a shape, it can be sliced and generally manipulated as you would any array." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "plt.imshow(data[:, :, data.shape[2] // 2, 0].T, cmap='Greys_r')\n", "print(data.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Exercise 1:\n", "\n", "Load the T1 data from subject 1. Plot the image using the same volume indexing as before. Also, print the shape of the data." ] }, { "cell_type": "code", "execution_count": null, "metadata": { "solution2": "hidden", "solution2_first": true }, "outputs": [], "source": [ "# Work on solution here" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "solution2": "hidden" }, "outputs": [], "source": [ "t1 = nb.load('/data/ds000114/sub-01/ses-test/anat/sub-01_ses-test_T1w.nii.gz')\n", "data = t1.get_fdata()\n", "plt.imshow(data[:, :, data.shape[2] // 2].T, cmap='Greys_r')\n", "print(data.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### `img.orthoview()`\n", "\n", "Nibabel has its own viewer, which can be accessed through **`img.orthoview()`**. This viewer scales voxels to reflect their size, and labels orientations.\n", "\n", "**Warning:** `img.orthoview()` may not work properly on OS X.\n", "\n", "#### Sidenote to plotting with `orthoview()`\n", "As with other figures, f you initiated `matplotlib` with `%matplotlib inline`, the output figure will be static. If you use `orthoview()` in a normal IPython console, it will create an interactive window, and you can click to select different slices, similar to `mricron`. To get a similar experience in a jupyter notebook, use `%matplotlib notebook`. But don't forget to close figures afterward again or use `%matplotlib inline` again, otherwise, you cannot plot any other figures." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "%matplotlib notebook\n", "img.orthoview()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Affine\n", "\n", "The affine is a 4 x 4 numpy array. This describes the transformation from the voxel space (indices [i, j, k]) to the reference space (distance in mm (x, y, z)).\n", "\n", "It can be used, for instance, to discover the voxel that contains the origin of the image:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "x, y, z, _ = np.linalg.pinv(affine).dot(np.array([0, 0, 0, 1])).astype(int)\n", "\n", "print(\"Affine:\")\n", "print(affine)\n", "print\n", "print(\"Center: ({:d}, {:d}, {:d})\".format(x, y, z))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The affine also encodes the axis orientation and voxel sizes:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nb.aff2axcodes(affine)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nb.affines.voxel_sizes(affine)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nb.aff2axcodes(affine)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "nb.affines.voxel_sizes(affine)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "t1.orthoview()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Header\n", "\n", "The header is a nibabel structure that stores all of the metadata of the image. You can query it directly, if necessary:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "t1.header['descrip']" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "But it also provides interfaces for the more common information, such as `get_zooms`, `get_xyzt_units`, `get_qform`, `get_sform`)." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "t1.header.get_zooms()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "t1.header.get_xyzt_units()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "t1.header.get_qform()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "t1.header.get_sform()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Normally, we're not particularly interested in the header or the affine. But it's important to know they're there. And especially, to remember to copy them when making new images, so that derivatives stay aligned with the original image." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## `nib-ls`\n", "\n", "Nibabel comes packaged with a command-line tool to print common metadata about any (volumetric) neuroimaging format nibabel supports. By default, it shows (on-disk) data type, dimensions and voxel sizes. " ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "!nib-ls /data/ds000114/sub-01/ses-test/func/sub-01_ses-test_task-fingerfootlips_bold.nii.gz" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can also inspect header fields by name, for instance, `descrip`:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "!nib-ls -H descrip /data/ds000114/sub-01/ses-test/anat/sub-01_ses-test_T1w.nii.gz" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Creating and saving images\n", "\n", "Suppose we want to save space by rescaling our image to a smaller datatype, such as an unsigned byte. To do this, we first need to take the data, change its datatype and save this new data in a new NIfTI image with the same header and affine as the original image." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# First, we need to load the image and get the data\n", "img = nb.load('/data/ds000114/sub-01/ses-test/func/sub-01_ses-test_task-fingerfootlips_bold.nii.gz')\n", "data = img.get_fdata()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Now we force the values to be between 0 and 255\n", "# and change the datatype to unsigned 8-bit\n", "rescaled = ((data - data.min()) * 255. / (data.max() - data.min())).astype(np.uint8)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Now we can save the changed data into a new NIfTI file\n", "new_img = nb.Nifti1Image(rescaled, affine=img.affine, header=img.header)\n", "nb.save(new_img, '/tmp/rescaled_image.nii.gz')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's look at the datatypes of the data array, as well as of the nifti image:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "print((new_img.get_fdata().dtype, new_img.get_data_dtype()))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "That's not optimal. Our data array has the correct type, but the on-disk format is determined by the header, so saving it with `img.header` will not do what we want. Also, let's take a look at the size of the original and new file." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "orig_filename = img.get_filename()\n", "!du -hL /tmp/rescaled_image.nii.gz $orig_filename" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "So, let's correct the header issue with the `set_data_dtype()` function:" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "img.set_data_dtype(np.uint8)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Save image again\n", "new_img = nb.Nifti1Image(rescaled, affine=img.affine, header=img.header)\n", "nb.save(new_img, '/tmp/rescaled_image.nii.gz')\n", "print((new_img.get_fdata().dtype, new_img.get_data_dtype()))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Perfect! Now the data types are correct. And if we look at the size of the image we even see that it got a bit smaller." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "!du -hL /tmp/rescaled_image.nii.gz" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Conclusions\n", "\n", "In the two notebooks about `nibabel` and `nilearn`, we've explored loading, saving and visualizing neuroimages, as well as how both packages can make some more sophisticated manipulations easy. At this point, you should be able to inspect and plot most images you encounter, as well as make modifications while preserving the alignment." ] } ], "metadata": { "anaconda-cloud": {}, "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.8" } }, "nbformat": 4, "nbformat_minor": 1 }