NeoRS: A Neonatal Resting State fMRI Data Preprocessing Pipeline

Resting state functional MRI (rsfMRI) has been shown to be a promising tool to study intrinsic brain functional connectivity and assess its integrity in cerebral development. In neonates, where functional MRI is limited to very few paradigms, rsfMRI was shown to be a relevant tool to explore regional interactions of brain networks. However, to identify the resting state networks, data needs to be carefully processed to reduce artifacts compromising the interpretation of results. Because of the non-collaborative nature of the neonates, the differences in brain size and the reversed contrast compared to adults due to myelination, neonates can’t be processed with the existing adult pipelines, as they are not adapted. Therefore, we developed NeoRS, a rsfMRI pipeline for neonates. The pipeline relies on popular neuroimaging tools (FSL, AFNI, and SPM) and is optimized for the neonatal brain. The main processing steps include image registration to an atlas, skull stripping, tissue segmentation, slice timing and head motion correction and regression of confounds which compromise functional data interpretation. To address the specificity of neonatal brain imaging, particular attention was given to registration including neonatal atlas type and parameters, such as brain size variations, and contrast differences compared to adults. Furthermore, head motion was scrutinized, and motion management optimized, as it is a major issue when processing neonatal rsfMRI data. The pipeline includes quality control using visual assessment checkpoints. To assess the effectiveness of NeoRS processing steps we used the neonatal data from the Baby Connectome Project dataset including a total of 10 neonates. NeoRS was designed to work on both multi-band and single-band acquisitions and is applicable on smaller datasets. NeoRS also includes popular functional connectivity analysis features such as seed-to-seed or seed-to-voxel correlations. Language, default mode, dorsal attention, visual, ventral attention, motor and fronto-parietal networks were evaluated. Topology found the different analyzed networks were in agreement with previously published studies in the neonate. NeoRS is coded in Matlab and allows parallel computing to reduce computational times; it is open-source and available on GitHub (https://github.com/venguix/NeoRS). NeoRS allows robust image processing of the neonatal rsfMRI data that can be readily customized to different datasets.

Introduction

The analysis of resting-state functional connectivity (RS-FC) constitutes a promising tool as it provides complementary information to structural imaging related to brain physiology. Indeed, since its discovery in 1995 (Biswal et al., 1995) resting state functional MRI (rsfMRI) studies have provided new insights in the understanding of brain architecture and cerebral development (Power et al., 2010; Smyser and Neil, 2015; Gao et al., 2017; Grayson and Fair, 2017; Keunen et al., 2017; Zhang et al., 2019). Smyser et al. (2013) demonstrated the feasibility of using rsfMRI to explore the alterations in resting state networks (RSN) associated with preterm birth and white matter injury. Alterations of the default mode and ventral attention networks at birth, are associated with behavioral inhibition at age of two years, (Sylvester et al., 2018) which suggests early alterations of the RSN present a correlation with clinical manifestations, and opens the opportunity of early diagnostics and treatment. Additionally, neonatal RSN are consistently identifiable and present with high similarities to older populations (Fransson et al., 2007, 2009; Gao et al., 2009). RS-FC is based on low frequency regional fluctuations (<0.1 Hz) in the Blood-Oxygen-Level-Dependent (BOLD) (Ogawa et al., 1990, 1993) signal while the participant is not performing any task, a useful feature when evaluating neonates (Smyser and Neil, 2015). RSN signal is very stable across subjects (Lee et al., 2013), but vulnerable to several artifacts such as head-motion (Maknojia et al., 2019), susceptibility distortions and or white matter (WM) and cerebrospinal fluid (CSF) signals (Jo et al., 2013; Power et al., 2014). Robust rsfMRI data processing is key to reduce the nuisance effects of the non-neural signals in the data to identify reliable resting state activity (Lund et al., 2006; Giove et al., 2009). Its clinical potential and implementation present several methodological challenges that need to be addressed before considering its use to develop a new generation of biomarkers. For this reason, straightforward to use and open-source tools for the neonatal rsfMRI data processing need to be readily available. Tools for mature brains already exist to process rsfMRI data, but analyzing the neonatal brain presents challenges that need to be addressed with new approaches (Smyser and Neil, 2015). There are several straightforward rsfMRI data processing pipelines developed for adults such as Conn toolbox (Whitfield-Gabrieli and Nieto-Castanon, 2012),fmriprep (Esteban et al., 2019), the Human Connectome Pipeline (HCP) (Glasser et al., 2013), the Resting-State Analysis Toolkit (REST) (Rubinov and Sporns, 2010), or the Connectome Computation System (CCS) (Xu et al., 2015), however, those are not adapted to the newborn brain which presents additional challenges, such as different contrast due to myelination (Enguix et al., 2018). T2-weighted images are usually needed for tissue segmentation in place of T1-weighted images. Further, varying brain sizes between subjects (Smyser and Neil, 2015) makes adult skull stripping less robust on the neonatal brain. Additionally, different age specific atlases and tissue probability maps are required for accurate segmentations, common space normalization and seed-based analysis.

To the best of our knowledge the only existing open-access pipeline to process neonatal rsfMRI data is the one developed by the developing Human Connectome Project (dHCP) (Fitzgibbon et al., 2020). While this pipeline has proven to provide excellent results with the dHCP data, its implementation on smaller or clinical datasets remains challenging, as it requires large datasets for independent component (IC) denoising. Furthermore, the dHCP pipeline can be difficult to set up for cohorts acquired at other centers, because the pipeline was developed/optimized from the dHCP database specifically. For example, the dHCP denoising step is based on spatial independent component analysis (sICA), which separates independent correlating signals that can be classified as neural or non-neural signal. This denoising technique has been shown to provide superior results in adults and infants when the dimensionality is accurately set (Griffanti et al., 2017; Alfaro-Almagro et al., 2018). However, the identified signals need to be classified as neural signal or structured noise, which in most cases is performed manually and a difficult process to automate. To overcome this limitation, the dHCP pipeline uses a machine learning approach (ICA-based Xnoiseifier) (Salimi-Khorshidi et al., 2014) to classify the independent components as neural signals or noise. The machine learning algorithm requires a minimum of 35 manually labeled subjects to be trained, which is not always possible in smaller cohorts and requires specialists to manually classify the independent components (Fitzgibbon et al., 2020).

To overcome the aforementioned challenges, we developed NeoRS, with the goal of creating a robust open-source pipeline containing the necessary tools to preprocess rsfMRI data. The main advantages of NeoRS are it has been developed specifically for neonates, simple to implement and flexible to process different datasets. Additionally, it can process single subject data, utilizes parallelizable environment and includes visual quality control checkpoints at each step.

The data processing steps include T2-weighted image alignment to a common space, slice timing correction and segmentation, and rsfMRI procedures are slice timing, distortion correction using reversed phase encoding polarity acquisitions, alignment in a common space, motion correction, removal of nuisance confounds and noise compromising functional data interpretation. Further, simple resting state functional connectivity cross-correlations based on seed-to-voxel and seed-to-seed approaches are incorporated.

Materials and Methods

Data

NeoRS has been evaluated on neonates (7 ± 1.4 weeks old) from the Baby Connectome Project (BCP) (Howell et al., 2019) dataset. For this study only participants scanned at 9 weeks old or less and contained T2-weighted images and rsfMRI were used (N = 10). Participants were naturally sleeping and were scanned on a 3.0 T MRI Prisma from Siemens using a 32-channel head coil. This study included a T2-weighted structural image (TE = 564 ms, TR = 3,200 ms, matrix = 320 mm × 320 mm, FOV = 256 mm × 256 mm, resolution = 0.8 × 0.8 × 0.8 mm, flip angle = variable, in-plane acceleration factor = 2, acquisition time = 5 min 57 s), two gradient-echo (GRE) echo-planar imaging (EPI) blip-up/blip-down (TE = 37 ms, TR = 800 ms, matrix = 104 mm × 91 mm, FOV = 208 mm × 182 mm, resolution = 2 mm × 2 mm × 2 mm, flip angle = 52°, multiband acceleration factor = 8, acquisition time = 5 min 47 s, 420 volumes), and two spin-echo (SE) EPI blip-up/blip-down for distortion correction purpose (TE = 66 ms, TR = 8,000 ms, matrix = 104 mm × 91 mm, FOV = 208 mm × 182 mm, resolution = 2 mm × 2 mm × 2 mm, flip angle = 52°, multiband acceleration factor = 1, acquisition time = 33 s, 3 volumes).

Data Structure

To facilitate collaborations, NeoRS uses the Brain Imaging Data Structure (BIDS) format as described in https://bids.neuroimaging.io/. See Figure 1 for an example of data naming and organization for NeoRS.

FIGURE 1

Figure 1. Example of data naming and organization for NeoRS.

Pipeline Overview

NeoRS is a neonatal rsfMRI data processing pipeline developed on Matlab and calls for commands developed on well-known open-source neuroimaging tools, such as FSL 6.0.3.1 (Smith et al., 2004; Woolrich et al., 2009; Jenkinson et al., 2012), AFNI 20.2.10 (Cox, 1996) and SPM 12.¹ It runs on both MacOS and Linux operating systems. NeoRS has been tested on a MacBook pro 2015 with operating system High Sierra and on a Linux computer using Ubuntu 18.04.5 by running the pipeline from the beginning to end on different subjects on both computers. NeoRS is built to accommodate MRI data acquired with different manufacturers. The pipeline was tested using the aforementioned BCP data, as well as the not publicly available data from CHU Sainte-Justine acquired on a GE 3T MR750, but this manuscript focuses only on BCP results. Furthermore, the pipeline has single subject capabilities. NeoRS has been developed to accommodate a parallelizable environment, allowing several subjects to be simultaneously processed depending on the number of selected cores. To investigate, two subjects were processed on an early 2015 MacBook pro with 2.7 GHz Intel Core i5 processor and 8 GB 1,867 MHz DDR3 memory by using a single core vs. 2 parallel cores and found a reduction of computing time of 1.8 times when using the 2 parallel cores. This function is optional and requires the Matlab parallel toolbox.

See NeoRS workflow in Figure 2.

FIGURE 2

Figure 2. NeoRS workflow for neonatal resting state functional connectivity processing and denoising.

To ensure images align to the orientation of the standard template a reorientation to standard is performed in both structural and functional data prior to other data processing procedures by employing fslreorient2std from FSL. Furthermore, to guarantee the accurate performance of NeoRS, output files for each processing step are saved in a folder called Output_files. Processing steps are dependent from previous outputs and should be inspected carefully. In case of fail, the user can parameterize the specific function, as specified later on. Various operations are not mandatory, such as slice-timing correction or distortion correction, and can be manually turned-off by setting the function parameter to 0 in the main file. See Figure 3 for an example of inputs configuration.

FIGURE 3

Figure 3. Example of NeoRS input parameters.

Data Processing

Structural

T2-Weighted Image Registration

NeoRS uses the term age stereotaxic space (Smyser et al., 2010) from Washington University – School of Medicine. The template is available in Talairach space (Talairach, 1988) 1 and 3 mm isotropic resolutions. Image registration in NeoRS is performed using FSL flirt and is implemented in a single step with 12 degrees of freedom and not applying the resampling blur when down sampling. These parameters can be modified by the user in the function anat2std.m. High resolution T2-weighted images are registered to a 1 and 3 mm isotropic template.

Skull Stripping

Skull stripping plays an important role in image processing, as it is mandatory for different processing functionalities, such as tissue segmentation, and requires special attention to avoid further complications in the process. NeoRS skull stripping step utilizes the FSL (Jenkinson et al., 2012) function bet2 and has been optimized for term neonatal brains. Skull stripping is performed after image registration to obtain consistent results independent of brain size. Furthermore, visual quality control is available in a file containing the brain with the skull and the overlay of the contour of the intracranial cavity. If the user is not satisfied with the results, modify the fractional intensity threshold, “-f”, and vertical gradient in fractional intensity threshold, “-g”, to properly adjust skull stripping in the Matlab function skull_stripping.m.

Segmentation

Extracted T2-weighted intracranial content is then segmented to create different tissue probability maps corresponding to each brain structure. Tissue segmentation is crucial in image processing as the outputs will be used for regression of confounds. For brain segmentation NeoRS applies Morphologically Adaptive Neonatal Tissue Segmentation: Mantis (Beare et al., 2016). Mantis is an SPM based toolbox and allows T2-weighted image segmentation based on template adaptation via topological filters and morphological segmentation tools, resulting in eight different tissue probability maps. The segmentation process is fully integrated in the NeoRS pipeline and has been tested on three different datasets (BCP, and CHU Sainte-Justine). After segmentation, the eight tissue probability maps are automatically combined, thresholded and binarized to create three different binary masks needed to run downstream processing. The masks correspond to WM, gray matter (GM), and CSF. The masks are resampled to 3-mm isotropic to match functional image space resolution.

Functional

Slice Timing Correction

A common acquisition technique for rsfMRI is the single-shot Gradient-Echo (GRE) Echo-Planar Imaging (EPI). In this acquisition sequence, slices are acquired at varying intervals, which need to be addressed. The NeoRS function for slice timing correction is FSL slicetimer and can automatically read the slice order from the.json file, if available. If the.json file is not available, the user can manually define the slice order or use one of the predefined options from fsl (i.e., interleaved ascending) in the configuration file.

Functional Cross Realignment

To correct for head movement, it is necessary to obtain a motion estimation based on 6 movement parameters (three rotation and three translation parameters). This is done by rigid-body registration (six degrees of freedom) between the different volumes with respect to a reference, in NeoRS, the reference is the first volume from the rsfMRI, but can be easily altered by the user if desired. This NeoRS function is performed using FSL mcflirt (Jenkinson et al., 2002) and works the same as for adults, however we set the smoothness level to 0, as smoothing occurs later in the pipeline, and used sinc interpolation. Parameters for cross realignment can be customized in the Matlab function cross_realign2.m. For quality control purposes, NeoRS creates a.png file where the total rotations, translations and framewise displacement (FD) for each volume can be evaluated. Framewise displacement is calculated as previously described by Power et al. (2012). To take into account head size differences, the calculations were done using a 35 mm radius sphere instead of 50 mm which approximately corresponds to the mean distance from the cerebral cortex to the center of the head in neonates. After motion correction, the motion parameters are saved in a text file that will be further used for denoising purposes.

Functional Best Resting State Section Selection

NeoRS incorporates the possibility to analyze sub sections of long time series (i.e., 20 min). This tool is deactivated by default but can be activated by setting options.best_volumes = 1 in the configuration file. The sectional analysis tool automatically identifies a section of the time-series (i.e., 5 min) with the lowest average FD, which is recommended to use on very long acquisitions that present a higher average FD than the threshold. The length of the section can be modified by the user, but it is recommended the duration remain above 5 min (Sylvester et al., 2018). To note, the average FD threshold has a default of 0.25 mm, but can also be tailored as needed by altering options.FDaverage in the configuration file. NeoRS chooses only the best sections of the time-series, which reduces computational times drastically.

Functional Distortion Correction

The EPI sequence is considerably sensitive to off-resonance fields due to susceptibility variations of participants. To address these distortions a typical approach is to use two SE-EPI reversed polarity acquisitions (reversed phase encoding direction) to estimate the distortion field. This field is implemented to correct for distortions in the original GRE EPI images. Directly using two reversed polarity GRE EPI to estimate the distortion field instead of two reversed polarity SE EPI is possible, but not recommended because GRE-EPI sequences are hampered by signal dropouts caused by intravoxel dephasing. Providing reversed phase encoding polarity SE EPI images and activating the distortion correction option (options.fmap = 1) in NeoRS, allows users to estimate these distortions utilizing FSL topup (Andersson et al., 2003; Smith et al., 2004; Graham et al., 2017). Such as specified in the topup documentation, a text file containing the encoding directions and total readout time needs to be included in the fmap folder to perform distortion corrections. If the.json file is found in the fmap folder, the text file will be automatically created by NeoRS. Distortion estimates are rectified with FSL applytopup for each EPI volume by applying the output from topup.

Functional Image Registration

Functional images are registered to the same stereotaxic space template from the Washington University – School of Medicine as the T2-weighted image registration procedure. Initially, a two-step registration was performed. First, a rigid body registration between the mean rsfMRI volume and T2-weighted images was calculated, followed by affine transformation between T2-weighted and the template. Finally, the output affine transformation matrices from each process were applied to the rsfMRI images to align to the 3 mm isotropic template. Additionally, a single step registration approach was investigated, where the rsfMRI images were aligned directly to the 3 mm isotropic template using 12 degrees of freedom. This registration approach was comparable to the 2-step registration process and was chosen as it was accurate and faster. Down-sampling blur, by default is set to off in NeoRS, however, if needed, the parameters can be customized in the Matlab function epi2std2.m.

Denoising

Motion Censoring

Before the regression of the confounding signals, volumes with excessive motion are removed based on the framewise displacement metric described by Power et al. (2012). NeoRS performs linear detrending and computes framewise displacement after functional cross realignment based on 6 motion parameters in radians (3 rotation parameters + 3 translation parameters). This step also automatically removes the first five volumes.

$$\mathrm{FD}=|\mathrm{rot}\mathrm{\_}\mathrm{x}|+|\mathrm{rot}\mathrm{\_}\mathrm{y}|+|\mathrm{rot}\mathrm{\_}\mathrm{z}|+|\mathrm{trans}\mathrm{\_}\mathrm{x}|+|\mathrm{trans}\mathrm{\_}\mathrm{y}|+|\mathrm{trans}\mathrm{\_}\mathrm{z}|$$

Where rot_x/rot_y/rot_z are rotations converted from radians to mm and trans_x/trans_y/trans_z are translations in mm. Once the FD is computed, a text file containing the information of volumes exceeding the FD threshold plus the first 5 frames, is created and head motion plots are saved and can be reviewed. The FD threshold is automatically set to FD < 0.25 mm (Smyser et al., 2016), so volumes with FD higher than or equal to 0.25 mm are excluded. Excluded volumes are set to zero value and no interpolation is applied to avoid artificial correlations.

High Motion Subjects

High motion acquisitions are source of artifacts and may confound neural correlations with non-neural signals. Apart from single frame motion censoring based on FD NeoRS evaluates the average FD for every BOLD run. By default, NeoRS was set to discard acquisitions with an average FD higher than 0.25 mm. The average FD threshold can be altered in the configuration file by defining options.FDaverage.

Regression of Confounds

Variables identified as potential confounders of the estimated BOLD signal are merged in a single file. To avoid frequency mismatch in the regression process, the file is used to compute linear regression in a single step with frequency filtering. This process is performed utilizing AFNI 3dTproject.

Motion Parameters

Head motion is considered as a rigid body moving in a 3D space with 6 degrees of freedom. In cartesian coordinates we can describe it with 3 translations x- (left/right), y- (anterior/posterior), and z-axes (inferior/superior), and 3 rotations around the x-axis (pitch), y-axis (yaw), and z-axis (roll). To address the residual motion related signal variance after a suboptimal rigid body registration, NeoRS uses a linear regression strategy based on the 6 aforementioned estimated motion parameters. Those parameters are considered as nuisance effects of the signal and are then removed. NeoRS allows various options including: 6 motion parameters, 12 (including temporal derivatives) (Power et al., 2012) or 24 (including temporal derivatives and their squares) (Friston et al., 1996) which can be defined in the configuration file parameter options.motion.

White Matter, Cerebrospinal Fluid, and Global Signals

White matter and cerebrospinal fluid signals are highly confounding and need to be removed from the rsfMRI (Smyser et al., 2016). Signals for regression of confounds are extracted from the WM and CSF masks generated previously in the pipeline from the segmentation. Masks are created in a conservative way by selecting voxels from the tissue probability masks with higher probability than 0.5. The voxels of the white matter are eroded by one voxel to ensure the mask doesn’t include any GM. Two files containing the average signal of the WM and CSF masks are created. Finally, global signal is approximated by averaging the signal in a GM mask (Vos de Wael et al., 2017) and used by default in NeoRS, as it improves data quality by reducing motion artifacts (cardiac, respiration, head motion) (Smyser et al., 2016).

Frequency Filter

Temporal frequencies outside the frequency range of [0.01–0.1] Hz are removed from the BOLD signal to correct for slow frequency drifts, reduce motion artifacts and other physiological noises while preserving the frequencies of resting state networks (Power et al., 2014). The use of a low-pass filter could drastically reduce the degrees of freedom of the time series in acquisitions with very short TR. For those cases, it is recommend to set the value of options.BPF = [HPH, LPF], to [0.01, 999] in the configuration file.

Smoothing

Functional Smoothing

Functional smoothing is the last processing procedure. After denoising the rsfMRI signal is convolved with a gaussian kernel. This reduces the effect of misregistration between functional regions and slightly increases the signal to noise ratio. Gaussian smoothing is performed implementing fslmaths from FSL. The size of the gaussian kernel is customizable in the NeoRS pipeline by modifying options.fwhm in the configuration file, which is 6 mm by default.

Data Analysis—Functional Connectivity

Prior to further data analysis, like ROI to ROI (region of interest) correlations, all the processed BOLD runs are merged together into a single 4D-file. NeoRS offers basic single subject data analysis, including seed based and seed to seed correlations, so the user can further assess data has been correctly processed.

Seed-Based Correlations

Seed-based functional connectivity identifies correlation between a defined ROI, also called a seed, and the rest of the brain. This metric facilitates the observation of simultaneously activated regions with the pre-defined ROI. The NeoRS pipeline provides 31 template seeds representing some of the most common resting state networks including: language, default mode, dorsal attention, visual, ventral attention, motor and fronto-parietal networks. An excel file (Perceptron_ROI_list.xlsx) can be found in the documentation with all the information related to seed positioning (Smyser et al., 2016).

Seed-To-Seed Measurements

Seed-to-seed data analysis provides measurements of functional connectivity between all the different pairs of seeds demonstrating a more global perspective about networks compared to seed-based functional connectivity. NeoRS performs Pearson correlation between the different ROIs to create a correlation matrix.

Results

Image Registration

Image registration results of the T2-weighted and BOLD images to the template for a representative subject are demonstrated in Figure 4, an example of a user checkpoint. Yellow lines represent the segmented cerebrospinal fluid, and is added as an overlay to the T2-weighted images, BOLD and template. After visual inspection, a correct alignment within the template for both registrations was observed for each test participant.

FIGURE 4

Figure 4. T2-weighted and BOLD image registration to stereotaxic space. Gray scale images represent the T2-weighted, average BOLD and template images; the yellow lines correspond to the cerebrospinal fluid contours obtained from the T2-weighted image segmentations.

When comparing single step registration versus a 2-step registration approach for rsfMRI there were no discernible differences between both registrations and final functional connectivity results presented the same correlation strengths and topology. The main difference between the two approaches was computational times, which were higher for the 2-step registration method.

Skull Stripping

Figure 5 illustrates the default skull stripping segmentations versus NeoRS adapted parameters for neonates. With the default settings we observed skull stripping was failing for some of the subjects with different brain sizes, in contrast, when using NeoRS parameters, skull stripping remained robust for all the processed subjects.

FIGURE 5

Figure 5. Skull stripping parameters comparing default bet2 settings in neonates and NeoRS settings optimized for neonates.

Segmentation and Mask Creation

Figure 6 displays the 1 mm isotropic binary masks created by NeoRS from Mantis tissue probability maps. The output contains three different binary files corresponding to white matter, SCF and GM. Figure 7 demonstrates the 3 mm isotropic masks for the regression of confounds process.

FIGURE 6

Figure 6. One millimeter isotropic masks created from the tissue probability maps obtained with Mantis. White matter (red), CSF (yellow), gray matter (blue).

FIGURE 7

Figure 7. White matter, cerebrospinal fluid, and gray matter masks for regression of confounds.

Functional Distortion Correction

Functional susceptibility-induced magnetic field inhomogeneity correction for a representative subject is shown in Figure 8. GRE EPI images, independent of brain size are distorted in the phase encoding direction whether they are acquired AP or PA but present those distortions in both areas of the brain. After susceptibility-induced magnetic field inhomogeneity distortion correction, the two acquisitions (AP and PA) present a similar morphology and a more accurate brain shape with respect to the undistorted T2-weighted image.

FIGURE 8

Figure 8. Susceptibility-induced magnetic field inhomogeneity causing geometric distortions along the phase encoding direction. (A) Original T2-weighted image without distortion, shown as reference; (B) original GRE-EPI acquired in anterior-posterior phase encoding direction (AP); (C) original GRE-EPI acquired in posterior-anterior phase encoding direction (PA); (D) corrected GRE-EPI AP; (E) corrected GRE-EPI PA.

Head Motion

After functional cross-realignment, an output graph is provided by NeoRS containing information concerning rotations and translations applied to cross-realign for each volume of the rsfMRI, as well as the computed framewise displacement (Figure 9).

FIGURE 9

Figure 9. Example of head motion plots from a single subject. Plots are generated for each bold run and contain three different graphs per run: estimated rotation in radians; estimated translation in millimeters; Framewise displacement in millimeters.

Figure 10 is an example of a single subject with 2 different rsfMRI acquisitions with different amounts of motion. In the seed-based functional connectivity results of the motor network, the correlation differences between an acquisition with an average framewise displacement higher than 0.25 mm (run 1) and an acquisition with an average framewise displacement lower than 0.25 mm (run 2). High motion acquisition presented increased amounts of noise and the network topology was difficult to identify.

FIGURE 10

Figure 10. Two acquisitions of a high motion subject, run 1 excluded for having an average FD ≥ 0.25 mm, run 2 kept with an average FD < 0.25 mm.

Resting State Networks—Seed-Based Correlations

Figure 11 illustrates, seven of the most common resting state networks after NeoRS processing employing an seed-based correlations (SBC) approach.

FIGURE 11

Figure 11. Example resting state networks obtained by seed-based functional connectivity after image processing with NeoRS.

Figure 12 is a single subject example of the 31 seeds included in NeoRS.

FIGURE 12

Figure 12. Representative subject resting state network example seed-to-seed functional connectivity correlations.

Discussion

NeoRS is an open-source image processing pipeline dedicated to neonatal rsfMRI. It includes seed-based and seed-to-seed 1st level analysis. NeoRS has neonatal brain templates for term in 1mm and 3mm Talairach space, as well as a set of 31 seed regions defining seven common resting state networks. NeoRS relies on the open-source neuroimaging pipelines: SPM, FSL and AFNI and encompasses robust methods to segment, register and denoise neonatal rsfMRI data. Each of the various processing steps were evaluated separately. Output was carefully inspected to ensure the best quality products, by optimizing skull stripping, image registration, head motion and denoising related results. Additionally, functional connectivity of the motor, visual, default mode, language, dorsal attention, ventral attention and fronto-parietal networks using seed-based functional connectivity analysis were demonstrated.

Image registration to the atlas was meticulously inspected for every subject and no significant misalignments were observed for T2-weighted or rsfMRI images for affine registration. We compared single step to two-step registration for accuracy and computational times. The single step registration was chosen, as this process required less computational time and demonstrated no substantial variation when compared to the two-step counterpart.

In contrast with some of the most common adult pipelines, image registration is implemented prior to skull stripping, as alignment quality results were identical. Further, performing image registration preceding skull stripping produced vastly robust skull stripping results across varying brain sizes. This was not the case employing traditional skull stripping before image registration. Skull stripping of the T2-weighted images was found to be a crucial step when working with neonates because poor stripping lead to misclassification of segmented data and ultimately unreliable representation of RSN because of misregistrations. Additionally, image registration prior to skull stripping facilitated brain extraction without any user intervention. Furthermore, if the brain was previously aligned to an atlas, bet2 was implemented on the subject specific atlas aligned data instead of applying an atlas brain mask to avoid subtle geometric inaccuracies. This procedure is critical and needs to be properly assessed. For this reason, an output for the skull stripping is provided which contains the image of the non-skull stripped mask with an overlay of the contour of the skull stripped brain.

After skull stripping, tissue segmentation is performed using Mantis, without the need of any further intervention. Mantis is integrated into the NeoRS pipeline as it provided robust results for disparate brain sizes using only T2-weighted images. Contrary to the adult brain that uses T1-weighted images for brain segmentation, it is fundamental in neonates to provide the pipeline with T2-weighted images as the water/cholesterol ratio is reversed with respect to adults due to lack of myelination. Neonatal T2-weighted images present a better contrast between brain structures (McArdle et al., 1987). While Mantis needs to be installed to use the NeoRS pipeline, no additional setup steps are required as it is fully assimilated in NeoRS. The results showed the binary masks created from Mantis tissue probability maps are perfectly aligned with their corresponding structures for various subjects without any manual intervention.

It is well known the GRE-EPI sequence for rsfMRI is prone to susceptibility-induced magnetic field inhomogeneity (Czervionke et al., 1988). The artifact primarily appears close to the extrema portions of the brain in the phase encoding direction (Andersson et al., 2018) and needs to be properly corrected. While several methods have been successfully used in these settings (Jezzard and Balaban, 1995; Cusack et al., 2003), NeoRS operates the standard topup/applytopup method. Two reversed phase encoding direction images are involved to correct for the deformations (Holland et al., 2010). The method is simple to implement in the acquisition protocol, provides high quality results and acquisition sequences require very short duration times.

Slice timing correction remains a controversial step when a very short repetition time (TR) is deployed (Parker and Razlighi, 2019). However, as it has been shown it can significantly improve z-scores (Parker and Razlighi, 2019) and as the aim was to make NeoRS work with the maximum number of datasets, it is included as an option. Slice timing correction can be deactivated for multi-band sequences with very low TR, as in this kind of low TR sequence, all slices in each volume are acquired very closer together (Glasser et al., 2013).

After data preprocessing, confounding signals and motion effects are removed. First, the framewise displacement threshold is defined as 0.25 mm, as performed by Smyser et al. on their neonatal study and shown to provide accurate results (Sylvester et al., 2018). Volumes with FD higher than 0.25 mm were removed from the time-series as described per Power et al. (2012). Motion censoring was performed prior to filtering to prevent spikes from passing through band-pass filtering, as this could introduce artifacts such as Gibb’s ringing and or skew correlation coefficients. Furthermore, extremely high motion acquisitions shouldn’t be taken into account as they can lead to inflated results. To do so, different metrics can be adopted, such as maximum framewise displacement, minimum number of low motion volumes or average framewise displacement. In NeoRS, acquisitions with average FD higher than 0.25 mm were introducing augmented correlations related to motion and improper denoising. Those acquisitions were completely removed. Therefore, the average framewise displacement was employed as the metric for exclusion.

To correct for nuisance variables NeoRS implements a traditional denoising strategy that performs global signal, white matter and cerebrospinal fluid signal regression, motion parameter regression and a band-pass filter. This simple approach provides robust results (Smyser et al., 2016), and doesn’t require manual intervention or large datasets for denoising purposes. In contrast to the aforementioned independent components denoising techniques, (Griffanti et al., 2017; Alfaro-Almagro et al., 2018) such as the one employed in the dHCP pipeline.

Finally, NeoRS incorporates seed-based functional connectivity analysis tools to assist the user in assessing initial results. Seed-based results across subjects showed patterns very similar to those observed in the literature for all the analyzed resting state networks (Fransson et al., 2009; Damaraju et al., 2010; Figure 11). Further data analysis can be carried out on the fully processed data, final_BOLD.nii, if desired.

Limitations and Future Directions

NeoRS was fully vetted on the BCP, and CHU (results not shown) cohorts for neonates less than 9 weeks old. Expanding NeoRS to a larger cohort in both size and neonatal age variation, i.e., preterm, would only further demonstrate its novelty and application, but the current number of subjects is a limitation. Obviously, this expansion would also necessitate additional age specific atlases. Additionally, the pipeline requires a usable T2-weighted structural image and would benefit from the adaptability of the option of using a T1-weighted image, especially for younger neonates older than 9 weeks where tissue boundaries may begin to vary. NeoRS has been developed based on a traditional denoising strategy which includes band-pass filtering. This is a limitation on data acquired with very low TR, such as the dHCP data, as the degrees of freedom might be highly reduced. In the Fourier domain, the maximum sample rate corresponds to the frequency of Nyquist, f_max = 1/2TR, and the frequency spacing Δf = 1/t_max, where t_max is the total scan duration. If there are two data sets with the same total acquisition time, but different TR, the one with the lower TR will lose more degrees of freedom when filtering (Bright et al., 2017). Another limitation of the pipeline is brains with significant malformations or injury. While a high range of brain sizes is accepted in the pipeline, anomalies such as neonatal hydrocephalus may require special attention, as part of the automated process could fail, such as cortical extraction. To overcome this limitation the user should perform rigorous individual image quality control of those brains and may adjust parameters, such as the head radius. Furthermore, it is known that a fundamental factor for measuring interindividual differences is reliability (Zuo et al., 2014) which is highly related to reproducibility. For this reason, performing a test-retest reliability approach should have been considered (Noble et al., 2019). However, measuring the test-retest reliability is a challenge in this population. As mentioned by Zuo et al. (2019), at least 20–30 min of data are needed to perform functional connectivity reliability measurements which can be challenging. Moreover, performing multiple sessions is also a challenge because of the fast brain development. Additionally, neonates also present higher levels of motion, which compromises reliability (Zuo and Xing, 2014). To overcome this issue, we think that we should focus on acquiring high quality data by using real-time head motion monitoring tools such as FIRMM (Dosenbach et al., 2017). Finally, we are currently working on the inclusion of second-level analysis to allow the users to make group-level inference about networks. The possibility to perform Independent Component Analysis will also be available in the next release. Further functionalities such as the characterization of the networks and lifespan developmental trajectories of cortical thickness and surface area, as described by Xu et al. (2015) might be of particular interest as it would allow the study of the maturation process of the neonatal brain and the resting state networks. The characterization of the resting state networks over time is especially important on newborns as a deeper knowledge on networks development might help detecting possible outliers in a given population, bringing neonatal rsfMRI one step closer of being a brain integrity biomarker.

Conclusion

NeoRS² is an open-source, straightforward to use rsfMRI data processing pipeline for the neonatal brain which relies on the open-source neuroimaging pipelines FSL, AFNI and SPM. NeoRS works with neuroimaging nifti format, BIDS folder structures and has been developed to work with different MRI vendors and diverse acquisition parameters with minimal user implication. After image processing with NeoRS, we observed resting state networks were in agreement with previously published studies at term age. Each processing step is easy to inspect to ensure consistent results through quality control checkpoint figures.

An open-source, rudimentary to use pipeline for neonatal resting state image processing will allow the community to process their data immediately after scanning sessions implementing a simple computational infrastructure. The democratization of rsfMRI processing will allow a higher number of centers to collaborate and process their datasets and optimistically bring the clinical biomarker application one step closer.

Data Availability Statement

The original contributions presented in this study are included in the article/supplementary material, further inquiries can be directed to the corresponding authors.

Author Contributions

VE developed NeoRS, processed and analyzed the data, and wrote the manuscript. JK provided guidance and helped comparing results with a non-open-source pipeline. DL contributed the writing and reviewing the first draft and provided guidance. JC-A provided guidance and reviewed the manuscript. GL supervised the project, contributed writing and reviewing all drafts, and provided guidance. All authors contributed to the manuscript revision, read, and approved the submitted version.

Funding

This work was funded by the Canada Research Chair in Quantitative Magnetic Resonance Imaging (950-230815), the Canadian Institute of Health Research (CIHR FDN-143263), the Canada Foundation for Innovation (32454 and 34824), the Fonds de Recherche du Québec - Santé (28826r), the Natural Sciences and Engineering Research Council of Canada (RGPIN-2019-07244), the Canada First Research Excellence Fund (IVADO and TransMedTech), the Courtois NeuroMod project and the Quebec BioImaging Network (5886 and 35450), and INSPIRED (Spinal Research, United Kingdom; Wings for Life, Austria; Craig H. Neilsen Foundation, United States).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We thank the Washington University – School of Medicine for sharing their templates for this project, and Jed Elison for letting us use the neonatal data from the Baby Connectome Project. We also thank the Québec Bio-Imaging Network for supporting VE with studentship and funding for travels.

Footnotes

1. ^ https://www.fil.ion.ucl.ac.uk/spm/
2. ^ https://github.com/venguix/NeoRS

References

Alfaro-Almagro, F., Jenkinson, M., Bangerter, N. K., Andersson, J. L. R., Griffanti, L., Douaud, G., et al. (2018). Image processing and quality control for the first 10,000 brain imaging datasets from UK Biobank. Neuroimage 166, 400–424. doi: 10.1016/j.neuroimage.2017.10.034

PubMed Abstract | CrossRef Full Text | Google Scholar

Andersson, J. L. R., Graham, M. S., Drobnjak, I., Zhang, H., and Campbell, J. (2018). Susceptibility-induced distortion that varies due to motion: correction in diffusion MR without acquiring additional data. Neuroimage 171, 277–295. doi: 10.1016/j.neuroimage.2017.12.040

PubMed Abstract | CrossRef Full Text | Google Scholar

Andersson, J. L., Skare, S., and Ashburner, J. (2003). How to correct susceptibility distortions in spin-echo echo-planar images: application to diffusion tensor imaging. Neuroimage 20, 870–888. doi: 10.1016/S1053-8119(03)00336-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Beare, R. J., Chen, J., Kelly, C. E., Alexopoulos, D., Smyser, C. D., Rogers, C. E., et al. (2016). Neonatal brain tissue classification with morphological adaptation and unified segmentation. Front. Neuroinform. 10:12. doi: 10.3389/fninf.2016.00012

PubMed Abstract | CrossRef Full Text | Google Scholar

Biswal, B., Yetkin, F. Z., Haughton, V. M., and Hyde, J. S. (1995). Functional connectivity in the motor cortex of resting human brain using echo-planar MRI. Magn. Reson. Med. 34, 537–541. doi: 10.1002/mrm.1910340409

PubMed Abstract | CrossRef Full Text | Google Scholar

Bright, M. G., Tench, C. R., and Murphy, K. (2017). Potential pitfalls when denoising resting state fMRI data using nuisance regression. Neuroimage 154, 159–168. doi: 10.1016/j.neuroimage.2016.12.027

PubMed Abstract | CrossRef Full Text | Google Scholar

Cox, R. W. (1996). AFNI: software for analysis and visualization of functional magnetic resonance neuroimages. Comput. Biomed. Res. 29, 162–173.

Google Scholar

Cusack, R., Brett, M., and Osswald, K. (2003). An evaluation of the use of magnetic field maps to undistort echo-planar images. Neuroimage 18, 127–142. doi: 10.1006/nimg.2002.1281

PubMed Abstract | CrossRef Full Text | Google Scholar

Czervionke, L. F., Daniels, D. L., Wehrli, F. W., Mark, L. P., Hendrix, L. E., Strandt, J. A., et al. (1988). Magnetic susceptibility artifacts in gradient-recalled echo MR imaging. AJNR Am. J. Neuroradiol. 9, 1149–1155.

PubMed Abstract | Google Scholar

Damaraju, E., Phillips, J. R., Lowe, J. R., Ohls, R., Calhoun, V. D., and Caprihan, A. (2010). Resting-state functional connectivity differences in premature children. Front. Syst. Neurosci. 4:23. doi: 10.3389/fnsys.2010.00023

PubMed Abstract | CrossRef Full Text | Google Scholar

Dosenbach, N. U. F., Koller, J. M., Earl, E. A., Miranda-Dominguez, O., Klein, R. L., Van, A. N., et al. (2017). Real-time motion analytics during brain MRI improve data quality and reduce costs. Neuroimage 161, 80–93. doi: 10.1016/j.neuroimage.2017.08.025

PubMed Abstract | CrossRef Full Text | Google Scholar

Enguix, V., Ding, Y., and Lodygensky, G. A. (2018). Recent advances in preclinical and clinical multimodal MR in the newborn brain. J. Magn. Reson. 292, 149–154. doi: 10.1016/j.jmr.2018.04.017

PubMed Abstract | CrossRef Full Text | Google Scholar

Esteban, O., Markiewicz, C. J., Blair, R. W., Moodie, C. A., Isik, A. I., Erramuzpe, A., et al. (2019). fMRIPrep: a robust preprocessing pipeline for functional MRI. Nat. Methods 16, 111–116. doi: 10.1038/s41592-018-0235-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Fitzgibbon, S. P., Harrison, S. J., Jenkinson, M., Baxter, L., Robinson, E. C., Bastiani, M., et al. (2020). The developing Human Connectome Project (dHCP) automated resting-state functional processing framework for newborn infants. Neuroimage 223:117303. doi: 10.1016/j.neuroimage.2020.117303

PubMed Abstract | CrossRef Full Text | Google Scholar

Fransson, P., Skiold, B., Engstrom, M., Hallberg, B., Mosskin, M., Aden, U., et al. (2009). Spontaneous brain activity in the newborn brain during natural sleep–an fMRI study in infants born at full term. Pediatr. Res. 66, 301–305. doi: 10.1203/PDR.0b013e3181b1bd84

PubMed Abstract | CrossRef Full Text | Google Scholar

Fransson, P., Skiold, B., Horsch, S., Nordell, A., Blennow, M., Lagercrantz, H., et al. (2007). Resting-state networks in the infant brain. Proc. Natl. Acad. Sci. U.S.A. 104, 15531–15536.

Google Scholar

Friston, K. J., Williams, S., Howard, R., Frackowiak, R. S., and Turner, R. (1996). Movement-related effects in fMRI time-series. Magn. Reson. Med. 35, 346–355. doi: 10.1002/mrm.1910350312

PubMed Abstract | CrossRef Full Text | Google Scholar

Gao, W., Lin, W., Grewen, K., and Gilmore, J. H. (2017). Functional connectivity of the infant human brain: plastic and modifiable. Neuroscientist 23, 169–184. doi: 10.1177/1073858416635986

PubMed Abstract | CrossRef Full Text | Google Scholar

Gao, W., Zhu, H., Giovanello, K. S., Smith, J. K., Shen, D., Gilmore, J. H., et al. (2009). Evidence on the emergence of the brain’s default network from 2-week-old to 2-year-old healthy pediatric subjects. Proc. Natl. Acad. Sci. U.S.A. 106, 6790–6795. doi: 10.1073/pnas.0811221106

PubMed Abstract | CrossRef Full Text | Google Scholar

Giove, F., Gili, T., Iacovella, V., Macaluso, E., and Maraviglia, B. (2009). Images-based suppression of unwanted global signals in resting-state functional connectivity studies. Magn. Reson. Imaging 27, 1058–1064. doi: 10.1016/j.mri.2009.06.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Glasser, M. F., Sotiropoulos, S. N., Wilson, J. A., Coalson, T. S., Fischl, B., Andersson, J. L., et al. (2013). The minimal preprocessing pipelines for the human connectome project. Neuroimage 80, 105–124. doi: 10.1016/j.neuroimage.2013.04.127

PubMed Abstract | CrossRef Full Text | Google Scholar

Graham, M. S., Drobnjak, I., Jenkinson, M., and Zhang, H. (2017). Quantitative assessment of the susceptibility artefact and its interaction with motion in diffusion MRI. PLoS One 12:e0185647. doi: 10.1371/journal.pone.0185647

PubMed Abstract | CrossRef Full Text | Google Scholar

Grayson, D. S., and Fair, D. A. (2017). Development of large-scale functional networks from birth to adulthood: a guide to the neuroimaging literature. Neuroimage 160, 15–31. doi: 10.1016/j.neuroimage.2017.01.079

PubMed Abstract | CrossRef Full Text | Google Scholar

Griffanti, L., Douaud, G., Bijsterbosch, J., Evangelisti, S., Alfaro-Almagro, F., Glasser, M. F., et al. (2017). Hand classification of fMRI ICA noise components. Neuroimage 154, 188–205. doi: 10.1016/j.neuroimage.2016.12.036

PubMed Abstract | CrossRef Full Text | Google Scholar

Holland, D., Kuperman, J. M., and Dale, A. M. (2010). Efficient correction of inhomogeneous static magnetic field-induced distortion in Echo Planar Imaging. Neuroimage 50, 175–183. doi: 10.1016/j.neuroimage.2009.11.044

PubMed Abstract | CrossRef Full Text | Google Scholar

Howell, B. R., Styner, M. A., Gao, W., Yap, P. T., Wang, L., Baluyot, K., et al. (2019). The UNC/UMN Baby Connectome Project (BCP): An overview of the study design and protocol development. Neuroimage 185, 891–905. doi: 10.1016/j.neuroimage.2018.03.049

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenkinson, M., Bannister, P., Brady, M., and Smith, S. (2002). Improved optimization for the robust and accurate linear registration and motion correction of brain images. Neuroimage 17, 825–841. doi: 10.1016/s1053-8119(02)91132-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Jenkinson, M., Beckmann, C. F., Behrens, T. E., Woolrich, M. W., and Smith, S. M. (2012). Fsl. Neuroimage 62, 782–790.

Google Scholar

Jezzard, P., and Balaban, R. S. (1995). Correction for geometric distortion in echo planar images from B0 field variations. Magn. Reson. Med. 34, 65–73. doi: 10.1002/mrm.1910340111

PubMed Abstract | CrossRef Full Text | Google Scholar

Jo, H. J., Gotts, S. J., Reynolds, R. C., Bandettini, P. A., Martin, A., Cox, R. W., et al. (2013). Effective preprocessing procedures virtually eliminate distance-dependent motion artifacts in resting state FMRI. J. Appl. Math 2013, doi: 10.1155/2013/935154

PubMed Abstract | CrossRef Full Text | Google Scholar

Keunen, K., Counsell, S. J., and Benders, M. (2017). The emergence of functional architecture during early brain development. Neuroimage 160, 2–14. doi: 10.1016/j.neuroimage.2017.01.047

PubMed Abstract | CrossRef Full Text | Google Scholar

Lee, M. H., Smyser, C. D., and Shimony, J. S. (2013). Resting-state fMRI: a review of methods and clinical applications. AJNR Am. J. Neuroradiol. 34, 1866–1872. doi: 10.3174/ajnr.A3263

PubMed Abstract | CrossRef Full Text | Google Scholar

Lund, T. E., Madsen, K. H., Sidaros, K., Luo, W. L., and Nichols, T. E. (2006). Non-white noise in fMRI: does modelling have an impact? Neuroimage 29, 54–66. doi: 10.1016/j.neuroimage.2005.07.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Maknojia, S., Churchill, N. W., Schweizer, T. A., and Graham, S. J. (2019). Resting state fMRI: going through the motions. Front. Neurosci. 13:825. doi: 10.3389/fnins.2019.00825

PubMed Abstract | CrossRef Full Text | Google Scholar

McArdle, C. B., Richardson, C. J., Nicholas, D. A., Mirfakhraee, M., Hayden, C. K., and Amparo, E. G. (1987). Developmental features of the neonatal brain: MR imaging. Part I. Gray-white matter differentiation and myelination. Radiology 162(1 Pt .1), 223–229. doi: 10.1148/radiology.162.1.3786767

PubMed Abstract | CrossRef Full Text | Google Scholar

Noble, S., Scheinost, D., and Constable, R. T. (2019). A decade of test-retest reliability of functional connectivity: a systematic review and meta-analysis. Neuroimage 203:116157. doi: 10.1016/j.neuroimage.2019.116157

PubMed Abstract | CrossRef Full Text | Google Scholar

Ogawa, S., Lee, T. M., Kay, A. R., and Tank, D. W. (1990). Brain magnetic resonance imaging with contrast dependent on blood oxygenation. Proc. Natl. Acad. Sci. U.S.A. 87, 9868–9872.

Google Scholar

Ogawa, S., Menon, R. S., Tank, D. W., Kim, S. G., Merkle, H., Ellermann, J. M., et al. (1993). Functional brain mapping by blood oxygenation level-dependent contrast magnetic resonance imaging. A comparison of signal characteristics with a biophysical model. Biophys. J. 64, 803–812. doi: 10.1016/S0006-3495(93)81441-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Parker, D. B., and Razlighi, Q. R. (2019). The benefit of slice timing correction in common fMRI preprocessing pipelines. Front. Neurosci. 13:821. doi: 10.3389/fnins.2019.00821

PubMed Abstract | CrossRef Full Text | Google Scholar

Power, J. D., Barnes, K. A., Snyder, A. Z., Schlaggar, B. L., and Petersen, S. E. (2012). Spurious but systematic correlations in functional connectivity MRI networks arise from subject motion. Neuroimage 59, 2142–2154. doi: 10.1016/j.neuroimage.2011.10.018

PubMed Abstract | CrossRef Full Text | Google Scholar

Power, J. D., Fair, D. A., Schlaggar, B. L., and Petersen, S. E. (2010). The development of human functional brain networks. Neuron 67, 735–748.

Google Scholar

Power, J. D., Mitra, A., Laumann, T. O., Snyder, A. Z., Schlaggar, B. L., and Petersen, S. E. (2014). Methods to detect, characterize, and remove motion artifact in resting state fMRI. Neuroimage 84, 320–341. doi: 10.1016/j.neuroimage.2013.08.048

PubMed Abstract | CrossRef Full Text | Google Scholar

Rubinov, M., and Sporns, O. (2010). Complex network measures of brain connectivity: uses and interpretations. Neuroimage 52, 1059–1069. doi: 10.1016/j.neuroimage.2009.10.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Salimi-Khorshidi, G., Douaud, G., Beckmann, C. F., Glasser, M. F., Griffanti, L., and Smith, S. M. (2014). Automatic denoising of functional MRI data: combining independent component analysis and hierarchical fusion of classifiers. Neuroimage 90, 449–468. doi: 10.1016/j.neuroimage.2013.11.046

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, S. M., Jenkinson, M., Woolrich, M. W., Beckmann, C. F., Behrens, T. E., Johansen-Berg, H., et al. (2004). Advances in functional and structural MR image analysis and implementation as FSL. Neuroimage 23(Suppl. 1), S208–S219. doi: 10.1016/j.neuroimage.2004.07.051

PubMed Abstract | CrossRef Full Text | Google Scholar

Smyser, C. D., and Neil, J. J. (2015). Use of resting-state functional MRI to study brain development and injury in neonates. Semin. Perinatol. 39, 130–140. doi: 10.1053/j.semperi.2015.01.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Smyser, C. D., Inder, T. E., Shimony, J. S., Hill, J. E., Degnan, A. J., Snyder, A. Z., et al. (2010). Longitudinal analysis of neural network development in preterm infants. Cereb. Cortex 20, 2852–2862. doi: 10.1093/cercor/bhq035

PubMed Abstract | CrossRef Full Text | Google Scholar

Smyser, C. D., Snyder, A. Z., Shimony, J. S., Blazey, T. M., Inder, T. E., and Neil, J. J. (2013). Effects of white matter injury on resting state fMRI measures in prematurely born infants. PLoS One 8:e68098. doi: 10.1371/journal.pone.0068098

PubMed Abstract | CrossRef Full Text | Google Scholar

Smyser, C. D., Snyder, A. Z., Shimony, J. S., Mitra, A., Inder, T. E., and Neil, J. J. (2016). Resting-state network complexity and magnitude are reduced in prematurely born infants. Cereb. Cortex 26, 322–333. doi: 10.1093/cercor/bhu251

PubMed Abstract | CrossRef Full Text | Google Scholar

Sylvester, C. M., Smyser, C. D., Smyser, T., Kenley, J., Ackerman, J. J. Jr., Shimony, J. S., et al. (2018). Cortical functional connectivity evident after birth and behavioral inhibition at age 2. Am. J. Psychiatry 175, 180–187. doi: 10.1176/appi.ajp.2017.17010018

PubMed Abstract | CrossRef Full Text | Google Scholar

Talairach, J. (1988). Co-planar Stereotaxic Atlas of the Human Brain: 3-Dimensional Proportional System:An Approach to Cerebral Imaging. New York, NY: Thieme Medical Publishers, Inc.

Google Scholar

Vos de Wael, R., Hyder, F., and Thompson, G. J. (2017). Effects of tissue-specific functional magnetic resonance imaging signal regression on resting-state functional connectivity. Brain Connect. 7, 482–490. doi: 10.1089/brain.2016.0465

PubMed Abstract | CrossRef Full Text | Google Scholar

Whitfield-Gabrieli, S., and Nieto-Castanon, A. (2012). Conn: a functional connectivity toolbox for correlated and anticorrelated brain networks. Brain Connect. 2, 125–141. doi: 10.1089/brain.2012.0073

PubMed Abstract | CrossRef Full Text | Google Scholar

Woolrich, M. W., Jbabdi, S., Patenaude, B., Chappell, M., Makni, S., Behrens, T., et al. (2009). Bayesian analysis of neuroimaging data in FSL. Neuroimage 45(1 Suppl.), S173–S186. doi: 10.1016/j.neuroimage.2008.10.055

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, T., Yang, Z., Jiang, L., Xing, X.-X., and Zuo, X.-N. (2015). A connectome computation system for discovery science of brain. Sci. Bull. 60, 86–95.

Google Scholar

Zhang, H., Shen, D., and Lin, W. (2019). Resting-state functional MRI studies on infant brains: a decade of gap-filling efforts. Neuroimage 185, 664–684. doi: 10.1016/j.neuroimage.2018.07.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Zuo, X. N., and Xing, X. X. (2014). Test-retest reliabilities of resting-state FMRI measurements in human brain functional connectomics: a systems neuroscience perspective. Neurosci. Biobehav. Rev. 45, 100–118. doi: 10.1016/j.neubiorev.2014.05.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Zuo, X. N., Anderson, J. S., Bellec, P., Birn, R. M., Biswal, B. B., Blautzik, J., et al. (2014). An open science resource for establishing reliability and reproducibility in functional connectomics. Sci. Data 1:140049. doi: 10.1038/sdata.2014.49

PubMed Abstract | CrossRef Full Text | Google Scholar

Zuo, X. N., Xu, T., and Milham, M. P. (2019). Harnessing reliability for neuroscience research. Nat. Hum. Behav. 3, 768–771. doi: 10.1038/s41562-019-0655-x

PubMed Abstract | CrossRef Full Text | Google Scholar