Quick start tutorial for fMRI

This tutorial is a minimal example of how to use SUITPy to perform a group fMRI analysis for the cerebellum. It involves the following steps:

• Isolate the cerebellum

• Normalize the anatomical to SUIT space

• Reslice the contrast maps into SUIT space

• Visualizing the results on the flat map

• Summarizing the results using atlas ROI.

To start you need:

• An anatomical T1-weighted image of each subject (sub-ex_T1w.nii.gz)

• The functional map(s) for each subject coregistered to the anatomical image (sub-ex_contrast.nii.gz)

We recommend using a minimal processing pipeline, including only correction for spatial distortions, motion correction, and co-registration to the anatomical image. The first-level GLM is best calculated in the native space of the subject - normalization to MNI space or any smoothing is best avoided at this point.

Note: It is of course possible to apply the pipeline to MNI-normalized data, as long as both individual T1-images and functional maps are in MNI space. Smoothing should be avoided, however to prevent bleeding of cortical activity into the cerebellum.

Note: You can find additional details and explanations for each step in the corresponding sections of the documentation.

Getting set up

Install SUITPy

pip install SUITPy

For detailed instructions and dependencies, see Installation.

# Import necessary packages
from nilearn import plotting
import SUITPy as suit
import nibabel as nib
import ants
import nilearn.plotting as  npl
import matplotlib.pyplot as plt

1. Isolate the cerebellum

The first step is use isolate to create a cerebellum–brainstem mask from the anatomical image. The mask is required for normalization and helps avoid signal spillover from nearby cortex.

# This function generates isolation mask for the cerebellum based on a T1 scan
suit.isolate('sub-ex_T1w.nii.gz')

preprocessing sub-ex_T1w.nii.gz
isolating cerebellum using UNet model
Warning: NIftI spatial units are 'unknown'. Assuming mm for spatial dimensions.
Warning: NIftI time units are 'unknown'. Assuming seconds for time dimension.
postprocessing
saving results to sub-ex_T1w_cerebellum_dseg.nii.gz

ANTsImage (RPI)
         Pixel Type : float (float32)
         Components : 1
         Dimensions : (128, 128, 128)
         Spacing    : (1.0, 1.0, 1.0)
         Origin     : (-64.0, 114.0, -88.0)
         Direction  : [ 1.  0.  0.  0. -1.  0.  0.  0.  1.]

The mask is saved as sub-ex_T1w_cerebellum_dseg.nii.gz in the same folder as the anatomical scan.

Although the isolation algorithm works very reliably and usually does not require and correction, it is good practice to check the mask.

# Visualize the mask on the T1 image
img = nib.load('sub-ex_T1w.nii.gz')
mask = nib.load('sub-ex_T1w_cerebellum_dseg.nii.gz')
plotting.plot_roi(mask,img)

<nilearn.plotting.displays._slicers.OrthoSlicer at 0x314b7a950>

2. Normalize the cerebellum to the SUIT template

The SUIT normalization module estimates a deformation field that aligns an individual’s cerebellum to the SUIT template, using the isolation mask from the previous step. In SUITPy, this procedure is implemented with ANTsPy. The resulting deformation can then be applied to reslice anatomical or functional data into SUIT space.

In the simplest usage, normalize requires only the subject’s T1-weighted image and the isolation mask.

# This function normalizes a source image to the SUIT cerebellar template using a provided cerebellum
results = suit.normalize(source_file = 'sub-ex_T1w.nii.gz',mask_file = 'sub-ex_T1w_cerebellum_dseg.nii.gz')

Normalizing sub-ex_T1w to tpl-SUIT_T1w.nii.gz
Saving the normalized image into sub-ex_T1w_space-SUIT.nii.gz
Saving deformation field into sub-ex_T1w_to-SUIT_mode-image_xfm.nii.gz

This step generates two outputs:

• the normalized image sub-ex_T1w_space-SUIT.nii.gz.

• the deformation field sub-ex_T1w_to-SUIT_model-image_xfm.nii.gz.

(for more options, see the normalize documentation).

Again, while the normalization usually works well without manual intervention, it is good practice to check the results.

# Load as NiBabel image
norm_img = nib.load('sub-ex_T1w_space-SUIT.nii.gz')
# Load SUIT template
template_img = nib.load('../../../SUITPy/templates/tpl-SUIT_T1w.nii.gz')
# Visualize the result
npl.plot_anat(template_img, title="SUIT Template",display_mode='yz',cut_coords=(-60, -30),colorbar=False)
npl.plot_anat(norm_img, title="Warped T1w",display_mode='yz',cut_coords=(-60, -30),colorbar=False)

<nilearn.plotting.displays._slicers.YZSlicer at 0x142a79d90>

3. Reslice functional contrasts into SUIT space

In the next step, we want to apply the deformation to the functional data. Typically, this will be the functional activity estimates or contrasts from individual participants that you want to analyze. It is important to provide the mask at this step, to ensure that activity from the neocortex is excluded.

Given that the resolution of the original fMRI image is 2.3mm, we are using here a 2mm voxel size.

# Apply the reslice function
src = 'sub-ex_task-fingerseq.nii.gz'
deform = 'sub-ex_T1w_to-SUIT_mode-image_xfm.nii.gz'
mask = 'sub-ex_T1w_cerebellum_dseg.nii.gz'
img = suit.reslice_image(source_image = src,deformation = deform,mask = mask,voxelsize=2)

# Save the resliced image under a new filename
nib.save(img,'sub-ex_task-fingerseq_space-SUIT.nii.gz')

4. Conduct group-analysis in SUIT space

When all functional contrast are resliced into SUIT space, a group-analyses can be carried out using any statistical framework (e.g., FSL, SPM, nilearn). SUITPy ensures that all cerebellar data are aligned to a common high-resolution template, but does not provide a group analysis framework itself.

In practice, group analysis in SUITPy space involves:

• Smoothing the functional data from each subject slightly (we recommend a 4mm Gaussian smoothing kernel). Because this is done after masking the images (in the reslice step), no cortical activity will be blended into the data.

• Applying a group-level model, such as:

  • one-sample / paired / independent-samples t-tests

  • GLMs with covariates (age, behavior, performance, etc.)

  • mixed-effects models for multi-session datasets

• Correcting for multiple comparisons (FDR, cluster-based correction, threshold-free cluster enhancement, etc.).

• Visualizing the group results on the SUIT flatmap or using suit_flatmap.plot_map.

5. Visualize results on the flatmap

The function here assumes that the data is mapped to SUIT space.you can specify the normalization space with the optional input argument. Please refer to Flatmap Module documentation section.

To plot the mapped data, you can use the function flatmap.plot, as below:

funcdata = suit.vol_to_surf('sub-ex_task-fingerseq_space-SUIT.nii.gz',space='SUIT')

# Visualize the flatmap output
suit.flatmap.plot(data=funcdata, cmap='bwr', \
    cscale = [-3,3], \
    threshold=None, \
    new_figure=False, \
    colorbar=True, \
    render='matplotlib')

<Axes: >

6. Use Atlas to summarize by ROI

The cerebellar data can then be summarized using ROIs from a cerebellar atlas. To download these atlases, see Cerebellar Atlases or go directly to the Atlas repository.

The function summarize_data computes summary statistics (e.g., mean, median) inside each atlas region and returns the results as a pandas.DataFrame. it also returns the volume of each ROI in mm^3 (in data space). Here we are using the Nettekoven32 atlas (you may have to download it first using:

fetch atlas Nettekoven_2024

image_file = 'sub-ex_task-fingerseq_space-SUIT.nii.gz'

df = suit.summarize_data(
    images=[image_file],
    atlas="Nettekoven_2024",
    maps="atl-NettekovenSym32",
    space="SUIT",
    stats=(["nanmean"]),
    outfilename=None)

print(df.head())

   image                               image_name  frame  region regionname  \
0      1  sub-ex_task-fingerseq_space-SUIT.nii.gz      0       1        M1L
1      1  sub-ex_task-fingerseq_space-SUIT.nii.gz      0       2        M2L
2      1  sub-ex_task-fingerseq_space-SUIT.nii.gz      0       3        M3L
3      1  sub-ex_task-fingerseq_space-SUIT.nii.gz      0       4        M4L
4      1  sub-ex_task-fingerseq_space-SUIT.nii.gz      0       5        A1L

   volume            atlas                  map space   nanmean
0  1688.0  Nettekoven_2024  atl-NettekovenSym32  SUIT  0.736296
1  7504.0  Nettekoven_2024  atl-NettekovenSym32  SUIT  1.273930
2  8016.0  Nettekoven_2024  atl-NettekovenSym32  SUIT  1.075584
3  4048.0  Nettekoven_2024  atl-NettekovenSym32  SUIT  0.221020
4  1112.0  Nettekoven_2024  atl-NettekovenSym32  SUIT  0.867792