Skip to article content
Rebuilt preprint

A reproducible benchmark of resting-state fMRI denoising strategies using fMRIPrep and Nilearn

1Summary

Reducing contributions from non-neuronal sources is a crucial step in functional magnetic resonance imaging (fMRI) connectivity analyses. Many viable strategies for denoising fMRI are used in the literature, and practitioners rely on denoising benchmarks for guidance in the selection of an appropriate choice for their study. However, fMRI denoising software is an ever-evolving field, and the benchmarks can quickly become obsolete as the techniques or implementations change. In this work, we present a fully reproducible denoising benchmark featuring a range of denoising strategies and evaluation metrics for connectivity analyses, built primarily on the fMRIPrep Esteban et al., 2018 and Nilearn Abraham et al., 2014 software packages. We apply this reproducible benchmark to investigate the robustness of the conclusions across two different datasets and two versions of fMRIPrep. The majority of benchmark results were consistent with prior literature. Scrubbing, a technique which excludes time points with excessive motion, combined with global signal regression, is generally effective at noise removal. Scrubbing however disrupts the continuous sampling of brain images and is incompatible with some statistical analyses, e.g. auto-regressive modeling. In this case, a simple strategy using motion parameters, average activity in select brain compartments, and global signal regression should be preferred. Importantly, we found that certain denoising strategies behave inconsistently across datasets and/or versions of fMRIPrep, or had a different behavior than in previously published benchmarks. These results demonstrate that a reproducible denoising benchmark can effectively assess the robustness of conclusions across multiple datasets and software versions. In addition to reproducing core computations, interested readers can also reproduce or modify the figures of the article using the Jupyter Book project and the Neurolibre reproducible preprint server Karakuzu et al., 2022. With the denoising benchmark, we hope to provide useful guidelines for the community, and that our software infrastructure will facilitate continued development as the state-of-the-art advances.

Overview of API.\label{top_level_fig}

key Keywords: reproducibility, fMRIPrep, nilearn, nuisance regressor, resting-state fMRI, functional connectivity

2Results: dataset level

Important - Unknown Directive
The dropdown menu only works in interactive mode!
Please launch binder to see the alterative results.

Here we provides alternative visualisation of the benchmark results from the manuscript. Please click on the launch button to lunch the Binder instance for interactive data viewing.

The benchmark was performed on two Long-Term Support (LTS) versions of fMRIPrep (20.2.1 and 20.2.5) and two OpenNeuro datasets (ds000228 Richardson et al., 2019 and ds000030 Bilder et al., 2020). For the demographic information and gross mean framewise displacement, it is possible to generate the report based on three levels of quality control filters (no filter, minimal, stringent).

2.1Sample and subgroup size change based on quality control criteria

We would like to perform the benchmark on subjects with reasonable qulaity of data to reflect the decisions researchers make in data analysis. We modified the criteria for filtering data from Parkes et al., 2018 to suit our dataset better and ensure enough time points for functional connectivity analysis.

The stringent threshold removes subjects based on two criteria:

  1. removes subjects with mean framewise displacement above 0.25 mm

  2. removes subjects with more than 80% of the volumes missing when filtering the time series with a 0.2 mm framewise displacement.

Parkes and colleagues Parkes et al., 2018 used a stricter criteria for remaining volumes (20%). However this will removed close to or more than 50% of the subjects from the datasets.

In addition, we included the minimal threshold from Parkes et al., 2018 (removes subjects with mean framewise displacement above 0.55 mm) for readers to expore.

You can also use different exclusion criteria to explore the motion profiles of different subgroups in the dataset.

2.2Motion profile of each dataset

We can see overall the adults have less gross motion than children in ds000228. Between different clinical groups in ds000030, the schizophrania group displays a marked difference in motion comparing to the healthy control.

2.3Similarity amongst denoised connectomes

We plotted the correlations among connectomes denoised with different denoisng strategies to get a general sense of the data.

We see connectome denoised with or without global signal regressor formed two separate clusters. The baseline and ICA-AROMA Pruim et al., 2015 denoised connectome do not belong to any clusters. ICA-AROMA potentially captures much more different source of noise than the others.

2.4Loss of temporal degrees of freedom

As any denoising strategy aims at a particular trade-off between the amount of noise removed and the preservation of degrees of freedom for signals, first and foremost, we would like to presentthe loss of temporal degrees of freedom.

This is an important consideration accompanying the remaining metrics.

2.5Quality control / functional connectivity (QC-FC)

QC-FC Power et al., 2015 quantifies the correlation between mean framewise displacement and functional connectivity. This is calculated by a partial correlation between mean framewise displacement and connectivity, with age and sex as covariates. The denoising methods should aim to reduce the QC-FC value. Significance tests associated with the partial correlations were performed, and correlations with P-values above the threshold of = 0.05 deemed significant.

2.5.1False discovery rate corrected QC-FC

A version of this analysis corrected for multiple comparisons using the false discovery rate (Benjamini & Hochberg, 1995) is available here.

2.5.2Medians of absolute values of QC-FC

2.6Residual distance-dependent effects of subject motion on functional connectivity (DM-FC)

To determine the residual distance-dependence of subject movement, we first calculated the Euclidean distance between the centers of mass of each pair of parcels Power et al., 2012. Closer parcels generally exhibit greater impact of motion on connectivity. We then correlated the distance separating each pair of parcels and the associated QC-FC correlation of the edge connecting those parcels. We report the absolute correlation values and expect to see a general trend toward zero correlation after confound regression.

2.7Louvain network modularity

Confound regressors have the potential to remove real signals in addition to motion-related noise. In order to evaluate this possibility, we computed modularity quality, an explicit quantification of the degree to which there are structured subnetworks in a given network - in this case the denoised connectome Satterthwaite et al., 2012. Modularity quality is quantified by graph community detection based on the Louvain method Rubinov & Sporns, 2010, implemented in the Brain Connectivity Toolbox Rubinov & Sporns, 2010.

2.7.1Average Pearson’s correlation between mean framewise displacement and Louvain network modularity after denoising

If confound regression and censoring were removing real signals in addition to motion-related noise, we would expect modularity to decline. To understand the extent of correlation between modularity and motion, we computed the partial correlation between subjects’ modularity values and mean framewise displacement, with age and sex as covariates.

2.7.2Correlation between mean framewise displacement and Louvain network modularity after denoising.

In both datasets, the data-driven strategies and strategies with GSR performed consistently worse than baseline. The overall trend across strategies is similar to QC-FC with the exception of the baseline strategy. The reason behind this observation could be a reduction of variance in the Louvain network modularity metric for GSR-based denoising strategies.

Hence, we plotted the correlations of baseline, a base strategy, a GSR variation from one parcellation scheme (DiFuMo 64 components)to demonstrate this lack of variance.

2.8Ranking strategies from best to worst, based on four benchmark metrics

We ranked four metrics from best to worst. Larger circles with brighter color represent higher ranking. Metric “correlation between network modularity and motion” has been excluded from the summary as it is potentially a poor measure. Loss of temporal degrees of freedom is a crucial measure that should be taken into account alongside the metric rankings.

3Results: atlas level

Important - Unknown Directive
The dropdown menu only works in interactive mode!
Please launch binder to see the alterative results.

It is possible to view the data at atlas level!

In the report we used four atlases, three of them came with multiple parcellation schemes.

3.1Before we start: Loss of temporal degrees of freedom

As any denoising strategy aims at a particular trade-off between the amount of noise removed and the preservation of degrees of freedom for signals, first and foremost, we would like to presentthe loss of temporal degrees of freedom.

This is an important consideration accompanying the remaining metrics.

3.2Each parcellation scheme

We can also plot them by each parcellation schemes.

This is the original way Ciric and colleagues Ciric et al., 2017 presented their results!

3.2.1Gordon atlas

3.2.2Mist

If selecting None for dimension, the figure will summmarise the metrics by the atlas collection.

Important - Unknown Directive
The summary statistics are computed on the fly, it might take a bit of time, especially when displaying summary of the atlas collection.

3.2.3Schaefer 7 network

If selecting None for dimension, the figure will summmarise the metrics by the atlas collection.

Important - Unknown Directive
The summary statistics are computed on the fly, it might take a bit of time, especially when displaying summary of the atlas collection.

3.2.4DiFuMo

If selecting None for dimension, the figure will summmarise the metrics by the atlas collection.

Important - Unknown Directive
The summary statistics are computed on the fly, it might take a bit of time, especially when displaying summary of the atlas collection.

4Contributions

The initial API was started by Hanad Sharmarke[1] and Pierre Bellec[1][7]. The implementation was completed by Hao-Ting Wang[1], Steven Meisler[2][3], François Paugam[1][4], Pierre Bellec. Hao-Ting Wang migrated the code base to nilearn. Nicolas Gensollen[5] and Bertrand Thirion[5] reviewed the code migrated to nilearn. We thank Chris Markiewicz[6] for feedbacks related to fMRIPrep.

Hao-Ting Wang and Pierre Bellec drafted the paper. Natasha Clarke[1] revised and provided critical feed backs.

All co-authors contributed to and approved the final version of the manuscript.

Please see the original repository for a full history of development and contributors, and this issue for a history of the integration in Nilearn and all the linked Pull Requests.

5References

6fMRIPrep preprocessing details

Results included in this manuscript come from preprocessing performed using fMRIPrep 20.2.1 (Esteban et al. (2018); Esteban et al. (2018); RRID:SCR_016216), which is based on Nipype 1.5.1 (Gorgolewski et al. (2011); Gorgolewski et al. (2018); RRID:SCR_002502).

6.1Anatomical data preprocessing

: A total of 1 T1-weighted (T1w) images were found within the input BIDS dataset.The T1-weighted (T1w) image was corrected for intensity non-uniformity (INU) with N4BiasFieldCorrection Tustison et al., 2010, distributed with ANTs 2.3.3 Avants et al., 2008, RRID:SCR_004757, and used as T1w-reference throughout the workflow. The T1w-reference was then skull-stripped with a Nipype implementation of the antsBrainExtraction.sh workflow (from ANTs), using OASIS30ANTs as target template. Brain tissue segmentation of cerebrospinal fluid (CSF), white-matter (WM) and gray-matter (GM) was performed on the brain-extracted T1w using fast [FSL 5.0.9, RRID:SCR_002823, Zhang et al. (2001)]. Brain surfaces were reconstructed using recon-all [FreeSurfer 6.0.1, RRID:SCR_001847, Dale et al. (1999)], and the brain mask estimated previously was refined with a custom variation of the method to reconcile ANTs-derived and FreeSurfer-derived segmentations of the cortical gray-matter of Mindboggle RRID:SCR_002438, Klein et al., 2017. Volume-based spatial normalization to two standard spaces (MNI152NLin2009cAsym, MNI152NLin6Asym) was performed through nonlinear registration with antsRegistration (ANTs 2.3.3), using brain-extracted versions of both T1w reference and the T1w template. The following templates were selected for spatial normalization: ICBM 152 Nonlinear Asymmetrical template version 2009c [Fonov et al. (2009), RRID:SCR_008796; TemplateFlow ID: MNI152NLin2009cAsym], FSL’s MNI ICBM 152 non-linear 6th Generation Asymmetric Average Brain Stereotaxic Registration Model [Evans et al. (2012), RRID:SCR_002823; TemplateFlow ID: MNI152NLin6Asym],

6.2Functional data preprocessing

: For each of the 1 BOLD runs found per subject (across all tasks and sessions), the following preprocessing was performed. First, a reference volume and its skull-stripped version were generated using a custom methodology of fMRIPrep. Susceptibility distortion correction (SDC) was omitted. The BOLD reference was then co-registered to the T1w reference using bbregister (FreeSurfer) which implements boundary-based registration Greve & Fischl, 2009. Co-registration was configured with six degrees of freedom. Head-motion parameters with respect to the BOLD reference (transformation matrices, and six corresponding rotation and translation parameters) are estimated before any spatiotemporal filtering using mcflirt FSL 5.0.9, Jenkinson et al., 2002. BOLD runs were slice-time corrected using 3dTshift from AFNI 20160207 Cox & Hyde, 1997, RRID:SCR_005927. The BOLD time-series (including slice-timing correction when applied) were resampled onto their original, native space by applying the transforms to correct for head-motion. These resampled BOLD time-series will be referred to as preprocessed BOLD in original space, or just preprocessed BOLD. The BOLD time-series were resampled into several standard spaces, correspondingly generating the following spatially-normalized, preprocessed BOLD runs: MNI152NLin2009cAsym, MNI152NLin6Asym. First, a reference volume and its skull-stripped version were generated using a custom methodology of fMRIPrep. Automatic removal of motion artifacts using independent component analysis ICA-AROMA, Pruim et al., 2015 was performed on the preprocessed BOLD on MNI space time-series after removal of non-steady state volumes and spatial smoothing with an isotropic, Gaussian kernel of 6mm FWHM (full-width half-maximum). Corresponding “non-aggresively” denoised runs were produced after such smoothing. Additionally, the “aggressive” noise-regressors were collected and placed in the corresponding confounds file. Several confounding time-series were calculated based on the preprocessed BOLD: framewise displacement (FD), DVARS and three region-wise global signals. FD was computed using two formulations following Power (absolute sum of relative motions, Power et al. (2014)) and Jenkinson (relative root mean square displacement between affines, Jenkinson et al. (2002)). FD and DVARS are calculated for each functional run, both using their implementations in Nipype following the definitions by Power et al., 2014. The three global signals are extracted within the CSF, the WM, and the whole-brain masks. Additionally, a set of physiological regressors were extracted to allow for component-based noise correction *CompCor*, Behzadi et al., 2007. Principal components are estimated after high-pass filtering the preprocessed BOLD time-series (using a discrete cosine filter with 128s cut-off) for the two CompCor variants: temporal (tCompCor) and anatomical (aCompCor). tCompCor components are then calculated from the top 2% variable voxels within the brain mask. For aCompCor, three probabilistic masks (CSF, WM and combined CSF+WM) are generated in anatomical space. The implementation differs from that of Behzadi et al. in that instead of eroding the masks by 2 pixels on BOLD space, the aCompCor masks are subtracted a mask of pixels that likely contain a volume fraction of GM. This mask is obtained by dilating a GM mask extracted from the FreeSurfer’s aseg segmentation, and it ensures components are not extracted from voxels containing a minimal fraction of GM. Finally, these masks are resampled into BOLD space and binarized by thresholding at 0.99 (as in the original implementation). Components are also calculated separately within the WM and CSF masks. For each CompCor decomposition, the k components with the largest singular values are retained, such that the retained components’ time series are sufficient to explain 50 percent of variance across the nuisance mask (CSF, WM, combined, or temporal). The remaining components are dropped from consideration. The head-motion estimates calculated in the correction step were also placed within the corresponding confounds file. The confound time series derived from head motion estimates and global signals were expanded with the inclusion of temporal derivatives and quadratic terms for each Satterthwaite et al., 2013. Frames that exceeded a threshold of 0.5 mm FD or 1.5 standardised DVARS were annotated as motion outliers. All resamplings can be performed with a single interpolation step by composing all the pertinent transformations (i.e. head-motion transform matrices, susceptibility distortion correction when available, and co-registrations to anatomical and output spaces). Gridded (volumetric) resamplings were performed using antsApplyTransforms (ANTs), configured with Lanczos interpolation to minimize the smoothing effects of other kernels Lanczos, 1964. Non-gridded (surface) resamplings were performed using mri_vol2surf (FreeSurfer).

Many internal operations of fMRIPrep use Nilearn 0.6.2 Abraham et al., 2014, RRID:SCR_001362, mostly within the functional processing workflow. For more details of the pipeline, see the section corresponding to workflows in fMRIPrep’s documentation.

The above boilerplate text was automatically generated by fMRIPrep with the express intention that users should copy and paste this text into their manuscripts unchanged. It is released under the CC0 license.

6.4References

7OHBM 2022 abstract - submitted text

The preliminary results will be presented at OHBM 2022 as a poster. Please find poster number WTh570.

Find the presenter at the virtual poster session on Wednesday, June 8, 2022, 8:30 PM GMT + 1.

At Glasgow, please contact the presenter on Twitter to schedule a time to chat, or come to see the presenter on Wednesday, June 22, 2022, 12:45 PM at the poster hall.

poster

7.1Impact of confound removal strategies on functional connectivity generated from fMRIPrep outputs

H-T Wang[1], S L Meisler[2][3], H Shamarke, F Paugam[1][4], N Gensollen[5], B Thirion[5], C Markiewicz[6], P Bellec[1][7]

7.1.1Introduction

Selecting an optimal denoising strategy is a key issue when processing fMRI data. The popular software fMRIPrep Esteban et al., 2018 aims to standardize fMRI preprocessing, but users are still offered a wide range of confound regressors to choose from to denoise data. Without a good understanding of the literature or the fMRIPrep documentation, users can select suboptimal strategies. Current work aims to provide a useful reference for fMRIPrep users by systematically evaluating the impact of different confound regression strategies, and by contrasting the results with past literature based on alternative preprocessing software.

7.1.2Methods

We selected dataset ds000228 Richardson et al., 2018 on OpenNeuro, which we preprocessed with fMRIPrep LTS20.2.1 using option --use-aroma. Time series were extracted using the Schaefer 7 network atlas with 400 ROIs Schaefer et al., 2017. We applied the denoising strategies listed in the table below using fMRIPrep-generated confounds. Subjects with less than 80% of remaining volumes after scrubbing with a 0.5 mm threshold were excluded from all analysis. We also calculated the connectome from high-pass filtered time series as a comparison baseline.

strategyimagehigh_passmotionwm_csfglobal_signalscrubfd_threshcompcorn_compcorica_aromademean
baselinedesc-preproc_boldTrueN/AN/AN/AN/AN/AN/AN/AN/ATrue
simpledesc-preproc_boldTruefullbasicN/AN/AN/AN/AN/AN/ATrue
simple+gsrdesc-preproc_boldTruefullbasicbasicN/AN/AN/AN/AN/ATrue
scrubbingdesc-preproc_boldTruefullfullN/A50.5N/AN/AN/ATrue
scrubbing+gsrdesc-preproc_boldTruefullfullbasic50.5N/AN/AN/ATrue
compcordesc-preproc_boldTruefullN/AN/AN/AN/Aanat_combinedallN/ATrue
compcor6desc-preproc_boldTruefullN/AN/AN/AN/Aanat_combined6N/ATrue
aromadesc-smoothAROMAnonaggr_boldTrueN/AbasicN/AN/AN/AN/AN/AfullTrue
aroma+gsrdesc-smoothAROMAnonaggr_boldTrueN/AbasicbasicN/AN/AN/AN/AfullTrue

We used three metrics Ciric et al., 2017, Parkes et al., 2018 to evaluate the denoising results:

  1. Quality control / functional connectivity (QCFC Power et al., 2015): partial correlation between motion and connectivity with age and sex as covariates. We control for multiple comparisons with false positive rate correction.

  2. Distance-dependent effects of motion on connectivity Power et al., 2012: correlation between node-wise Euclidean distance and QC-FC.

  3. Network modularity Satterthwaite et al., 2012: graph community detection based on Louvain method, implemented in the Brain Connectome Toolbox.

7.1.3Results

7.1.3.1QC-FC and distance-dependent effect

No denoise strategy removed the correlation with motion captured by mean framewise displacement. aroma, compcor6, and simple reduced the correlation between connectivity edges and mean framewise displacement. scrubbing and scrubbing+gsr performed the best, as seen in previous work Power et al., 2015. compcor, which applies compcor-based regressors covering 50% of the variance, performs worse than the connectome baseline created with high-pass filtered time series. Surprisingly, all strategies with global signal regression underperform, contradicting the existing literature Ciric et al., 2017 Parkes et al., 2018.

7.1.3.2Distance-dependent effects of motion on connectivity

Consistent with the literature, aroma reduces the distance dependency of motion on connectivity.

7.1.3.3Network modularity

All strategies increased the overall network modularity compared to the baseline, with scrubbing based methods performing the best out of all. GSR-based strategies improved the network modularity compared to their conunterparts. The correlation between modularity quality and motion for each denoising approach shows that compcor-based and ICA-AROMA strategies are the best at eliminating correlations between motion and modularity.

7.1.4Conclusions

We replicated previous findings demonstrating the usefulness of standard denoising strategies (compcor, aroma, etc.). However, results involving global signal regression methods systematically contradict the literatureCiric et al., 2017 Parkes et al., 2018. This evaluation is implemented in a fully reproducible jupyter book framework, and it can be applied to evaluate denoising strategies for future fMRIPrep release. This software may also be useful for researchers to select the most suitable strategy and produce denoising benchmarks for their own dataset.

8Software implementation

For the full report, please see our preprint.

The information here aims to provide context to understand the documentations.

8.1Denoising workflow

The denoising workflow is implemented through nilearn. fig-fmriprep-nilearn-denoise presents the graphic summary of the workflow.

The Python-based workflow describes the basic procedure to generate functional connectomes from fMRIPrep outputs with a Nilearn data loading routine (e.g., NiftiMapsMasker or NiftiLabelsMasker), fMRIPrep confounds output retrieval function (e.g., load_confounds_strategy), and connectome generation routine (ConnectivityMeasure). Path to the preprocessed image data is passed to load_confounds_strategy and the function fetches the associated confounds from the .tsv file. The path of an atlas and the path of the preprocessed image file is then passed to the masker, along with the confounds, for time series extraction. The time series are then passed to ConnectivityMeasure for generating connectomes.

8.2Benchmark workflow

fig-denoise-benchmark-workflow presents the graphic summary of the benchmark workflow.

The denoising benchmark workflow expands on the workflow in fig-fmriprep-nilearn-denoise (represented by the Nilearn logo in this figure). We retrieved the datasets from OpenNeuro through DataLad and all steps indicated with the arrows are implemented with bash scripts written for the SLURM scheduler. Atlases were either retrieved from the TemplateFlow archive or reformatted to fit the TemplateFlow format. The extracted time series, denoising metrics, and all metadata for generating the report are available on Zenodo.

9Quick start

I (HTW) aim to provide some documentation for people who would like to rebuild all material of this book. For those who wish to generate all the time series, connectomes, and metric, please see Setup to run all customised scripts for preprocessing.

If you only want to rebuild the content of the book for the review process, or help me edit the manuscript, the instruction is as followed:

10Setup to run all customised scripts for preprocessing

The dataset and features requires HPC for the computation. Here are scripts we used to process the current project presented in the study. If not otherwise specified, the scripts were executive from the root of the project. If you are a member of the SIMEXP team, you can run the script as is; otherwise, please install other software and modify the path and variables accordingly.

We recommend to understand how to create and detach a tmux session before further reading.

In most of the slurm scripts, we send SLURM job status email to slack bot, please see here to set-up your slack bot email address. Please modify the scripts accordingly.

10.1Atlases

After cloning and setting up the project, we will need to get all the atlases for the analysis.

10.1.1Option 1: Download the atlas

You can download the customised tempalteflow directory, untar into data/. For the purpose of separating processing and book building, we will copy the templates to inputs/custome_templateflow

make data
cp -r data/fmriprep-denoise-benchmark/custome_templateflow inputs/custome_templateflow

10.1.2Option 2: Generate the atlas

Alternatively, you can re-run the scripts. We will need to download from templateflow and nilearn, and organise them in templateflow standard. This is a shared step for those who wants to rerun the whole workflow. I separate the template fetching from template flowset up as the computing node of the HCP doesn’t have network access.

  1. Download all the atlases (require network access):

    python scripts/fetch_templates.py
  2. Process the atlases and organise the files to templateflow standard. One of the measure require the node-wise distances, so we need to separate the Difumo atlas dimensions into parcels.

    bash scripts/setup_templateflow.sh

    If this step is computational too intensive, you can run it in on a computing node.

  3. Calculate the centroid of the customised DiFuMo

    python scripts/calculate_centroids.py

11Preprocessing fMRI data

The steps here are required if you want to fully rerun the analysis. If you are not a SIMEXP member on Aliance Canada, you will need to clone fMRIPrep-slurm and modify paths accordingly. All scripts are under scripts/preprocessing. If you wish to use the produced metrics to generate the book, feel free to skip this section.

  1. preprocessing/get_fmriprep.sh:

    We need to build the fMRIPrep singularity containers version 20.2.1 (fMRIPrep-slurm default) and 20.2.5. For details, please find the documentation on niprep. Alternatively, you can execute the script prepared get_fmriprep.sh. We recommend building the container in a tmux session. This will take a few hours! Run it before going home.

  2. preprocessing/get_datasets.sh:

    While running Step 1 you can get the OpenNeuro datasets in parallel. Run the script in a tmux session. This will take an hour or so.

  3. preprocessing/create_fmriprep_slurm_scripts.sh:

    After the containers are built and datasets downloaded, we generate the slurm files using fMRIPrep-slurm.

    Remember to check your script to make sure everything runs correctly. This will take the time of a coffee break.

  4. Run fmriprep

    fMRIPrep-slurm will give you the exact commands you need to run. It should be something looking like this:

    find /scratch/${USER}/ds000228/UNIXTIME/.slurm/smriprep_sub-*.sh -type f | while read file; do sbatch "$file"; done

    This process will take up a day.

12Generate the reports

12.1timeseries and metrics

  1. Create a symbolic link of the metrics to the inputs/ directory

    mkdir $SCRATCH/fmriprep-denoise-benchmark/ 
  2. generate_timeseries_slurm_scripts.py:

    Extract time series with different atlases. The scripts generates slurm to extract time series with different atlases. Here’s the docs.

    usage: generate_timeseries_slurm_scripts.py [-h]
                                                [--slurm-account SLURM_ACCOUNT]
                                                scratch_path fmriprep_output
                                                participants_tsv virtualenv
    
    create timeseries extraction scripts
    
    positional arguments:
    scratch_path          Path to scratch space.
    fmriprep_output       Path to fMRIPrep output directory.
    participants_tsv      Path to participants.tsv in the original BIDS dataset.
    virtualenv            Path to virtual environment of this project.
    
    optional arguments:
    -h, --help            show this help message and exit
    --slurm-account SLURM_ACCOUNT
                            SLURM account for job submission (default: rrg-pbellec)

    We created two separate scripts for discrete and probability atlas due to different memory requirement. You will find the output under:

    /scratch/${USER}/fmriprep-denoise-benchmark/giga_timeseries/{DATASET_NAME}/{FMRIPREP_VERSION}/{UNIXTIME}/.slurm

    Similar to fmriprep-slurm, it will give you the exact commands you need to run. It should be something looking like this:

    find /scratch/${USER}/ds000228/UNIXTIME/.slurm/smriprep_sub-*.sh -type f | while read file; do sbatch "$file"; done

    This process will take a few hours.

  3. slurm_meta_confounds.sh:

    Create files to determine which subject will enter the next stage for metric generation.

    sbatch slurm_metric/slurm_meta_confounds.sh
  4. generate_metrics/*/slurm/metrics*.sh:

    Calculate metrics on denoising quality per atlas. Use this line to submit all jobs at once.

    find scripts/generate_metrics/slurm/metrics*.sh -type f | while read file; do sbatch $file; done

    The extra scripts in generate_metrics/ are for running directly on a computing node.

12.2Build the book

To improve build time, we need to summarise the metrics further. If you generated the data from scratch, you will need to run the following command.

usage: summarise_metadata [-h]
                          [--fmriprep_version {fmriprep-20.2.1lts,fmriprep-20.2.5lts}]
                          [--dataset_name {ds000228,ds000030}]
                          [--qc {stringent,minimal,None}]
                          output_root

Summarise denoising metrics for visualization and save at the top level of the denoise metric outputs directory.

positional arguments:
  output_root           Output root path data.

optional arguments:
  -h, --help            show this help message and exit
  --fmriprep_version {fmriprep-20.2.1lts,fmriprep-20.2.5lts}
                        Path to a fmriprep dataset.
  --dataset_name {ds000228,ds000030}
                        Dataset name.
  --qc {stringent,minimal,None}
                        Automatic motion QC thresholds.

Now you can build the book:

make book

Acknowledgments

The project is funded by IVADO PRF3, CCNA and J Courtois, the neuromind collaboration. HTW and NC funded by Institut de valorisation des données (IVADO) postdoctoral research funding. SLM was funded by the National Institute on Deafness and Other Communication Disorders (NIDCD; Grant 5T32DC000038). CJM funded by NIMH 5R24MH117179. FP funded by Courtois Neuromod. PB funded by Fonds de Recherche du Québec - Santé (FRQ-S).

Footnotes
  1. Centre de recherche de l’institut Universitaire de gériatrie de Montréal (CRIUGM), Montréal, Québec, Canada

  2. Psychology Department, Université de Montréal, Montréal, Québec, Canada

  3. Harvard University, MA, USA

  4. Massachusetts Institute of Technology, MA, USA

  5. Computer Science and Operations Research Department, Université de Montréal, Montréal, Québec, Canada

  6. Inria, CEA, Université Paris-Saclay, Paris, France

  7. Department of Psychology, Stanford University, Stanford, United States

References
  1. Esteban, O., Markiewicz, C., Blair, R. W., Moodie, C., Isik, A. I., Erramuzpe Aliaga, A., Kent, J., Goncalves, M., DuPre, E., Snyder, M., Oya, H., Ghosh, S., Wright, J., Durnez, J., Poldrack, R., & Gorgolewski, K. J. (2018). fMRIPrep: a robust preprocessing pipeline for functional MRI. Nature Methods. 10.1038/s41592-018-0235-4
  2. Abraham, A., Pedregosa, F., Eickenberg, M., Gervais, P., Mueller, A., Kossaifi, J., Gramfort, A., Thirion, B., & Varoquaux, G. (2014). Machine learning for neuroimaging with scikit-learn. Frontiers in Neuroinformatics, 8. 10.3389/fninf.2014.00014
  3. Karakuzu, A., DuPre, E., Tetrel, L., Bermudez, P., Boudreau, M., Chin, M., Poline, J.-B., Das, S., Bellec, P., & Stikov, N. (2022). NeuroLibre : A preprint server for full-fledged reproducible neuroscience. OSF Preprints. 10.31219/osf.io/h89js
  4. Richardson, H., Lisandrelli, G., Riobueno-Naylor, A., & Saxe, R. (2019). “MRI data of 3-12 year old children and adults during viewing of a short animated film.” OpenNeuro. 10.18112/openneuro.ds000228.v1.1.0
  5. Bilder, R., Poldrack, R., Cannon, T., London, E., Freimer, N., Congdon, E., Karlsgodt, K., & Sabb, F. (2020). “UCLA Consortium for Neuropsychiatric Phenomics LA5c Study.” OpenNeuro. 10.18112/openneuro.ds000030.v1.0.0
  6. Parkes, L., Fulcher, B., Yücel, M., & Fornito, A. (2018). An evaluation of the efficacy, reliability, and sensitivity of motion correction strategies for resting-state functional MRI. NeuroImage, 171, 415–436. 10.1016/j.neuroimage.2017.12.073
  7. Pruim, R. H. R., Mennes, M., van Rooij, D., Llera, A., Buitelaar, J. K., & Beckmann, C. F. (2015). ICA-AROMA: A robust ICA-based strategy for removing motion artifacts from fMRI data. NeuroImage, 112(Supplement C), 267–277. 10.1016/j.neuroimage.2015.02.064
  8. Power, J. D., Schlaggar, B. L., & Petersen, S. E. (2015). Recent progress and outstanding issues in motion correction in resting state fMRI. NeuroImage, 105, 536–551. 10.1016/j.neuroimage.2014.10.044
  9. Power, J. D., Barnes, K. A., Snyder, A. Z., Schlaggar, B. L., & Petersen, S. E. (2012). Spurious but systematic correlations in functional connectivity MRI networks arise from subject motion. NeuroImage, 59(3), 2142–2154. 10.1016/j.neuroimage.2011.10.018
  10. Satterthwaite, T. D., Wolf, D. H., Loughead, J., Ruparel, K., Elliott, M. A., Hakon, H., Gur, R. C., & Gur, R. E. (2012). Impact of In-Scanner Head Motion on Multiple Measures of Functional Connectivity: Relevance for Studies of Neurodevelopment in Youth. NeuroImage, 60(1), 623–632. 10.1016/j.neuroimage.2011.12.063
  11. Rubinov, M., & Sporns, O. (2010). Complex network measures of brain connectivity: Uses and interpretations. NeuroImage, 52(3), 1059–1069. 10.1016/j.neuroimage.2009.10.003
  12. Gordon, E. M., Laumann, T. O., Adeyemo, B., Huckins, J. F., Kelley, W. M., & Petersen, S. E. (2014). Generation and Evaluation of a Cortical Area Parcellation from Resting-State Correlations. Cerebral Cortex, 26(1), 288–303. 10.1093/cercor/bhu239
  13. Schaefer, A., Kong, R., Gordon, E. M., Laumann, T. O., Zuo, X.-N., Holmes, A. J., Eickhoff, S. B., & Yeo, B. T. T. (2017). Local-Global Parcellation of the Human Cerebral Cortex from Intrinsic Functional Connectivity MRI. Cerebral Cortex, 1–20. 10.1093/cercor/bhx179
  14. Urchs, S., Armoza, J., Moreau, C., Benhajali, Y., St-Aubin, J., Orban, P., & Bellec, P. (2019). MIST: A multi-resolution parcellation of functional brain networks. MNI Open Research, 1(3), 3.
  15. Dadi, K., Varoquaux, G., Machlouzarides-Shalit, A., Gorgolewski, K. J., Wassermann, D., Thirion, B., & Mensch, A. (2020). Fine-grain atlases of functional modes for fMRI analysis. NeuroImage, 221, 117126. 10.1016/j.neuroimage.2020.117126