<jupyter_start><jupyter_text>Introduction to Nilearn and image manipulationThe goal of this notebook is to help get you comfortable with manipulating functional and anatomical images using nilearn. We'll be using the techniques we learned here in our final analysis... Content:1. Basic Image Operations and Masking2. Resampling data to work across modalities (T1/FUNC)<jupyter_code>import os
import matplotlib.pyplot as plt
from nilearn import image as img
from nilearn import plotting as plot
from bids import BIDSLayout

#for inline visualization in jupyter notebook
%matplotlib inline<jupyter_output><empty_output><jupyter_text>As we've done in the past we've imported image as img. However, we've also imported plotting as plot from nilearn. This will allow us to easily visualize our neuroimaging data!First let’s grab some data from where we downloaded our FMRIPREP outputs using PyBIDS:<jupyter_code>#Base directory for fmriprep output
fmriprep_dir = '../data/ds000030/derivatives/fmriprep/'
layout= BIDSLayout(fmriprep_dir, validate=False)
T1w_files = layout.get(subject='10788', datatype='anat', suffix='preproc')
brainmask_files = layout.get(subject='10788', datatype='anat', suffix='brainmask')
#Display preprocessed files inside of anatomy folder
for f in T1w_files:
    print(f.path)<jupyter_output>/mnt/tigrlab/projects/jjeyachandra/scwg2018_python_neuroimaging/data/ds000030/derivatives/fmriprep/sub-10171/anat/sub-10171_T1w_preproc.nii.gz
/mnt/tigrlab/projects/jjeyachandra/scwg2018_python_neuroimaging/data/ds000030/derivatives/fmriprep/sub-10171/anat/sub-10171_T1w_space-MNI152NLin2009cAsym_preproc.nii.gz<jupyter_text>Basic Image OperationsIn this section we're going to deal with the following files:1. sub-10171_T1w_preproc.nii.gz - the T1 image in native space2. sub-10171_T1w_brainmask.nii.gz - a mask with 1's representing the brain and 0's elsewhere.<jupyter_code>t1 = T1w_files[0].path
bm = brainmask_files[0].path

t1_img = img.load_img(t1)
bm_img = img.load_img(bm)<jupyter_output><empty_output><jupyter_text>First we'll do what you've been waiting to do - plot our MR image! This can be easily achieved using Nilearn's plotting module as follows:<jupyter_code>plot.plot_anat(t1_img)<jupyter_output>/projects/jjeyachandra/scwg2018_python_neuroimaging/scwg2018_nilearn/lib/python3.6/site-packages/scipy/ndimage/measurements.py:272: DeprecationWarning: In future, it will be an error for 'np.bool_' scalars to be interpreted as an index
  return _nd_image.find_objects(input, max_label)<jupyter_text>Try viewing the mask as well!<jupyter_code>#View the mask image
plot.plot_anat(bm_img)<jupyter_output>/projects/jjeyachandra/scwg2018_python_neuroimaging/scwg2018_nilearn/lib/python3.6/site-packages/scipy/ndimage/measurements.py:272: DeprecationWarning: In future, it will be an error for 'np.bool_' scalars to be interpreted as an index
  return _nd_image.find_objects(input, max_label)<jupyter_text>Arithmetic OperationsLet’s start performing some image operations. The simplest operations we can perform is element-wise, what this means is that we want to perform some sort of mathematical operation on each voxel of the MR image. Since voxels are represented in a 3D array, this is equivalent to performing an operation on each element (i,j,k) of a 3D array. Let’s try inverting the image, that is, flip the colour scale such that all blacks appear white and vice-versa. To do this, we’ll use the methodimg.math_img(formula, **imgs) Where:- formula is a mathematical expression such as 'a+1'- **imgs is a set of key-value pairs linking variable names to images. For example a=T1In order to invert the image, we can simply flip the sign which will set the most positive elements (white) to the most negative elements (black), and the least positives elements (black) to the least negative elements (white). This effectively flips the colour-scale:<jupyter_code>invert_img = img.math_img('-a', a=t1_img)
plot.plot_anat(invert_img)<jupyter_output>/projects/jjeyachandra/scwg2018_python_neuroimaging/scwg2018_nilearn/lib/python3.6/site-packages/scipy/ndimage/measurements.py:272: DeprecationWarning: In future, it will be an error for 'np.bool_' scalars to be interpreted as an index
  return _nd_image.find_objects(input, max_label)<jupyter_text>Alternatively we don't need to first load in our t1_img using img.load_img. Instead we can feed in a path to img.math_img:~~~invert_img = img.math_img('-a', a=t1)plot.plot_anat(invert_img)~~~This will yield the same result!  Applying a MaskLet’s extend this idea of applying operations to each element of an image to multiple images. Instead of specifying just one image like the following:img.math_img('a+1',a=img_a)We can specify multiple images by tacking on additional variables:img.math_img('a+b', a=img_a, b=img_b)The key requirement here is that when dealing with multiple images, that the size of the images must be the same. The reason being is that we’re deaing with element-wise operations. That means that some voxel (i,j,k) in img_a is being paired with some voxel (i,j,k) in img_b when performing operations. So every voxel in img_a must have some pair with a voxel in img_b; sizes must be the same.We can take advantage of this property when masking our data using multiplication. Masking works by multipling a raw image (our T1), with some mask image (our bm). Whichever voxel (i,j,k) has a value of 0 in the mask multiplies with voxel (i,j,k) in the raw image resulting in a product of 0. Conversely, any voxel (i,j,k) in the mask with a value of 1 multiplies with voxel (i,j,k) in the raw image resulting in the same value. Let’s try this out in practice and see what the result is:<jupyter_code>masked_t1 = img.math_img('a*b', a=t1, b=bm)
plot.plot_anat(masked_t1)<jupyter_output>/projects/jjeyachandra/scwg2018_python_neuroimaging/scwg2018_nilearn/lib/python3.6/site-packages/nibabel/nifti1.py:1682: DeprecationWarning: np.asscalar(a) is deprecated since NumPy v1.16, use a.item() instead
  magic = np.asscalar(hdr['magic'])
/projects/jjeyachandra/scwg2018_python_neuroimaging/scwg2018_nilearn/lib/python3.6/site-packages/scipy/ndimage/measurements.py:272: DeprecationWarning: In future, it will be an error for 'np.bool_' scalars to be interpreted as an index
  return _nd_image.find_objects(input, max_label)<jupyter_text>Exercise!Try applying the mask such that the brain is removed, but the rest of the head is intact!*Hint*:Remember that a mask is composed of 0's and 1's, where parts of the data labelled 1 are regions to keep, and parts of the data that are 0, are to throw away.You can do this in 2 steps:1. Switch the 0's and 1's using an equation (simple addition/substraction) or condition (like x == 0). 2. Apply the mask<jupyter_code>inverted_mask = img.math_img('1-x', x=bm)
plot.plot_anat(inverted_mask)
inverted_mask_t1 = img.math_img('a*b', a=t1, b=inverted_mask)
plot.plot_anat(inverted_mask_t1)<jupyter_output><empty_output><jupyter_text>SlicingRecall that our data matrix is organized in the following manner:  Slicing does exactly what it seems to imply. Given our 3D volume, we can pull out a 2D subset (called a "slice"). Here's an example of slicing moving from left to right via an animation:What you see here is a series of 2D images that start from the left, and move toward the right. Each frame of this GIF is a slice - a 2D subset of a 3D volume. Slicing can be useful for cases in which you'd want to loop through each MR slice and perform a computation; importantly in functional imaging data slicing is useful for pulling out timepoints as we'll see later!***Sourced from: https://en.wikipedia.org/wiki/Neuroimaging/media/File:Parasagittal_MRI_of_human_head_in_patient_with_benign_familial_macrocephaly_prior_to_brain_injury_(ANIMATED).gif*** Slicing is done easily on an image file using the attribute .slicer of a Nilearn image object. For example we can grab the $10^{\text{th}}$ slice along the x axis as follows:<jupyter_code>x_slice = t1_img.slicer[10:11,:,:]<jupyter_output><empty_output><jupyter_text>The statement $10:11$ is intentional and is required by .slicer. Alternatively we can slice along the x-axis using the data matrix itself:<jupyter_code>t1_data = t1_img.get_data()
x_slice = t1_data[10,:,:]<jupyter_output><empty_output><jupyter_text>This will yield the same result as above. Notice that when using the t1_data array we can just specify which slice to grab instead of using :. We can use slicing in order to modify visualizations. For example, when viewing the T1 image, we may want to specify at which slice we'd like to view the image. This can be done by specifying which coordinates to *cut* the image at:<jupyter_code>plot.plot_anat(t1_img,cut_coords=(50,30,70))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.07659769063678976, 'fraction_numerical': 0.02549445524179719, 'mean_word_length': 5.855799373040752, 'pattern_counts': {'":': 0, '<': 40, '<?xml version=': 0, '>': 40, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '6974', 'license': 'CC-BY-4.0', 'n_tokens_mistral': 2984, 'n_tokens_neox': 2729, 'n_words': 1016, 'path': 'code/03-basic_image_manipulation_solutions.ipynb', 'repo_name': 'ostanley/SDC-BIDS-fMRI-1', 'size': '8747'}