Benchmarking confound regression strategies for the control of motion artifact in studies of functional connectivity
Rastko Ciric, Daniel H. Wolf, Jonathan D. Power, David R. Roalf, Graham Baum, Kosha Ruparel, Russell T. Shinohara, Mark A. Elliott, Simon B. Eickhoff, Christos Davatzikos, Ruben C. Gur, Raquel E. Gur, Danielle S. Bassett, Theodore D. Satterthwaite
BBenchmarking confound regression strategies for the control of motion artifact instudies of functional connectivity
Rastko Ciric a , Daniel H. Wolf a , Jonathan D. Power b , David R. Roalf a , Graham Baum a , Kosha Ruparel a , Russell T.Shinohara c , Mark A. Elliott d , Simon B. Eickhoff e,f , Christos Davatzikos d , Ruben C. Gur a , Raquel E. Gur a , Danielle S.Bassett g,h , Theodore D. Satterthwaite a, ∗ a Department of Psychiatry, Perelman School of Medicine, University of Pennsylvania, Philadelphia PA, USA b Sackler Institute for Developmental Psychobiology, Weill Medical College of Cornell University, New York, NY, USA c Department of Biostatistics and Epidemiology, Perelman School of Medicine, University of Pennsylvania, Philadelphia, PA, USA d Department of Radiology, Perelman School of Medicine, University of Pennsylvania, Philadelphia PA, USA e Department of Clinical Neuroscience and Medical Psychology, Heinrich-Heine University D¨usseldorf f Institute of Neuroscience and Medicine (INM-1), Research Center J¨ulich g Department of Bioengineering, University of Pennsylvania, Philadelphia PA, USA h Department of Electrical and Systems Engineering, University of Pennsylvania, Philadelphia PA, USA
Abstract
Since initial reports regarding the impact of motion artifact on measures of functional connectivity, there has been aproliferation of confound regression methods to limit its impact. However, recent techniques have not been systematicallyevaluated using consistent outcome measures. Here, we provide a systematic evaluation of 12 commonly used confoundregression methods in 193 young adults. Specifically, we compare methods according to three benchmarks, includingthe residual relationship between motion and connectivity, distance-dependent effects of motion on connectivity, andadditional degrees of freedom lost in confound regression. Our results delineate two clear trade-offs among methods. First,methods that include global signal regression minimize the relationship between connectivity and motion, but unmaskdistance-dependent artifact. In contrast, censoring methods mitigate both motion artifact and distance-dependence,but use additional degrees of freedom. Taken together, these results emphasize the heterogeneous efficacy of proposedmethods, and suggest that different confound regression strategies may be appropriate in the context of specific scientificgoals.
Keywords: fMRI, functional connectivity, artifact, confound, motion, noise
Introduction
Resting-state (intrinsic) functional connectivity (rsfc-MRI) has evolved to become one of the most commonbrain imaging phenotypes (Craddock et al., 2013; Fox andRaichle, 2007; Power et al., 2014b; Smith et al., 2013; VanDijk et al., 2010), and has been critical for understand-ing fundamental properties of brain organization (Damoi-seaux et al., 2006; Fox et al., 2005; Power et al., 2011;Yeo et al., 2011), brain development over the lifespan (Di-Martino et al., 2014; Dosenbach et al., 2011; Fair et al.,2008), and abnormalities associated with diverse clinicalconditions (Baker et al., 2014; Buckner et al., 2008; Fairet al., 2010). rsfc-MRI has numerous advantages, includ-ing ease of acquisition and suitability for a wide and ex-panding array of analysis techniques. However, despiteknowledge that in-scanner motion can influence measuresof activation from task-related fMRI (Friston et al., 1996), ∗ Corresponding author
Email address: [email protected] (Theodore D.Satterthwaite ) the impact of in-scanner motion on measures of functionalconnectivity was not explored for 16 years after its initialdiscovery. However, since the near-simultaneous publica-tion of three independent reports in early 2012 (Van Dijket al., 2012; Power et al., 2012; Satterthwaite et al., 2012),it has been increasingly recognized that motion can have alarge impact on rsfc-MRI measurements, and can system-atically bias inference. This bias is particularly problem-atic in developmental or clinical populations where mo-tion is correlated with the independent variable of inter-est (age, diagnosis) (Satterthwaite et al., 2012; Fair et al.,2012), and has resulted in the re-evaluation of numerouspublished findings.In response to this challenge, there has been a recentproliferation of methods aimed at mitigating the impactof motion on functional connectivity (Yan et al., 2013a;Power et al., 2015). These methods can be broadlygrouped into several categories. First, high-parameter con-found regression strategies use expansions of realignmentparameters or tissue-compartment signals, often includingderivative and quadratic regressors (Friston et al., 1996;Satterthwaite et al., 2013; Yan et al., 2013a). Second, prin-
Preprint submitted to arXiv August 15, 2016 a r X i v : . [ q - b i o . N C ] A ug ipal component analysis (PCA) based methods (Comp-Cor; Behzadi et al. (2007); Muschelli et al. (2014)) find theprimary directions of variation within high-noise areas de-fined by anatomy (e.g., aCompCor) or temporal variance(tCompCor). Third, whole-brain independent componentanalysis (ICA; Beckmann et al. (2005)) of single-subjecttime series has increasingly been used for de-noising, withnoise components selected either by a trained classifier(ICA-FIX; Griffanti et al. (2014); Salimi-Khorshidi et al.(2014)) or using a priori heuristics (ICA-AROMA; Pruimet al. (2015b,a)). Fourth, temporal censoring techniquesidentify and remove (or de-weight) specific volumes con-taminated by motion artifact, often followed by interpo-lation. These techniques include scrubbing (Power et al.,2012, 2014a, 2015), spike regression (Satterthwaite et al.,2013), and de-spiking (Jo et al., 2013; Patel et al., 2014).Censoring techniques have been reported to attenuate mo-tion artifact, but at the cost of a shorter time series andvariably reduced degrees of freedom. Fifth, one recent re-port emphasized the relative merits of spatially-tailoredconfound regression using local white matter signals (wm-Local; Jo et al. (2013)). Finally, the inclusion of globalsignal regression (GSR) (Macey et al., 2004) in confoundregression models remains a source of controversy (Foxet al., 2009; Murphy et al., 2009; Chai et al., 2012; Saadet al., 2012; Yan et al., 2013c). While several studies havesuggested its utility in de-noising (Fox et al., 2009; Poweret al., 2015; Satterthwaite et al., 2013; Yan et al., 2013a),other studies have emphasized the risk of removing a valu-able signal (Yang et al., 2014; Hahamy et al., 2014), po-tentially biasing group differences (Gotts et al., 2013; Saadet al., 2012), or exacerbating distance-dependent motionartifact. Distance-dependent artifact (Power et al., 2012;Satterthwaite et al., 2012) manifests as increased connec-tivity in short-range connections, and reduced connectiv-ity in long-range connections, which has the potential toimpact measures of network topology (Yan et al., 2013b).This recent proliferation of de-noising techniques hasprompted excitement but also sowed confusion. Unsurpris-ingly, new de-noising pipelines have often tended to em-phasize outcome measures that suggest their relative supe-riority. As a result, investigators often anecdotally reportsubstantial uncertainty regarding which pipeline should beused. Such uncertainty has been exacerbated by the lack ofcommon outcome measures used across studies, which hashampered direct comparison among pipelines. While onereview paper has summarized recent developments in thisrapidly-evolving sub-field (Power et al., 2015), systematicevaluation of de-noising pipelines according to a range ofbenchmarks remains lacking.Accordingly, in this report we compare a dozen of themost commonly used confound regression strategies in alarge ( N = 193) dataset of young adults. Pipelines eval-uated include standard techniques, high-parameter con-found regression, PCA-based techniques such as aComp-Cor and tCompcor, ICA-based approaches such as ICA-AROMA, spatially-tailored local white matter regression, and three different censoring techniques (spike regression,de-spiking, and scrubbing); GSR is included in manypipelines as well. Critically, we compare these pipelines ac-cording to three intuitive benchmarks, including the resid-ual relationship between functional connectivity and sub-ject motion, the degree of distance-dependent artifact, andthe loss of temporal degrees of freedom. As described be-low, results underscore the relative strengths and weak-nesses among these methods, and outline clear trade-offsamong commonly used confound regression approaches. Materials and methods
Participants and data acquisition
The task-free BOLD data used in this study ( N = 193)were selected from the Philadelphia NeurodevelopmentalCohort (PNC) (Satterthwaite et al., 2014, 2016) on the ba-sis of age, health, and data quality. In order to minimizethe potential impact of developmental effects, subjects se-lected for inclusion were aged at least 18 years at the dateof scan. All participants were free of medical problemsthat could impact brain function (Merikangas et al., 2010),lacked gross structural brain abnormalities (Gur et al.,2013), and were not taking psychotropic medication. Fur-thermore, subjects were excluded from the study if theyfailed to satisfy any of three criteria for functional imagequality: if subject movement (relative root mean squaredisplacement) averaged over all frames exceeded 0.2mm;if over 20 individual frames featured movement in excessof 0.25mm (Satterthwaite et al., 2012); or if brain cov-erage during acquisition was incomplete. As a result ofthese criteria, participants with gross in-scanner motionwere excluded, allowing us to evaluate the utility of con-found regression strategies for the mitigation of artifactdue to micro-movements.Structural and functional subject data were acquiredon a 3T Siemens Tim Trio scanner with a 32-channelhead coil (Erlangen, Germany), as previously described(Satterthwaite et al., 2014, 2016). High-resolution struc-tural images were acquired in order to facilitate align-ment of individual subject images into a common space.Structural images were acquired using a magnetization-prepared, rapid-acquisition gradient-echo (MPRAGE) T1-weighted sequence ( T R = 1810ms; T E = 3 . × T R = 3000ms; T E = 32ms;FoV = 192 × Structural image processing
A study-specific template was generated from a sam-ple of 120 PNC subjects balanced across sex, race, andage bins using the buildTemplateParallel procedure inANTs (Avants et al., 2011a). Study-specific tissue priorswere created using a multi-atlas segmentation procedure(Wang et al., 2014). Next, each subject’s high-resolutionstructural image was processed using the ANTs CorticalThickness Pipeline (Tustison et al., 2014). Following biasfield correction (Tustison et al., 2010), each structural im-age was diffeomorphically registered to the study-specificPNC template using the top-performing SyN deformation(Klein et al., 2009). Study-specific tissue priors were usedto guide brain extraction and segmentation of the subject’sstructural image (Avants et al., 2011b).
BOLD time series processing
Task-free functional images were processed using theXCP Engine (Ciric et al., In Preparation), which was con-figured to support the 12 pipelines evaluated in this study(see
Figure 1 ). Each pipeline was based on de-noisingstrategies previously described in the neuroimaging lit-erature. A number of preprocessing procedures were in-cluded across all de-noising pipelines. Common elementsof preprocessing included (1) correction for distortions in-duced by magnetic field inhomogeneities, (2) removal ofthe 4 initial volumes of each acquisition, (3) realignment ofall volumes to a selected reference volume using mcflirt (Jenkinson et al., 2002), (4) demeaning and removal ofany linear or quadratic trends, (5) co-registration of func-tional data to the high-resolution structural image usingboundary-based registration (Greve and Fischl, 2009), and(6) temporal filtering using a first-order Butterworth fil-ter with a passband between 0.01 and 0.08 Hz. Thesepreliminary processing stages were then followed by theconfound regression procedures described below. In or-der to prevent frequency-dependent mismatch during con-found regression (Hallquist et al., 2013), all regressors wereband-pass filtered to retain the same frequency range asthe data. As in our prior work (Satterthwaite et al., 2012,2013), the primary summary metric of subject motion usedwas the mean relative root-mean-square displacement cal-culated during time series realignment using mcflirt . Overview of confound regression strategies
The primary objective of the current study was to evalu-ate the performance of common de-noising strategies. Weselected 12 de-noising models, labelled – below, forevaluation (Figure 1). Models – used nuisance parame-ters derived from 6 movement estimates and 3 physiolog-ical time series, as well as their temporal derivatives andquadratic expansions. • Model 1. (2P) Used only the 2 physiological time se-ries: mean signal in WM and mean signal in CSF, andfunctioned as a base model for comparison to othermore complex confound regression models. • Model 2. (6P) Used only the 6 motion estimates de-rived from mcflirt realignment as explanatory vari-ables. • Model 3. (9P) Combined the 6 motion estimatesand 2 physiological time series with global signal re-gression. This model has been widely applied to func-tional connectivity studies (Fox et al., 2005, 2009). • Model 4. (24P) Expansion of model that in-cludes 6 motion parameters, 6 temporal derivatives,6 quadratic terms, and 6 quadratic expansions of thederivatives of motion estimates for a total 24 regres-sors (Friston et al., 1996). • Model 5. (36P) Similar expansion of model :9 regressors, their derivatives, quadratic terms, andsquares of derivatives (Satterthwaite et al., 2013).Models – further expanded upon this maximal 36Pstrategy by incorporating censoring approaches. • Model 6. (36P+despike) Included 36 regressors aswell as despiking (Cox, 1996). • Model 7. (36P+spkreg) Included 36 regressors aswell as spike regression, as in Satterthwaite et al.(2013). • Model 8. (36P+scrub) Included 36 regressors as wellas motion scrubbing, as in Power et al. (2014a).Models and adapted variants of the PCA-based CompCor approach. • Model 9. (aCompCor) Used 5 principal componentseach from the WM and CSF, in addition to motionestimates and their temporal derivatives (Muschelliet al., 2014). • Model 10. (tCompCor) Used 6 principal compo-nents from high-variance voxels (Behzadi et al., 2007).The final two models evaluated used a regressor de-rived from local WM signals, and subject-specific ICA de-noising, respectively.3
Model 11. (wmLocal) Used a voxelwise, localisedWM regressor in addition to motion estimates andtheir temporal derivatives and despiking (Jo et al.,2013). • Model 12. (ICA-AROMA) Used a recently devel-oped ICA-based procedure for removal of motion-related variance from BOLD data, together with meanWM and CSF regressors (Pruim et al., 2015a,b).We explicitly limited our scope to models that did notrequire training a classifier, and did not evaluate confoundregression strategies that require extensive parameter op-timization (Salimi-Khorshidi et al., 2014; Griffanti et al.,2014; Patel et al., 2014). Furthermore, in order to con-strain the parameter space, we did not examine unpub-lished combinations of de-noising approaches.
Confound regression using realignment parameters
Time series of six realignment parameters (three transla-tional and three rotational) for each subject were returnedby mcflirt as part of time series realignment (motion cor-rection). Additionally, the temporal derivative, quadraticterms, and quadratic of the temporal derivative of each ofthe realignment parameters were calculated, yielding 24 re-alignment regressors in total. The original six realignmentparameters were included in confound regression models and . Models and included 12 realignment regressors– the 6 realignment parameters and their temporal deriva-tives – while the full set of 24 expanded realignment regres-sors were included as part of confound regression models – . Global signal regression
The mean global signal was computed by averagingacross all voxelwise time series located within a subject-specific mask covering the entire brain. The global signalwas included in model , while the expanded models – included 4 global signal regressors: the global signal, itsderivative, its square, and the derivative of its square. Tissue class regressors
Mean white matter (WM) and cerebrospinal fluid (CSF)signals were computed by averaging within masks derivedfrom the segmentation of each subject’s structural image;these masks were eroded using AFNI’s 3dmask tool (Cox,1996) to prevent inclusion of gray matter signal via partial-volume effects. The WM mask was eroded at the 2-voxellevel, while the CSF mask was eroded at the 1-voxel level.More liberal erosion often resulted in empty masks in ourdata. Temporal derivatives, quadratic terms, and squaresof the derivative were computed as above. Two tissue classregressors (WM and CSF) were included in models and , whereas their expansions (8 regressors) were includedin models – . Local white matter regression
Model used a local WM regressor (Jo et al., 2013).This was computed in AFNI using 3dLocalstat (Cox,1996). Unlike the regressors described above, which werevoxel-invariant, the value of the local WM regressor wascomputed separately at each voxel. For each voxel, asphere of radius 45mm was first centered on that voxel;this sphere defined that voxel’s local neighborhood. Next,this spherical neighborhood was intersected with an erodedWM mask to produce a local WM mask, which includedonly the fraction of the WM that was also in the voxel’sneighborhood. The mean signal within this new local WMmask was then used to model the local WM signal at thevoxel (Jo et al., 2013). This process was repeated at everyvoxel in order to generate the local WM regressor. Thislocal WM regressor was included in model along withrealignment parameters and their derivatives (12 total);this model also included voxelwise de-spiking. CompCor
Principal component analysis (PCA) can be used tomodel noise in BOLD time series (Behzadi et al., 2007;Muschelli et al., 2014). Broadly, the use of PCA-derivedregressors to model noise is called component-based cor-rection (CompCor). Numerous variants of CompCor havebeen developed; here, our focus will be on anatomicalCompCor (aCompCor, model ) and temporal CompCor(tCompCor, model ). In aCompCor, a PCA is per-formed within an anatomically defined tissue class of in-terest. We extracted 5 components for WM and CSFeach, yielding 10 compcor components (Muschelli et al.,2014). As part of model , as in Muschelli et al. (2014),these 10 aCompCor components were combined with 12re-alignment parameters (raw and temporal derivative).In tCompCor, the temporal variance of the BOLD signalis first computed at each voxel. Subsequently, a mask isgenerated from high-variance voxels, and principal compo-nents are extracted from the time series at these voxels.In confound regression model , tCompcor was imple-mented using ANTs, with 6 tCompCor components usedas confound regressors for each participant. ICA-AROMA
ICA-AROMA (automatic removal of motion artifact) isa recently-introduced, widely-used method for de-noisingusing single-subject ICA (Pruim et al., 2015a,b); we eval-uated ICA-AROMA in confound regression model . Incontrast to other ICA based methods (e.g., ICA-FIX:Salimi-Khorshidi et al. (2014)), it does not require dataset-specific training data. The input to ICA-AROMA is a vox-elwise time series that has been smoothed at 6mm FWHMusing a Gaussian kernel. After decomposing this time se-ries using FSL’s melodic (with model order estimated us-ing the Laplace approximation) (Beckmann et al., 2005),ICA-AROMA uses four features to determine whethereach component corresponds to signal or noise. The first 2features are spatial characteristics of the signal source: (1)4he fraction of the source that falls within a CSF compart-ment and (2) the fraction of the source that falls alongthe edge or periphery of the brain. The remaining fea-tures are derived from the time series of the source: (3) itsmaximal robust correlation with time series derived fromrealignment parameters and (4) its high-frequency spec-tral content. ICA-AROMA includes two de-noising steps.The first de-noising step occurs immediately after classi-fication. All component time series (signal and noise) areincluded as predictors in the linear model, and the resid-ual BOLD time series is obtained via partial regression ofonly the noise time series. A second confound regressionstep occurs after temporal filtering, wherein mean signalsfrom WM and CSF were regressed from the data. Temporal censoring: de-spiking, spike regression, andscrubbing
In addition to regression of nuisance time series, a num-ber of ‘temporal censoring’ approaches were used to iden-tify motion-contaminated volumes in the BOLD time se-ries and reduce their impact on further analysis. These ap-proaches included despiking, spike regression, and scrub-bing. Despiking is a procedure that identifies outliers inthe intensity of each voxel’s detrended BOLD time seriesand then interpolates over these outliers. Despiking wasimplemented in AFNI using the 3dDespike utility (Cox,1996) as part of confound regression model .Unlike despiking, which identifies outliers on a voxel-wise basis, spike regression and scrubbing censor completevolumes based on metrics of subject movement defined a priori . For spike regression, as in Satterthwaite et al.(2013), volumes were flagged for spike regression if theirvolume-to-volume RMS displacement exceeded 0.25mm.Next, as part of confound regression model , k ‘spike’regressors were included as predictor variables in the de-noising model, where k equalled the number of volumesflagged (Satterthwaite et al., 2013). For each flagged timepoint, a unit impulse function that had a value of 1 atthat time point and 0 elsewhere was included as a spikeregressor.For scrubbing, the framewise displacement (FD) (Poweret al., 2012) was computed at each time point as the sumof the absolute values of the derivatives of translationaland rotational motion estimates. If framewise displace-ment (FD) at any point in time exceeded 0.2mm, thenthat time point was flagged for scrubbing. It should benoted that the conversion of FD to RMS displacement isapproximately 2:1, and thus the published criterion forscrubbing has a lower threshold for flagging high-motionvolumes than does spike regression. Scrubbing of BOLDdata was performed iteratively (Power et al., 2014a) aspart of confound regression model . At any stage wherea linear model was applied to the data (for instance, dur-ing detrending procedures), high-motion epochs were tem-porally masked out of the model so as not to influencefit. During temporal filtering, a frequency transform wasused to generate surrogate data with the same phase and spectral characteristics as the unflagged data. This sur-rogate data was used to interpolate over flagged epochsprior to application of the filter. During confound regres-sion, flagged timepoints were excised from the time seriesso as not to contribute to the model fit. For scrubbing (butnot spike regression) if fewer than five contiguous volumeshad unscrubbed data, these volumes were scrubbed andinterpolated as well. Overview of outcome measures
We evaluated each de-noising pipeline according to threebenchmarks. Residual
QC-FC correlations and distance-dependence provided a metric of each pipeline’s efficacy,while loss of temporal DOF provided an estimate of eachpipeline’s efficiency.
QC-FC correlations
In order to estimate the residual relationship betweensubject movement and connectivity after de-noising, wecomputed
QC-FC correlations (quality control / func-tional connectivity) (Power et al., 2015; Satterthwaiteet al., 2012, 2013; Power et al., 2012). While other met-rics have been used in prior reports, including FD-DVARScorrelations, we favor
QC-FC as the most useful metric ofinterest as it directly quantifies the relationship betweenmotion and the primary outcome of interest (rather thantwo quality metrics, as in FD-DVARS). For an extendeddiscussion of the rationale for this measure, see Power et al.(2015).We evaluated
QC-FC relationships within twocommonly-used whole-brain networks, the first consistingof spherical nodes distributed across the brain (Poweret al., 2011) and the second comprising an areal par-cellation of the cerebral cortex (Gordon et al., 2016).For each network, the mean time series for each nodewas calculated from the denoised residual data, and thepairwise Pearson correlation coefficient between nodetime series was used as the network edge weight (Biswalet al., 1995). For each edge, we then computed thecorrelation between the weight of that edge and themean relative RMS motion. To eliminate the potentialinfluence of demographic factors,
QC-FC relationshipswere calculated as partial correlations that accountedfor participant age and sex. We thus obtained, for eachde-noising pipeline, a distribution of
QC-FC correlations.This distribution was used to obtain two measures of thepipeline’s ability to mitigate motion artifact, including:1) the number of edges significantly related to motion,which was computed after using the false discovery rate(FDR) to account for multiple comparisons; and 2) themedian absolute value of all
QC-FC correlations. Allgraphs were generated using ggplot2 in R version 3.2.3;brain renderings were prepared in BrainNet Viewer (Xiaet al., 2013).5 istance-dependent effects of motion
Early work on motion artifact demonstrated that in-scanner motion can bias connectivity estimates betweentwo nodes in a manner that is related to the distancebetween those nodes (Satterthwaite et al., 2012; Poweret al., 2012). Under certain processing conditions, subjectmovement enhances short-distance connections while re-ducing long-distance connections. To determine the resid-ual distance-dependence of subject movement, we firstused the center of mass of each node to obtain a distancematrix D where entry D ij indicates the Euclidean dis-tance between the centers of mass of nodes i and j . Wethen obtained the correlation between the distance sepa-rating each pair of nodes and the QC-FC correlation ofthe edge connecting those nodes; this correlation served asan estimate of the distance-dependence of motion artifact.
Additional degrees of freedom lost in confound regression
Including additional regressors in a confound model re-duces the impact of motion on future analyses, but it isnot without cost. Confound regressors and censoring bothreduce the temporal degrees of freedom (DOF) in data.This loss in temporal DOF may diminish the ability offunctional data acquisitions to sample a subject’s connec-tome or may introduce bias if the loss is variable acrosssubjects. Thus, de-noising strategies ideally limit the lossof temporal DOF, for instance by including fewer, moreefficacious regressors. In the present study, we also as-sessed the number of temporal DOF lost in each confoundregression approach.As in previous work (Pruim et al., 2015a), we assumedthat each time series regressed out and each volume ex-cised from the data constituted a single temporal DOF.Consequently, the loss of temporal DOF was estimatedas the sum of the number of regressors in each confoundmodel and the number of volumes flagged for excision un-der that model. It should be emphasized that the val-ues thus obtained are imperfect estimates. First, becausefunctional MR time series typically exhibit temporal au-tocorrelation, the actual loss in DOF will be less than theestimated loss in DOF. Accordingly, censoring adjacentvolumes does not remove the same number of DOF asdoes censoring volumes separated in time. Furthermore,a temporal bandpass filter was uniformly applied to alldata prior to confound regression; this filtering procedurewould itself have removed additional temporal DOF andelevated the autocorrelation of the data. Because this fil-ter was uniform across all de-noising strategies, it was notconsidered when estimating the loss of additional temporalDOF in each model.
Results
Heterogeneity in confound regression performance
Confound regression strategies typically remove some,but not all, of the artifactual variance that head motion introduces into the BOLD signal. The motion-related arti-fact that survives de-noising can be quantified to provide ametric of pipeline performance. Here, our primary bench-mark of confound regression efficacy was the residual rela-tionship between brain connectivity and subject motion, orthe
QC-FC correlation . We measured
QC-FC correlationsusing two metrics: the percentage of network connectionswhere a significant relationship with motion was present(
Figure 2 ), and the absolute median correlation (AMC)between connection strength and head movement acrossall connections (
Figure 3 ).No preprocessing strategy was completely effective inabolishing the relationship between head movement andconnectivity. However, different approaches exhibitedwidely varying degrees of efficacy. The top four con-found regression strategies included 36 parameters, com-prising an expansion of GSR, tissue-specific regressors(WM, CSF), and realignment parameters. Beyond thisbase 36-parameter model, all censoring techniques pro-vided provided some additional benefit, reducing the num-ber of edges that were significantly related to motion toless than 1%. Convergent results were present across both
QC-FC measures (% edges, AMC) and networks (Power,Gordon) that were evaluated.In contrast, many pipelines performed relatively poorly,leaving a majority of network edges with a residual re-lationship with motion. Specifically, 82% of edges wereimpacted by motion when the least effective method wasused (6 realignment parameters). The commonly used 24-parameter expansion of realignment parameters originallysuggested by Friston et al. (1996) did not provide muchof an improvement (79% edges). Similarly, the local WMregressor model (69% edges) and tCompCor (48% edges)also resulted in substantial residual
QC-FC correlations.In fact, these methods performed worse than a basic 2-parameter model composed of mean WM and CSF signals(33% edges). Finally, several methods demonstrated in-termediate performance, with 1-20% of edges impacted bymotion. This middle group included methods as disparateas aCompCor (5% edges), ICA-AROMA (12% edges), andthe classic 9-parameter confound regression model whichincluded GSR (4% edges).
Variability in distance-dependent motion artifact afterconfound regression
Our second benchmark quantified the distance-dependent motion artifact that was present in data pro-cessed by each pipeline (
Figure 4 ). We observed thatdistance-dependence was present even under conditionswhere artifact magnitude was attenuated. For example,though the 36-parameter model was among the most ef-fective in attenuating
QC-FC relationships, its applicationrevealed strongly distance-dependent artifact. Examina-tion of graphs that plot
QC-FC by Euclidean distance (see
Figure 4C ) revealed that this is due to effective mitiga-tion of motion artifact for long-range but not short-rangeconnections.6istance-dependence was most prominent in modelsthat included GSR, but did not include censoring (e.g.,9-parameter and 36-parameter models). However, despitethe lack of global signal in the aCompCor or tCompcormodels, data returned from both of these component-based approaches revealed substantial distance-dependentartifact. Notably, inclusion of censoring consistently re-duced distance-dependence, although scrubbing was some-what more effective than spike regression or voxelwise de-spiking.The top performing method according to this bench-mark was ICA-AROMA, which completely abolished anydistance-dependence of residual motion artifact. In otherwords, the motion artifact that was still present in the dataafter ICA-AROMA impacted all connections in a mannerthat was not dependent on the spatial separation betweennodes. There was similar lack of distance-dependence inthe wmLocal model, although as noted above this modeldid not provide effective de-noising according to
QC-FC benchmarks. Despite the presence of the global signal,the 36-parameter model that included scrubbing also per-formed quite well.
Effective preprocessing strategies use many additional de-grees of freedom
Perhaps unsurprisingly, the preprocessing strategiesthat consistently reduced both
QC-FC correlations anddistance-dependence were also among the costliest interms of loss of temporal degrees of freedom (
Figure 5 ).By definition, the 36-parameter models included a highfixed number of regressors. Furthermore, models thatadditionally included censoring resulted in a substantialadditional loss of data that varied across subjects. ICA-AROMA also had a variable loss of DOF, but of a lowermagnitude than censoring or high-parameter confound re-gression.
Discussion
In response to rapid evolution of confound regressionstrategies available for the mitigation of motion artifact,in this report we evaluated 12 commonly-used pipelines.Results indicate that there is substantial heterogeneity inthe performance of these confound regression techniquesacross all measures evaluated. The context, implications,and limitations of these findings are discussed below.
Confound regression techniques have substantial perfor-mance variability
We evaluated confound regression strategies accordingto three intuitive benchmarks that were selected to capturedifferent domains of effectiveness. These included
QC-FC associations, distance-dependence of motion artifact, andadditional degrees of freedom lost in confound regression.Across each benchmark, there was a striking heterogeneityin pipeline performance. Notably, five of the top six confound regression ap-proaches included GSR. This effect was consistent in bothnetworks we evaluated. The effectiveness of GSR is mostlikely due to the nature of motion artifact itself: in-scannerhead motion tends to induce widespread reductions in sig-nal intensity across the entire brain parenchyma (see Sat-terthwaite et al. (2013), Figure 4). As discussed in de-tail elsewhere (Power et al., Under Revision) this effectis highly reproducible across datasets, and is effectivelycaptured by time series regression of the global signal.Beyond GSR, a second strategy that clearly minimizes
QC-FC relationships is temporal censoring. We evaluatedthree censoring variants, including scrubbing, spike regres-sion, and de-spiking. Compared to spike regression andde-spiking, scrubbing appears to be more effective in re-moving distance-dependent artifact in this dataset. This ismost likely due to the explicit tension between data qual-ity and data quantity: because of the lower threshold forscrubbing than spike regression (due to differences in FDvs. RMS measures of motion; see Figure 9C in Yan et al.(2013a)), more low-quality data was excised during scrub-bing. Furthermore, scrubbing includes a criterion to notleave isolated epochs ( < QC-FC relationships, they exhibitedopposite effects on distance-dependence. While censor-ing techniques appear to consistently reduce the presenceof distance-dependence, GSR is associated with increaseddistance-dependence. However, it should be emphasizedthat the distance-dependence associated with GSR is notthe result of worsening associations with motion in certainconnections. Rather, the distance-dependence seen withGSR stems from differential de-noising efficacy, wherebymotion artifact is more effectively minimized for long-distance connections than for short-range connections.Certain models such as the local WM regression approach(Jo et al., 2013) thus have minimal distance-dependence,but this is a consequence of lack of efficacy across all dis-tances. In contrast, ICA-AROMA (Pruim et al., 2015b,a)reduced motion to a moderate degree over all connec-tion distances, resulting in almost no distance-dependence.However, while clearly an improvement over some othermethods, data processed using ICA-AROMA was nois-ier than other methods which included GSR or censor-ing, and resulting networks contained a substantial num-ber of edges impacted by motion. Somewhat to our sur-prise, benchmark results for aCompCor (Behzadi et al.,2007; Muschelli et al., 2014) were most similar to models7hat included GSR. Alone among models where GSR wasnot included, aCompCor was both relatively effective inthe mitigation of residual motion (5% of edges impacted)and also exhibited substantial distance-dependence (e.g., r = − . Trade-offs of confound regression approaches: implicationsfor investigators
The current results emphasize two clear trade-offs inthe choice of confound regression strategy. First, pipelinesthat include global signal regression tend to be more effec-tive at minimizing
QC-FC relationships, but at the costof some increase in distance-dependence. As noted above,for minimizing
QC-FC relationships, nearly all of the topstrategies (except aCompCor) included GSR. Conversely,the two techniques that had the most substantial distance-dependence (the 9-regressor and 36-regressor methods)both included GSR. Second, censoring techniques providea clear benefit in reducing
QC-FC relationships and addi-tionally tend to attenuate distance-dependence. However,by definition, removing contaminated volumes results inless data and loss of degrees of freedom.These two trade-offs suggest that a single confound re-gression strategy is unlikely to be optimal for every study.For example, in studies of network organization, the pres-ence of both anti-correlations and distance-dependent ar-tifact resulting from GSR may result in an undesirableimpact on nodal degree distribution (Yan et al., 2013b).In contrast, for studies of group or individual differ-ences, minimizing
QC-FC relationships is likely to be ofparamount importance so as to limit the possibility thatinference is driven by artifactual signals. This concern isparticularly relevant for studies of development or clinicalsub-groups where motion is systematically related to thesubject-level variable of interest (e.g., age, disease status).Co-varying for motion at the group level is unlikely to bea panacea for such studies when inadequate subject-leveltime series de-noising is employed, as prior work (Poweret al., 2014a) has suggested that motion effects at thegroup level may potentially both be nonlinear and varyacross sub-samples in a manner that is difficult to predict.However, aggressive volume censoring may be problematicin datasets with relatively brief acquisitions. In datasetswhere long time series are acquired, such as multi-band ac-quisitions (Feinberg et al., 2010) and intensive acquisitionsof single subjects (Laumann et al., 2015), loss of tempo-ral degrees of freedom is less likely to be a major concern.The 36-parameter models without volume censoring offeruniformity, as does randomly or systematically censoringadditional volumes until all subjects retain approximatelythe same degrees of freedom.
Limitations
Several limitations of the current approach should benoted. One of the principal challenges in evaluating the performance of de-noising approaches is the lack of a noise-free ground truth. Our primary benchmark of confoundregression performance assumes that mitigation of the re-lationship between QC (participant motion) and FC (i.e.,the imaging measurement) is desirable. To the degree thatin-scanner motion itself represents a biologically informa-tive phenotype, this approach will mistake signal for noise.Indeed, prior data suggests that this may sometimes bethe case. For example, Zeng et al. (2014) found specificchanges in connectivity for participants who had generallyhigh levels of motion, even on scans where motion was low.However, without multiple scans which allow such carefuldissociation, most studies are incapable of disambiguatingthe large confounding effects of motion on connectivity.Second, in place of
QC-FC relationships, one could focuson alternative benchmarks such as test-retest reliability(Zuo et al., 2014). Reliability is certainly of interest, butto the degree that motion tends to be highly correlatedwithin individuals across scanning sessions, there is a sub-stantial potential for the presence of consistent motion ar-tifact across sessions to artificially inflate reliability, anddiminish the biological relevance of observed results. Athird and related concern is that certain de-noising meth-ods could conceivably both minimize
QC-FC relationshipsand even enhance reliability by aggressively removing bothsignal and noise, but in the process diminish sensitivity tomeaningful individual differences. Indeed, one prior studydemonstrated the association between canonical restingstate networks and randomly generated confound parame-ters (Bright and Murphy, 2015). This concern is somewhatmitigated by other studies in independent datasets, whichsuggest that higher-order confound regressors improve theconfound regression model fit (Yan et al., 2013a; Satterth-waite et al., 2013), while random regressors do not (see Fig-ure 8 in Satterthwaite et al. (2013)). Nonetheless, the ten-sion between the goals of noise reduction and sensitivity toindividual differences remains an outstanding issue for thefield. Fourth, while our evaluation included many of themost commonly used techniques, other approaches whichrequire substantial training or parameter selection (i.e.,ICA-FIX (Salimi-Khorshidi et al., 2014; Griffanti et al.,2014), wavelet de-spiking (Patel et al., 2014)) may be valu-able and merit further consideration. Fifth and finally, itshould be noted that while improvements in image acqui-sition (including multi-echo techniques) may not salvageexisting motion-contaminated data, it is likely that theywill change the methodological landscape of connectivityresearch moving forward (Kundu et al., 2012, 2013; Brightand Murphy, 2013).
Conclusions
Taken together, the present results underline the per-formance heterogeneity of recently-introduced, commonly-used confound regression methods. In selecting amongthese methods, investigators should be aware of the rel-ative strengths and weaknesses of each approach, and un-derstand how processing strategy may impact inference.8learly, the relative merit of each approach will vary byresearch question and study design. Perhaps most im-portantly, as has been emphasized in nearly every otherstudy of motion artifact, the choice of confound regres-sion strategy is often dwarfed in importance by the needto transparently report and evaluate the impact of mo-tion in each dataset. At a minimum, this includes report-ing the relationship between motion artifact and not onlysubject phenotypes (e.g., group, age, symptom or cogni-tive score) but also the functional connectivity phenotypesbeing considered. In the context of such data, the distinc-tion between observed results and the impact of motionartifact can be understood. Such transparency bolstersconfidence in reported findings, but also will likely tendto emphasize the remaining challenges for de-noising go-ing forward. Especially when considered in the context ofthe the rapid evolution of available techniques since 2012,there is no doubt that innovations in post-processing con-found regression strategies will continue.
Acknowledgements
ReferencesReferences
Avants, B.B., Tustison, N.J., Song, G., Cook, P.A., Klein, A., Gee,J.C., 2011a. A reproducible evaluation of ANTs similarity metric performance in brain image registration. NeuroImage 54, 2033–2044. doi: .Avants, B.B., Tustison, N.J., Wu, J., Cook, P.A., Gee, J.C., 2011b.An open source multivariate framework for N-tissue segmenta-tion with evaluation on public data. Neuroinformatics 9, 381–400.doi: .Baker, J.T., Holmes, A.J., Masters, G.a., Yeo, B.T.T., Krienen, F.,Buckner, R.L., ¨Ong¨ur, D., 2014. Disruption of cortical associationnetworks in schizophrenia and psychotic bipolar disorder. JAMApsychiatry 71, 109–18. doi: , arXiv:15334406 .Beckmann, C.F., Deluca, M., Devlin, J.T., Smith, S.M., 2005. Inves-tigations into resting-state connectivity using independent compo-nent analysis. Philos Trans R Soc Lond B Biol Sci 360, 1001–13.doi: .Behzadi, Y., Restom, K., Liau, J., Liu, T.T., 2007. A compo-nent based noise correction method (CompCor) for BOLD andperfusion based fMRI. NeuroImage 37, 90–101. doi: , arXiv:NIHMS150003 .Biswal, B., Yetkin, F.Z., Haughton, V.M., Hyde, J.S., 1995. Func-tional connectivity in the motor cortex of resting human brainusing echo-planar MRI. Magnetic Resonance in Medicine 34, 537–541. doi: .Bright, M.G., Murphy, K., 2013. Removing motion and physiolog-ical artifacts from intrinsic BOLD fluctuations using short echodata. NeuroImage 64, 526–537. doi: .Bright, M.G., Murphy, K., 2015. Is fMRI ”noise” really noise?Resting state nuisance regressors remove variance with networkstructure. NeuroImage 114, 158–169. doi: .Buckner, R.L., Andrews-Hanna, J.R., Schacter, D.L., 2008. Thebrain’s default network: Anatomy, function, and relevance to dis-ease. Annals of the New York Academy of Sciences 1124, 1–38.doi: .Chai, X.J., Casta˜n´an, A.N., ¨Ong¨ur, D., Whitfield-Gabrieli, S., 2012.Anticorrelations in resting state networks without global signal re-gression. NeuroImage 59, 1420–1428. doi: .Cox, R., 1996. AFNI: Software for analysis and visualization offunctional magnetic resonance neuroimages. Comput Biomed Res29, 162–173.Craddock, R.C., Jbabdi, S., Yan, C.G., Vogelstein, J.T., Castel-lanos, F.X., Di Martino, A., Kelly, C., Heberlein, K., Colcombe,S., Milham, M.P., 2013. Imaging human connectomes at themacroscale. Nature Methods 10, 524–539. doi: , arXiv:NIHMS150003 .Damoiseaux, J.S., Rombouts, S.A.R.B., Barkhof, F., Scheltens, P.,Stam, C.J., Smith, S.M., Beckmann, C.F., 2006. Consistentresting-state networks across healthy subjects. Proceedings of theNational Academy of Sciences of the United States of America103, 13848–53. doi: .DiMartino, A., Fair, D.A., Kelly, C., Satterthwaite, T.D., Castel-lanos, F.X., Thomason, M.E., Craddock, R.C., Luna, B., Lev-enthal, B.L., Zuo, X.N., Milham, M.P., 2014. Unraveling themiswired connectome: A developmental perspective. Neuron 83,1335–1353. doi: .Dosenbach, N.U.F., Nardos, B., Cohen, A.L., Fair, D.A., Power,D., Church, J.a., Nelson, S.M., Wig, G.S., Vogel, A.C., Lessov-schlaggar, C.N., Barnes, K.A., Dubis, J.W., Feczko, E., Coalson,R.S., Jr, J.R.P., Barch, D.M., Petersen, S.E., Schlaggar, B.L.,2011. Prediction of individual brain maturity using fMRI. Science329, 1358–1361. doi: .Fair, D.A., Cohen, A.L., Dosenbach, N.U.F., Church, J.A., Miezin,F.M., Barch, D.M., Raichle, M.E., Petersen, S.E., Schlaggar, B.L.,2008. The maturing architecture of the brain’s default network.Proceedings of the National Academy of Sciences of the UnitedStates of America 105, 4028–4032. doi: .Fair, D.A., Nigg, J.T., Iyer, S., Bathula, D., Mills, K.L., Dosenbach,N.U.F., Schlaggar, B.L., Mennes, M., Gutman, D., Bangaru, S.,Buitelaar, J.K., Dickstein, D.P., Di Martino, A., Kennedy, D.N., elly, C., Luna, B., Schweitzer, J.B., Velanova, K., Wang, Y.F.,Mostofsky, S., Castellanos, F.X., Milham, M.P., 2012. Distinctneural signatures detected for ADHD subtypes after controllingfor micro-movements in resting state functional connectivity MRIdata. Frontiers in systems neuroscience 6, 80. doi: .Fair, D.A., Posner, J., Nagel, B.J., Bathula, D., Dias, T.G.C., Mills,K.L., Blythe, M.S., Giwa, A., Schmitt, C.F., Nigg, J.T., 2010.Atypical default network connectivity in youth with attention-deficit/hyperactivity disorder. Biological psychiatry 68, 1084–91.doi: .Feinberg, D.A., Moeller, S., Smith, S.M., Auerbach, E.J., Ramanna,S., Glasser, M.F., Miller, K.L., Ugurbil, K., Yacoub, E.S., 2010.Multiplexed echo planar imaging for sub-second whole brain fMRIand fast diffusion imaging. PLoS ONE 5, e15710. doi: .Fox, M.D., Raichle, M.E., 2007. Spontaneous fluctuations in brainactivity observed with functional magnetic resonance imaging.Nat Rev Neurosci 8, 700–711. doi: .Fox, M.D., Snyder, A.Z., Vincent, J.L., Corbetta, M., Van Es-sen, D.C., Raichle, M.E., 2005. The human brain is intrinsi-cally organized into dynamic, anticorrelated functional networks.Proceedings of the National Academy of Sciences of the UnitedStates of America 102, 9673–8. doi: , arXiv:NIHMS150003 .Fox, M.D., Zhang, D., Snyder, A.Z., Raichle, M.E., 2009. Theglobal signal and observed anticorrelated resting state brain net-works. Journal of neurophysiology 101, 3270–3283. doi: .Friston, K., Williams, S., Howard, R., Frackowiak, R.S.J., Turner,R., 1996. Movement-related effects in { fMRI } time-series. Mag-netic Resonance in Medicine 35, 346–355.Gordon, E.M., Laumann, T.O., Adeyemo, B., Huckins, J.F., Kel-ley, W.M., Petersen, S.E., 2016. Generation and evaluation of acortical area parcellation from resting-state correlations. CerebralCortex 26, 288–303. doi: .Gotts, S.J., Saad, Z.S., Jo, H.J., Wallace, G.L., Cox, R.W., Martin,A., 2013. The perils of global signal regression for group compar-isons: a case study of Autism Spectrum Disorders. Frontiers inHuman Neuroscience 7, 356. doi: .Greve, D.N., Fischl, B., 2009. Accurate and robust brain imagealignment using boundary-based registration. NeuroImage 48, 63–72. doi: .Griffanti, L., Salimi-Khorshidi, G., Beckmann, C.F., Auerbach, E.J.,Douaud, G., Sexton, C.E., Zsoldos, E., Ebmeier, K.P., Filippini,N., Mackay, C.E., Moeller, S., Xu, J., Yacoub, E., Baselli, G.,Ugurbil, K., Miller, K.L., Smith, S.M., 2014. ICA-based artefactremoval and accelerated fMRI acquisition for improved restingstate network imaging. NeuroImage 95, 232–247. doi: .Gur, R.E., Kaltman, D., Melhem, E.R., Ruparel, K., Prabhakaran,K., Riley, M., Yodh, E., Hakonarson, H., Satterthwaite, T., Gur,R.C., 2013. Incidental findings in youths volunteering for brainMRI research. American Journal of Neuroradiology 34, 2021–2025.doi: .Hahamy, A., Calhoun, V., Pearlson, G., Harel, M., Stern, N., Attar,F., Malach, R., Salomon, R., 2014. Save the global: global signalconnectivity as a tool for studying clinical populations with func-tional magnetic resonance imaging. Brain connectivity 4, 395–403.doi: .Hallquist, M.N., Hwang, K., Luna, B., 2013. The nuisance ofnuisance regression: Spectral misspecification in a common ap-proach to resting-state fMRI preprocessing reintroduces noiseand obscures functional connectivity. NeuroImage 82, 208–225.doi: , arXiv:NIHMS150003 .Jenkinson, M., Bannister, P., Brady, M., Smith, S., 2002.Improved optimization for the robust and accurate linearregistration and motion correction of brain images. Neu-roImage 17, 825–841. doi: , arXiv:arXiv:1011.1669v3 .Jo, H.J., Gotts, S.J., Reynolds, R.C., Bandettini, P.A., Martin, A., Cox, R.W., Saad, Z.S., 2013. Effective preprocessing proceduresvirtually eliminate distance-dependent motion artifacts in restingstate FMRI. Journal of Applied Mathematics 2013. doi: .Klein, A., Andersson, J., Ardekani, B.A., Ashburner, J., Avants,B., Chiang, M.C., Christensen, G.E., Collins, D.L., Gee, J., Hel-lier, P., Song, J.H., Jenkinson, M., Lepage, C., Rueckert, D.,Thompson, P., Vercauteren, T., Woods, R.P., Mann, J.J., Parsey,R.V., 2009. Evaluation of 14 nonlinear deformation algorithms ap-plied to human brain MRI registration. NeuroImage 46, 786–802.doi: .Kundu, P., Brenowitz, N.D., Voon, V., Worbe, Y., V´ertes, P.E.,Inati, S.J., Saad, Z.S., Bandettini, P.A., Bullmore, E.T., 2013.Integrated strategy for improving functional connectivity mappingusing multiecho fMRI. Proceedings of the National Academy ofSciences of the Unites States of America 110, 16187–16192. doi: , arXiv:arXiv:1408.1149 .Kundu, P., Inati, S.J., Evans, J.W., Luh, W.M., Bandettini, P.A.,2012. Differentiating BOLD and non-BOLD signals in fMRI timeseries using multi-echo EPI. NeuroImage 60, 1759–1770. doi: .Laumann, T.O., Gordon, E.M., Adeyemo, B., Snyder, A.Z., Joo,S.J., Chen, M.Y., Gilmore, A.W., McDermott, K.B., Nelson,S.M., Dosenbach, N.U.F., Schlaggar, B.L., Mumford, J.A., Pol-drack, R.A., Petersen, S.E., 2015. Functional system and arealorganization of a highly sampled individual human brain. Neuron87, 658–671. doi: .Macey, P.M., Macey, K.E., Kumar, R., Harper, R.M., 2004. Amethod for removal of global effects from fMRI time series. Neu-roImage 22, 360–366. doi: .Merikangas, K.R., He, J.P., Brody, D., Fisher, P., Bourdon, K., Ko-retz, D.S., 2010. Prevalence and treatment of mental disordersamong US children in the 20012004 NHANES. Pediatrics 125,75–81. doi: .Murphy, K., Birn, R.M., Handwerker, D.A., Jones, T.B., Bandet-tini, P.A., 2009. The impact of global signal regression on restingstate correlations: Are anti-correlated networks introduced? Neu-roImage 44, 893–905. doi: , arXiv:NIHMS150003 .Muschelli, J., Nebel, M.B., Caffo, B.S., Barber, A.D., Pekar, J.J.,Mostofsky, S.H., 2014. Reduction of motion-related artifacts inresting state fMRI using aCompCor. NeuroImage 96, 22–35.doi: .Patel, A.X., Kundu, P., Rubinov, M., Jones, P.S., V´ertes, P.E.,Ersche, K.D., Suckling, J., Bullmore, E.T., 2014. A waveletmethod for modeling and despiking motion artifacts from resting-state fMRI time series. NeuroImage 95, 287–304. doi: .Power, J.D., Barnes, K.A., Snyder, A.Z., Schlaggar, B.L., Petersen,S.E., 2012. Spurious but systematic correlations in functionalconnectivity MRI networks arise from subject motion. Neu-roImage 59, 2142–2154. doi: , arXiv:NIHMS150003 .Power, J.D., Cohen, A.L., Nelson, S.M., Wig, G.S., Barnes, K.A.,Church, J.A., Vogel, A.C., Laumann, T.O., Miezin, F.M., Schlag-gar, B.L., Petersen, S.E., 2011. Functional network organizationof the human brain. Neuron 72, 665–678. doi: .Power, J.D., Mitra, A., Laumann, T.O., Snyder, A.Z., Schlag-gar, B.L., Petersen, S.E., 2014a. Methods to detect, charac-terize, and remove motion artifact in resting state fMRI. Neu-roImage 84, 320–341. doi: , arXiv:NIHMS150003 .Power, J.D., Schlaggar, B.L., Petersen, S.E., 2014b. Studying brainorganization via spontaneous fMRI signal. Neuron 84, 681–696.doi: , arXiv:NIHMS150003 .Power, J.D., Schlaggar, B.L., Petersen, S.E., 2015. Recent progressand outstanding issues in motion correction in resting state fMRI.NeuroImage 105, 536–551. doi: , arXiv:NIHMS150003 .Pruim, R.H.R., Mennes, M., Buitelaar, J.K., Beckmann, C.F., 2015a. valuation of ICA-AROMA and alternative strategies for motionartifact removal in resting state fMRI. NeuroImage 112, 278–287.doi: .Pruim, R.H.R., Mennes, M., van Rooij, D., Llera, A., Buitelaar, J.K.,Beckmann, C.F., 2015b. ICA-AROMA: A robust ICA-based strat-egy for removing motion artifacts from fMRI data. NeuroImage112, 267–277. doi: .Saad, Z.S., Gotts, S.J., Murphy, K., Chen, G., Jo, H.J., Martin, A.,Cox, R.W., 2012. Trouble at rest: how correlation patterns andgroup differences become distorted after global signal regression.Brain connectivity 2, 25–32. doi: .Salimi-Khorshidi, G., Douaud, G., Beckmann, C.F., Glasser, M.F.,Griffanti, L., Smith, S.M., 2014. Automatic denoising of func-tional MRI data: Combining independent component analysisand hierarchical fusion of classifiers. NeuroImage 90, 449–468.doi: , arXiv:NIHMS150003 .Satterthwaite, T.D., Connolly, J.J., Ruparel, K., Calkins, M.E.,Jackson, C., Elliott, M.A., Roalf, D.R., Hopsona, R., Prab-hakaran, K., Behr, M., Qiu, H., Mentch, F.D., Chiavacci, R.,Sleiman, P.M.A., Gur, R.C., Hakonarson, H., Gur, R.E., 2016.The Philadelphia Neurodevelopmental Cohort: A publicly avail-able resource for the study of normal and abnormal brain devel-opment in youth. NeuroImage 124, 1115–1119. doi: .Satterthwaite, T.D., Elliott, M.A., Gerraty, R.T., Ruparel, K., Loug-head, J., Calkins, M.E., Eickhoff, S.B., Hakonarson, H., Gur,R.C., Gur, R.E., Wolf, D.H., 2013. An improved frameworkfor confound regression and filtering for control of motion arti-fact in the preprocessing of resting-state functional connectivitydata. NeuroImage 64, 240–256. doi: , arXiv:NIHMS150003 .Satterthwaite, T.D., Elliott, M.A., Ruparel, K., Loughead, J., Prab-hakaran, K., Calkins, M.E., Hopson, R., Jackson, C., Keefe, J.,Riley, M., Mentch, F.D., Sleiman, P., Verma, R., Davatzikos, C.,Hakonarson, H., Gur, R.C., Gur, R.E., 2014. Neuroimaging ofthe Philadelphia Neurodevelopmental Cohort. NeuroImage 86,544–553. doi: .Satterthwaite, T.D., Wolf, D.H., Loughead, J., Ruparel, K., Elliott,M.A., Hakonarson, H., Gur, R.C., Gur, R.E., 2012. Impact of in-scanner head motion on multiple measures of functional connec-tivity: Relevance for studies of neurodevelopment in youth. Neu-roImage 60, 623–632. doi: , arXiv:NIHMS150003 .Smith, S.M., Vidaurre, D., Beckmann, C.F., Glasser, M.F., Jenk-inson, M., Miller, K.L., Nichols, T.E., Robinson, E.C., Salimi-Khorshidi, G., Woolrich, M.W., Barch, D.M., U??urbil, K., VanEssen, D.C., 2013. Functional connectomics from resting-statefMRI. Trends in Cognitive Sciences 17, 666–682. doi: .Tustison, N.J., Avants, B.B., Cook, P.A., Zheng, Y., Egan, A.,Yushkevich, P.A., Gee, J.C., 2010. N4ITK: Improved N3 bias cor-rection. IEEE Transactions on Medical Imaging 29, 1310–1320.doi: .Tustison, N.J., Cook, P.A., Klein, A., Song, G., Das, S.R., Duda,J.T., Kandel, B.M., van Strien, N., Stone, J.R., Gee, J.C., Avants,B.B., 2014. Large-scale evaluation of ANTs and FreeSurfer corticalthickness measurements. NeuroImage 99, 166–179. doi: .Van Dijk, K.R.A., Hedden, T., Venkataraman, A., Evans, K.C.,Lazar, S.W., Buckner, R.L., 2010. Intrinsic functional con-nectivity as a tool for human connectomics: theory, proper-ties, and optimization. Journal of neurophysiology 103, 297–321.doi: .Van Dijk, K.R.A., Sabuncu, M.R., Buckner, R.L., 2012. The in-fluence of head motion on intrinsic functional connectivity MRI.NeuroImage 59, 431–438. doi: , arXiv:NIHMS150003 .Wang, H., Cao, Y., Syeda-mahmood, T., 2014. Multi-atlas segmen-tation with learning-based label fusion. Machine Learning in Med-ical Imaging 35, 256–263. doi: .Xia, M., Wang, J., He, Y., 2013. BrainNet Viewer: A network visualization tool for human brain connectomics. PLoS ONE 8.doi: .Yan, C.G., Cheung, B., Kelly, C., Colcombe, S., Craddock, R.C., DiMartino, A., Li, Q., Zuo, X.N., Castellanos, F.X., Milham, M.P.,2013a. A comprehensive assessment of regional variation in theimpact of head micromovements on functional connectomics. Neu-roImage 76, 183–201. doi: , arXiv:NIHMS150003 .Yan, C.G., Craddock, R.C., He, Y., Milham, M.P., 2013b. Ad-dressing head motion dependencies for small-world topologies infunctional connectomics. Frontiers in human neuroscience 7, 910.doi: .Yan, C.G., Craddock, R.C., Zuo, X.N., Zang, Y.F., Milham, M.P.,2013c. Standardizing the intrinsic brain: Towards robust measure-ment of inter-individual variation in 1000 functional connectomes.NeuroImage 80, 246–262. doi: , arXiv:NIHMS150003 .Yang, G.J., Murray, J.D., Repovs, G., Cole, M.W., Savic, A.,Glasser, M.F., Pittenger, C., Krystal, J.H., Wang, X.J., Pearl-son, G.D., Glahn, D.C., Anticevic, A., 2014. Altered global brainsignal in schizophrenia. Proceedings of the National Academy ofSciences of the United States of America 111, 7438–43. doi: .Yeo, B.T., Krienen, F.M., Sepulcre, J., Sabuncu, M.R., Lashkari, D.,Hollinshead, M., Roffman, J.L., Smoller, J.W., Zollei, L., Poli-meni, J.R., Fischl, B., Liu, H., Buckner, R.L., 2011. The organi-zation of the human cerebral cortex estimated by intrinsic func-tional connectivity. Journal of neurophysiology 106, 1125–1165.doi: Zeng, L.L., Wang, D., Fox, M.D., Sabuncu, M., Hu, D., Ge, M.,Buckner, R.L., Liu, H., 2014. Neurobiological basis of head motionin brain imaging. Proceedings of the National Academy of Sciencesof the United States of America 111, 6058–62. doi: .Zuo, X.N., Anderson, J.S., Bellec, P., Birn, R.M., Biswal, B.B.,Blautzik, J., Breitner, J.C.S., Buckner, R.L., Calhoun, V.D.,Castellanos, F.X., Chen, A., Chen, B., Chen, J., Chen, X., Col-combe, S.J., Courtney, W., Craddock, R.C., Di Martino, A.,Dong, H.M., Fu, X., Gong, Q., Gorgolewski, K.J., Han, Y., He,Y., He, Y., Ho, E., Holmes, A., Hou, X.H., Huckins, J., Jiang, T.,Jiang, Y., Kelley, W., Kelly, C., King, M., LaConte, S.M., Lain-hart, J.E., Lei, X., Li, H.J., Li, K., Li, K., Lin, Q., Liu, D., Liu, J.,Liu, X., Liu, Y., Lu, G., Lu, J., Luna, B., Luo, J., Lurie, D., Mao,Y., Margulies, D.S., Mayer, A.R., Meindl, T., Meyerand, M.E.,Nan, W., Nielsen, J.A., O’Connor, D., Paulsen, D., Prabhakaran,V., Qi, Z., Qiu, J., Shao, C., Shehzad, Z., Tang, W., Villringer,A., Wang, H., Wang, K., Wei, D., Wei, G.X., Weng, X.C., Wu,X., Xu, T., Yang, N., Yang, Z., Zang, Y.F., Zhang, L., Zhang, Q.,Zhang, Z., Zhang, Z., Zhao, K., Zhen, Z., Zhou, Y., Zhu, X.T.,Milham, M.P., 2014. An open science resource for establishing re-liability and reproducibility in functional connectomics. Scientificdata 1, 140049. doi: . igure 1: Schematic of the 12 de-noising models evaluated in the present study.
For each of the 12 models indexed at left,the table details what processing procedures and confound regressors were included in the model. De-noising models were selectedfrom the functional connectivity literature and represented a range of strategies. igure 2: Number of edges significantly related to motion after de-noising.
Successful de-noising strategies reduced therelationship between connectivity and motion. The number of edges (network connections) for which this relationship persists providesevidence of a pipeline’s efficacy. A, The percentage of edges significantly related to motion in a 264-node network defined by Poweret al. (2011). Fewer significant edges is indicative of better performance. B, The percentage of edges significantly related to motion ina second, 333-node network defined by Gordon et al. (2016). C, Renderings of significant edges for each de-noising strategy, rankedaccording to efficacy. Strategies that include regression of the mean global signal are framed in blue and consistently ranked as thebest performers. igure 3: Residual QC-FC correlations after de-noising.
The absolute median
QC-FC correlation is another measure of therelationship between connectivity and motion. A, The absolute median correlation between functional connectivity and motion in a 264-node network defined by Power et al. (2011). A lower absolute median correlation indicates better performance. B, The absolute mediancorrelation between functional connectivity and motion in a second, 333-node network defined by Gordon et al. (2016). C, Distributionsof all edgewise
QC-FC correlations after each de-noising strategy, ranked according to efficacy. Results largely recapitulated thosereported in Figure 2, with GSR-based approaches (blue frame) collectively exhibiting the best performance. Whereas approaches thatincluded more regressors generally yielded a narrower distribution, those approaches that included GSR tended to shift the distribution’scenter toward 0. igure 4: Distance-dependence of motion artifact after de-noising.
The magnitude of motion artifact varies with the Euclideandistance separating a pair of nodes, with closer nodes generally exhibiting greater impact of motion on connectivity. A, The residualdistance-dependence of motion artifact in a 264-node network defined by Power et al. (2011) following confound regression. B, Theresidual distance-dependence of motion artifact in a second, 333-node network defined by Gordon et al. (2016). C, Density plotsindicating the relationship between the Euclidean distance separating each pair of nodes (x-axis) and the
QC-FC correlation ofthe edge connecting those nodes (y-axis). The overall trend lines for each de-noising strategy, from which distance-dependence iscomputed, are indicated in red. The best performing models either excised high-motion volumes (36-parameter + scrubbing) or usedmore localized regressors (ICA-AROMA and wmLocal). In general, approaches that made use of GSR without censoring resulted insubstantial distance-dependence. This effect was driven by differential efficacy of de-noising, with effective de-noising for long rangeconnections but not short-range connections. igure 5: Estimated loss of temporal degrees of freedom for each pipeline evaluated.