Evaluation of Denoising Strategies to Address Motion-Correlated Artifacts in Resting-State Functional Magnetic Resonance Imaging Data from the Human Connectome Project

Abstract

Like all resting-state functional connectivity data, the data from the Human Connectome Project (HCP) are adversely affected by structured noise artifacts arising from head motion and physiological processes. Functional connectivity estimates (Pearson's correlation coefficients) were inflated for high-motion time points and for high-motion participants. This inflation occurred across the brain, suggesting the presence of globally distributed artifacts. The degree of inflation was further increased for connections between nearby regions compared with distant regions, suggesting the presence of distance-dependent spatially specific artifacts. We evaluated several denoising methods: censoring high-motion time points, motion regression, the FMRIB independent component analysis-based X-noiseifier (FIX), and mean grayordinate time series regression (MGTR; as a proxy for global signal regression). The results suggest that FIX denoising reduced both types of artifacts, but left substantial global artifacts behind. MGTR significantly reduced global artifacts, but left substantial spatially specific artifacts behind. Censoring high-motion time points resulted in a small reduction of distance-dependent and global artifacts, eliminating neither type. All denoising strategies left differences between high- and low-motion participants, but only MGTR substantially reduced those differences. Ultimately, functional connectivity estimates from HCP data showed spatially specific and globally distributed artifacts, and the most effective approach to address both types of motion-correlated artifacts was a combination of FIX and MGTR.

Introduction

The Human Connectome Project (HCP) endeavors to reveal variations in connectivity and their relationship to behavior, function, and genetics in 1200 healthy participants (Van Essen et al., 2013). HCP advanced cutting-edge pulse sequences to provide resting-state functional magnetic resonance imaging (rfMRI) data with high spatial and temporal resolution and whole-brain coverage (Ugurbil et al., 2013). Nonetheless, HCP rfMRI data, like most rfMRI data, are likely contaminated by artifacts resulting from a number of influences, including head motion, scanner-related issues, and physiological processes related to cardiac, respiratory, and pCO₂ fluctuations. The intention of the current work was to focus on artifacts correlated with head motion and the effectiveness of methods designed to reduce such artifacts.

Motion-correlated artifacts can bias our understanding of functional networks and their relationship with individual and group difference variables (Power et al., 2012; Van Dijk et al., 2012; Yan et al., 2013b). Prior studies using lower resolution rfMRI data provided approaches to address motion-correlated artifacts (Jo et al., 2013; Muschelli et al., 2014; Satterthwaite et al., 2013), reviewed in Power and coworkers (2015). These denoising strategies demonstrated varying degrees of efficacy, but it is unclear whether they will benefit higher resolution HCP rfMRI data to the same degree.

Indicators suggestive of motion-correlated artifact

In this article, we took four approaches to interrogate the relationship between head motion and artifacts in rfMRI data (c.f., Power et al., 2014): (1) intensity fluctuations in time series data, (2) distance-dependent artifacts, (3) elevated differences between low- and high-motion groups, and (4) relationships between head motion and resting-state functional connectivity (rsFC) estimates. As in the existing literature, we expect motion-correlated artifacts to take two forms: global effects and spatially specific or distance-dependent effects.

Fluctuations in time series data

Power and colleagues (2012, 2014) demonstrated that modest movements of the head are associated with large blood oxygen level-dependent (BOLD) signal changes across gray matter, white matter, and cerebrospinal fluid voxels. Motion-correlated fluctuations in BOLD signal appear quite complex: they may increase, decrease, or both before returning to baseline; they may fluctuate across the whole brain (i.e., globally distributed), in some regions more than others (i.e., spatially specific); and they may be brief (Satterthwaite et al., 2013) or temporally extended (10 or more seconds after motion ends) (Power et al., 2014). Some of these BOLD fluctuations may directly result from motion, and others may be motion-correlated artifacts due to physiological processes time-locked to motion (e.g., yawning not only moves the head but also changes heart rate and pCO₂ concentration). Either way, the concern is that the motion-correlated artifact influences rfMRI data for all individuals on average, but more for high-motion than low-motion individuals.

Distance-dependent artifact

Artifactual variance during motion tends to be more similar for nearby voxels than distant voxels (see discussion in Power et al., 2015). This results in correlations between head motion and rsFC estimates that are higher for short-distance connections and lower for long-distance connections (Satterthwaite et al., 2012). Censoring (i.e., removing) high-motion time points reduces motion-correlated artifact and usually decreases correlations between nearby parcels and increases correlations between more distant parcels (Power et al., 2012, 2014).

Motion-group differences

Head motion during rfMRI scans varies across individuals and may be confounded with factors of interest, such as age (Power et al., 2012), attention-deficit/hyperactivity disorder and impulsivity (Epstein et al., 2007; Kong et al., 2014), and bipolar disorder and schizophrenia (Mamah et al., 2013). Unfortunately, rsFC estimates also appear biased in higher motion versus lower motion individuals and groups (Van Dijk et al., 2012). Differences between groups that vary in head motion are apparent across the frequency spectrum (Satterthwaite et al., 2013) and across analysis strategies (Power et al., 2015; Satterthwaite et al., 2012; Yan et al., 2013b).

Quality control-rsFC plots

The relationship between individual differences in head motion and rsFC estimates can also be interrogated with quality control (QC)-rsFC plots: plots of the correlation across participants between QC measures of head motion during the scan and rsFC estimates. Previous investigations using QC-rsFC plots (Muschelli et al., 2014; Power et al., 2014; Satterthwaite et al., 2013) show higher rsFC estimates in individuals with greater head motion. These increased rsFC estimates exist across the whole brain (i.e., global), but to a greater extent for short-distance connections (i.e., spatially specific). Although methods such as aCompCor (Muschelli et al., 2014) and censoring (Power et al., 2014) reduce spatially specific effects, they fail to substantially reduce global effects.

Goals of the current study

In this study, we investigate the relationship between estimated head motion and resting-state correlations in HCP data. Using the procedures of Power and colleagues (2014), we investigate BOLD fluctuations and distance-dependent changes in correlations related to high-motion time points, and report how group and individual differences in head motion relate to differences in rsFC estimates. We investigated the efficacy of several denoising techniques, including motion regression, censoring of high-motion time points, FMRIB independent component analysis (ICA)-based X-noiseifier (FIX) denoising, and mean grayordinate time series regression (MGTR) (see the Mean grayordinate time series regression section for more information).

Materials and Methods

Participants

This investigation evaluated rfMRI data from the HCP 500 Subjects Public Data Release. HCP participants were between 22 and 35 years of age at the time of recruitment and did not have a documented history of psychiatric, neurological, or medical disorders known to influence brain function. For a more detailed description of inclusion and exclusion criteria for HCP, see Van Essen and colleagues (2013).

The sample of participants included 183 participants (mean age = 29.11, standard deviation [SD] = 3.55) divided into three motion groups with 26 men and 35 women in each. (Information about demographics and motion is presented in Table 1. See the Selecting Motion Groups section in the Supplementary Data for more details; Supplementary Data are available online at www.liebertpub.com/brain). The three motion groups did not differ with respect to age [F(2,180) = 0.37, p = 0.693], years of education [F(2,178) = 1.35, p = 0.263], race [χ ²(8) = 3.79, p = 0.878], or ethnicity [χ ²(2) = 2.26, p = 0.323]. By design, the motion groups did differ with respect to the proportion of high-motion time points in their scans [F(2,180) = 232.79, p = 1.20 × 10⁻⁵⁰] and mean framewise displacement (FD) [F(2,180) = 66.17, p = 3.06 × 10⁻²²].

Table 1.

Demographic and Head Motion Information Broken Down by Motion Group

	Motion groups
	Low-motion	Medium-motion	High-motion
Age
Mean	29.43	28.92	28.98
SD	3.62	3.74	3.31
Min.	22	22	23
Max.	36	36	35
Gender
Men	26	26	26
Women	35	35	35
Race
Asian/Hawaiian/Pacific Is.	1	2	2
Black or African American	14	12	15
More than one	0	1	1
Unknown or not reported	2	1	0
White	44	45	43
Ethnicity
Hispanic/Latino	8	9	4
Not Hispanic/Latino	53	52	57
Education
Mean	15.082	14.750	14.550
SD	1.9519	1.7527	1.6917
Min.	11	11	12
Max.	17	17	17
Proportion of time points censored
Mean	0.1602	0.2333	0.3384
SD	0.0274	0.0210	0.0715
Min.	0.0792	0.2000	0.2658
Max.	0.1992	0.2650	0.6517
Mean FD
Mean	0.1187	0.1460	0.2157
SD	0.0217	0.0340	0.0727
Min.	0.0800	0.0800	0.1000
Max.	0.1800	0.2100	0.5000

Head motion, measured by mean FD or proportion of time points censored, differs between groups. However, there are no statistically significant differences among groups in age, gender, race, or ethnicity.

FD, framewise displacement.

Image acquisition

Details of the MRI acquisition parameters for the HCP were described elsewhere (Ugurbil et al., 2013). Structural T1-weighted and T2-weighted images were collected at 0.7 mm isotropic resolution. Whole-brain EPI acquisitions were acquired on the 3T Siemens Connectome scanner: 32-channel head coil, TR = 720 msec, TE = 33.1 msec, in-plane FOV = 208 × 180 mm, 72 slices, 2.0 mm isotropic voxels, and multiband acceleration factor of 8 (Feinberg et al., 2010).

Overview of four types of preprocessed rfMRI data

We evaluated the reduction in motion-correlated artifacts after FIX denoising (e.g., removing variance classified as noise by FIX) versus without FIX denoising, and after MGTR versus without MGTR. Crossing these two factors of interest yielded four types of preprocessed rfMRI data: FIX, MPP, FIX+MGTR, and MPP+MGTR.

Preprocessing of HCP rfMRI data

The HCP FIX-denoising pipeline uses a gentle high-pass temporal filter (using fslmaths with 2000 sec cutoff), motion regression (i.e., regression of 24 movement parameters: six rigid -body motion parameters, their backward temporal derivatives, and squares of those 12 time series), and applies a nonaggressive regression based on ICA to remove variance in noise components that was orthogonal to signal components (Salimi-Khorshidi et al., 2014).

The inputs to our preprocessing stream were both minimally preprocessed (MPP) and FIX-denoised rfMRI data from the HCP 500 subject release. To ensure that comparisons between MPP and FIX-denoised data primarily reflected differences due to nonaggressive regression of ICA noise components, we preprocessed the MPP data using procedures similar to the FIX-denoising pipeline, including the lenient high-pass temporal filter (2000 sec cutoff) and motion regression, but excluding the regression of ICA noise component variance. Additional high-pass filtering (0.009 Hz) was conducted after regressing these confound time series (Carp, 2013). (Additional details about the HCP Minimal Preprocessing Procedures and HCP FIX-Denoising Procedures are provided in Supplementary Data.)

Mean grayordinate time series regression

Our analyses utilized the connectivity informatics technology initiative (CIFTI) dense time series files, which represent the rfMRI time series for 91282 grayordinates (i.e., surface-based vertices and subcortical voxels, constrained to gray matter) [see Glasser et al. (2013) for more details about CIFTI format]. Therefore, rather than global signal regression, we performed MGTR, computing the mean grayordinate time series from the MPP or FIX data, then regressing it and its backward derivative from each grayordinate. Although the mean grayordinate time series reflects gray matter only, MGTR seems a reasonable replacement for global signal regression due to the strong average correlation between the global signal and mean grayordinate time series (r = 0.94).

Using grayordinate plots to visualize fluctuations in rfMRI data

Residual grayordinate plots (shown in grayscale) display the time series after denoising. These reflect the nature of rfMRI data going into rsFC estimation after each denoising strategy. Difference grayordinate plots (shown in color) reflect BOLD fluctuations removed by each denoising strategy, computed as the difference between the current stage and a specified previous stage.

Both types of grayordinate plots display the time series data across time points (columns) and grayordinates (rows). Intensities are displayed as z-scores, standardized relative to the mean and SD for that grayordinate. White or black values in the residual grayordinate plots and red or blue values in the difference grayordinate plots reflect time points when the BOLD signal for that grayordinate is relatively extreme (greater than 2 SD from the mean). These values might influence correlation values strongly.

Censoring and its role in identifying motion-correlated artifact

In analyses involving censoring, we explicitly deleted high-motion time points from the time series before analysis. High-motion time points were identified using a combination of FD and DVARS (temporal Derivative, then RMS VARiance over elementS) thresholds (as defined in the Defining High-Motion Time Points: Framewise Displacement and DVARS section in Supplementary Data), were defined from the original MPP data, and were identical for each denoising strategy.

We censored high-motion time points not only as a strategy to reduce motion-correlated artifacts but also as a means to index any remaining artifact still correlated with head motion. Previous evaluations showed changes in rsFC estimates after censoring high-motion time points (Power et al., 2012, 2014; Yan et al., 2013a, 2013b). If motion-correlated artifact is present, high-motion time points will distort rsFC estimates. If denoising removes motion-correlated artifacts entirely, then rsFC estimates involving high-motion time points will no longer be distorted, and will not differ from rsFC estimates involving only low-motion time points. Therefore, comparing rsFC estimates for censored data (i.e., low-motion time points only) with rsFC estimates for uncensored data (i.e., low-motion and high-motion time points) will identify motion-correlated artifacts remaining after each denoising strategy.

Parcellated connectomes

Functional connectivity was evaluated in subsequent analyses using full correlations (Pearson's correlation coefficients) between parcel time series extracted from CIFTI grayordinates. We utilized 333 cortical parcels from the Gordon and colleagues (2014) parcellation because they have higher functional homogeneity than several other published parcellation schemes. We also added 19 subcortical parcels from the group-average CIFTI atlas for a total of 352 parcels. (Additional details about Parcellated Connectomes are included in Supplementary Data.)

ΔR plots display motion-correlated artifacts from high-motion time points

The influence of residual motion-correlated artifacts may be revealed with ΔR plots—that is, the difference between rsFC estimates derived from censored versus uncensored data. For each connection, rsFC estimates are computed using the censored time series and the uncensored time series. rsFC estimates are Fisher z transformed, averaged across participants, and converted back to Pearson's r. The difference between rsFC estimates (averaged across participants) from censored and uncensored data is then plotted as a function of the distance between parcels.

The slope of these ΔR plots may demonstrate distance-dependent effects, indicating that the residual motion-correlated artifact differs for short-distance versus long-distance connections. A shift in the mean of the ΔR plots may indicate global influences on rsFC estimates such that high-motion time points increase correlations across all connections regardless of distance.

The mean ΔR and the linear relationship of ΔR with distance from a general linear model (GLM) will estimate the global and distance-dependent effects, respectively. However, because each parcel contributes to multiple rsFC observations, rsFC observations are nonindependent. Statistical inference on the mean and slope is biased by nonindependent observations. However, the parameter estimates themselves are unbiased. After accounting or controlling for family structure, participants provide independent observations to conduct valid statistical inference.

Therefore, we fit a GLM for each participant, predicting ΔR with an intercept term and a distance term (i.e., Euclidean distance between parcels after mean centering). This yielded independent observations of the mean and slope. Then, one-sample t-tests determined whether those parameters were different from zero in the high-motion group (Table 2) and in the low-motion group (Supplementary Table S1).

Table 2.

Means and Slopes of the ΔR Plots as Estimated Using a General Linear Model

	MPP		FIX		Change due to FIX
Denoising strategy	Mean	Slope	Mean	Slope	Mean	Slope
Before MGTR	−6.659 × 10⁻⁴	5.677 × 10⁻³ ^***	−4.204 × 10⁻³ ^***	3.450 × 10⁻³ ^***	−3.538 × 10⁻³	−2.226 × 10⁻³
After MGTR	−8.820 × 10⁻⁴ ^***	5.683 × 10⁻³ ^***	−2.086 × 10⁻⁴ ^*	2.896 × 10⁻³ ^***	6.734 × 10⁻⁴ ^**	−2.787 × 10⁻³ ^**
Change due to MGTR	−2.161 × 10⁻⁴	6.694 × 10⁻⁶	3.995 × 10⁻³ ^***	−5.538 × 10⁻⁴

Differences between denoising strategies were determined by computing the difference in ΔR values for those strategies. Slopes are expressed as change in ΔR per 100 mm, roughly reflecting the difference between average short-distance and long-distance connections.

***

p < .0001, ^** p < 0.01, ^* p < 0.05.

FIX, FMRIB independent component analysis-based X-noiseifier; MGTR, mean grayordinate time series regression; MPP, minimally preprocessed.

To estimate differences between two denoising strategies, we computed the difference in ΔR values between two strategies, then predicted those values using a GLM for each participant with intercept and distance terms. The mean and slope reflected the difference in global and distance-dependent effects between the two strategies. Finally, one-sample t-tests across high-motion (Table 2) and low-motion (Supplementary Table S1) participants indicated whether the means and slopes differed from zero.

QC-rsFC plots reveal residual distance-dependent and global artifacts after censoring

For each pairwise connection, correlations were computed across 183 participants between the value of the rsFC estimate (after censoring) and the proportion of time points censored under the combined FD and DVARS criteria. These QC-rsFC correlations were plotted as a function of distance between the regions to examine whether distance-dependent artifacts might be present in the rfMRI time series. At the same time, QC-rsFC correlations that are elevated across all distances are suggestive of global artifacts remaining in the time series.

To determine the statistical significance of global and distance-dependent effects in the QC-rsFC plots, we used a similar approach to significance testing of the ΔR plots. First, for each participant, a GLM predicted rsFC estimates with intercept and distance terms, providing independent estimates of the global (mean) and distance-dependent (slope) effects. Second, to determine if individual differences in head motion modulated global and distance-dependent effects, group-level GLMs predicted the subject-level parameters with a QC term (i.e., head motion estimated by proportion of time points censored after mean centering) and an intercept. We estimated the significance of the global QC-rsFC relationship by QC predicting the mean and the significance of the distance-dependent QC-rsFC relationship by QC predicting the slope. We report these results separately for censored data (Table 3) and uncensored data (Supplementary Table S2).

Table 3.

Relationship of Head Motion (QC) with the Mean and Slope of the Resting-State Functional Connectivity–Distance Relationship

	MPP		FIX		Change due to FIX
Denoising strategy	Mean	Slope	Mean	Slope	Mean	Slope
Before MGTR	0.3272^***	−0.0549^*	0.1940^***	−0.0450^**	−0.1332^**	0.0099
After MGTR	0.0296^***	−0.0731^***	0.0084^*	−0.0246	−0.0212^***	0.0485^**
Change due to MGTR	−0.2975^***	−0.0182	−0.1856^***	0.0205^*

See the Materials and Methods section for more details regarding how these relationships were computed. Slopes are expressed as change in QC-rsFC relationship per 100 mm, roughly reflecting the difference between average short-distance and long-distance connections.

***

p < 0.001, ^** p < 0.01, ^* p < 0.05.

QC, quality control; rsFC, resting-state functional connectivity.

To estimate differences between denoising strategies, we computed the difference in rsFC values between two strategies, then predicted those values for each participant using a GLM with intercept and distance terms. The mean and slope reflect the difference in global and distance-dependent effects between those two strategies. Then, we estimated whether those parameters related to the degree of head motion (QC) across participants.

Computing motion-group differences

We computed the percentage of significant motion-group differences observed in the correlation matrices after each denoising strategy. For a given rsFC estimate, unpaired t-tests determined whether motion-group differences (e.g., between high- and low-motion groups) were significant (see the Establishing Alpha Level for Testing Motion-Group Differences section in the Supplementary Data for more details). We reported the percentage of significant motion-group differences across all connections, and separately for short-distance (4.8–58.6 mm), medium-distance (58.6–112.3 mm), and long-distance (112.3–166.1 mm) connections.

For each denoising strategy, we used permutation testing to determine whether the observed number of motion-group differences was significantly greater than expected by chance. The null distribution was estimated by permuting motion-group labels across participants 10,000 times.

Results

Grayordinate plots display spatially specific and global fluctuations

Figure 1 shows grayordinate plots from HCP participant 107422 from time points 600–1000 (4.8 min) of the rfMRI_REST1_RL scan. This participant had below average head motion, but features evident in these data are present across the full sample of HCP participants.

FIG. 1.

Nature of blood oxygen level-dependent fluctuations in HCP data, and aspects removed by each denoising stage: Five residual grayordinate plots (in grayscale on left) show rfMRI data after each denoising stage: (A) MPP before motion regression, (B) MPP after motion regression, (C) FIX, (D) MPP+MGTR, and (E) FIX+MGTR. Four difference grayordinate plots (in color on right) show variance removed by specific denoising steps (estimated by subtracting the current preprocessing stage from a specified prior stage): (F) motion regression (A, B); (G) ICA denoising (B, C); (H) MGTR (B–D); and (J) ICA denoising plus MGTR (B–E). For each grayordinate plot, columns reflect time points and rows reflect grayordinates. Intensities are z-scored (across time, separately for each vertex) and range from −2 to +2. The top panel on both sides shows FD (red), with horizontal lines marking FD = 0.2 mm (suggested as a censoring threshold by Power et al., 2014) and FD = 0.39 mm (current study threshold for FD). MGT (black lines) and DVARS (blue lines) are derived from data after each denoising strategy. The horizontal line in (A) corresponds to the DVARS censoring threshold of 4.9 arbitrary MR units. Green arrows indicate time periods displaying the global artifact, which manifests as similar effects across space and occurs across most grayordinates. Spatially specific artifacts, which manifest as dissimilar effects across space, are indicated by red arrows (instances that occur at few grayordinates) and blue arrows (instances that occur across most grayordinates). FD, framewise displacement; FIX, FMRIB ICA-based X-noiseifier; HCP, Human Connectome Project; ICA, independent component analysis; MGTR, mean grayordinate time series regression; MPP, minimally preprocessed; rfMRI, resting-state functional magnetic resonance imaging.

Global and spatially specific artifact is evident in HCP MPP rfMRI data before motion regression. The residual grayordinate plot (Fig. 1A) prominently displays global fluctuations as vertical bands that show similar sign (i.e., positive or negative) and magnitude across most grayordinates. Some example time periods are indicated by green arrows. Global fluctuations vary in duration, with some lasting 20 sec or more. Spatially specific fluctuations appear as horizontal bands that varied in sign and magnitude across grayordinates. Some examples are indicated by red and blue arrows.

Spatially specific artifact is reduced by motion regression and FIX denoising

In some cases, spatially specific noise was isolated to a small proportion of grayordinates, but extended over several time points. These manifest as rows of higher or lower intensity (red arrows in Fig. 1A). These are reduced by motion regression (red arrow in Fig. 1F) and FIX denoising (red arrows in Fig. 1G), resulting in reduced noise in the residual grayordinate plots (red arrows in Fig. 1B, C). In other cases, spatially specific artifact was evident at specific time points, but showed effects that varied in sign and magnitude across grayordinates. This spatially specific artifact is evident in the variance removed by motion regression (Fig. 1F, blue arrows) and FIX denoising (Fig. 1G, blue arrows). Despite the artifact being removed by these strategies, the residual grayordinate plots show clear global fluctuations left behind by both motion regression (Fig. 1B) and FIX denoising (Fig. 1C).

Global fluctuations are slightly reduced by FIX denoising, but dramatically reduced by MGTR

FIX denoising appeared to target some time points that exhibit global fluctuations (Fig. 1G, green arrows). Some of these global fluctuations occurred after spikes in DVARS or FD, but others did not. FIX denoising results in a reduction in the magnitude of global fluctuations, but FIX rfMRI data still show substantial global fluctuations (Fig. 1C, green arrows) in the grayordinate time series. However, MGTR dramatically reduced or eliminated the global fluctuations (Fig. 1D, H).

Combining FIX denoising and MGTR reduced both spatially specific and global fluctuations

The above results suggest that FIX denoising and MGTR have different types of effects on motion-correlated fluctuations in the rfMRI data. The differences are most easily observed in the difference grayordinate plots, where it appears that combining FIX denoising and MGTR controls both spatially specific fluctuations and global fluctuations in those data (Fig. 1J).

ΔR plots reveal distance-dependent and global artifacts

Within the high-motion group, ΔR plots show both distance-dependent (spatially specific) and global artifacts. Censoring high-motion time points elucidates distance-dependent artifact present in MPP rfMRI data (Fig. 2A), confirmed by a statistically significant slope as a function of distance (Table 2). In the FIX rfMRI data (Fig. 2B), censoring appears to reduce both global and distance-dependent artifacts, as indicated by a statistically significant mean and slope related to distance (Table 2). The reduction in distance-dependent artifact due to censoring is visually evident in ΔR plots from MPP+MGTR and FIX+MGTR data (Fig. 2C, D). However, both the mean and slope of the ΔR effect were statistically different from zero in MPP+MGTR and FIX+MGTR data (Table 2).

FIG. 2.

Censoring high-motion time points reveals spatially specific and global shift artifacts in ΔR plots: Red cloud (and white loess fit line) shows effects of censoring high-motion time points on rsFC estimates in the high-motion group, plotted as the function of distance between parcels being correlated. Black cloud (and gray loess fit) shows positive control (censoring equal number of randomized time points). Range of ΔR (y-axis) from 0.1 to −0.1, following Power and associates (2014). Panels show effects of censoring on average rsFC estimates from high-motion group for (A) MPP, (B) FIX, (C) MPP+MGTR, and (D) FIX+MGTR time series data. Analogous plots for the low-motion group are in Supplementary Figure S1. rsFC, resting-state functional connectivity.

We evaluated differences in ΔR between denoising strategies to determine whether those denoising strategies might reduce distance-dependent and global artifacts. Neither the difference between FIX and MPP rfMRI data nor the difference between MPP+MGTR and MPP rfMRI data showed significant distance-dependent or global ΔR effects (Table 2). One might assume that this indicates that neither denoising strategy removes motion-correlated artifacts. However, the high variability of the MPP rfMRI data might have reduced our sensitivity to the effects of FIX and MGTR. Consistent with this hypothesis, the difference between FIX+MGTR and FIX rfMRI data is significant for the mean term, and the difference between FIX+MGTR and MPP+MGTR is significant for both the mean and slope. These findings are consistent with our interpretation of grayordinate plots that MGTR reduces global fluctuations and that FIX reduces spatially specific fluctuations as well as global fluctuations. The results were similar for low-motion participants (Supplementary Fig. S1 and Supplementary Table S1).

QC-rsFC plots reveal residual distance-dependent and global artifacts after censoring

We investigated the relationship between rsFC estimates and individual differences in estimated head motion (Fig. 3 and Table 3). The MPP rfMRI data show both global and distance-dependent relationships between rsFC estimates and the amount of motion during the scan, even after censoring, supported by the fact that head motion is significantly related to both the mean and slope as a function of distance. FIX denoising resulted in a statistically significant reduction in the relationship of head motion with the mean, but head motion still modulated both the mean and slope in FIX data. MPP+MGTR data showed a larger reduction in the relationship between head motion and the mean compared with MPP data, but again head motion still modulated both the mean and slope in MPP+MGTR data. For FIX+MGTR data (after censoring), the global relationship with head motion was still statistically significant, although at its lowest level, but the distance-dependent relationship with head motion was eliminated. Global and distance-dependent effects were significantly reduced in FIX+MGTR data relative to both FIX data and MPP+MGTR data. These effects were similar for uncensored data (Supplementary Fig. S2 and Supplementary Table S2).

FIG. 3.

QC-rsFC plots show the correlation across participants between the rsFC estimates after censoring and degree of head motion (quantified by proportion of time points censored using the combined FD and DVARS criteria). The QC-rsFC relationship is plotted for each of the 61,776 connections as a function of the distance between parcels for (A) MPP, (B) FIX, (C) MPP+MGTR, and (D) FIX+MGTR time series data. Analogous plots for uncensored data are in Supplementary Figure S2. QC, quality control.

Substantial differences exist between motion groups and are primarily reduced by MGTR

Consistent with prior work (Power et al., 2014), the observed number of differences between the low- and high-motion groups was inflated strongly above chance in the MPP rfMRI data (Table 4). Both censoring high-motion time points and FIX denoising reduced the number of motion-group differences slightly, but the number of motion-group differences was still substantially elevated. In contrast, MGTR dramatically reduced the number of differences between high- and low-motion groups, although not to chance levels, for both MPP+MGTR and FIX+MGTR rfMRI data.

Table 4.

Percentage of Connections Showing Motion-Group Differences at α = 0.000242

		Uncensored (%)			Censored (%)
Denoising strategy	Distance bin	Low vs. High	Med. vs. High	Low vs. Med.	Low vs. High	Med. vs. High	Low vs. Med.
MPP	All	22.02^***	0.13^**	0.04^*	12.27^***	0.11^*	0.04^*
	Short	23.91^***	0.23^**	0.05	13.90^***	0.17^*	0.04
	Medium	21.35^***	0.09^**	0.04	11.55^***	0.08^*	0.05
	Long	21.23^***	0.09	0.04	12.21^***	0.11	0.03
FIX	All	14.23^***	0.19^*	0.04	12.07^***	0.12^*	0.08^*
	Short	16.77^***	0.24^**	0.04	14.00^***	0.11^*	0.06
	Medium	14.20^***	0.19^*	0.03	12.03^***	0.13^*	0.09^*
	Long	8.64^***	0.14	0.05	7.98^***	0.08	0.08
MPP+MGTR	All	0.44^***	0.03^**	0.01	0.20^***	0.02	0.03
	Short	0.90^***	0.04	0.01	0.28^**	0.02	0.03
	Medium	0.49^***	0.03	0.01	0.16^**	0.01	0.03
	Long	0.51^**	0.01	0.01	0.26	0.01	0.03
FIX+MGTR	All	0.38^***	0.05^*	0.02	0.28^***	0.04	0.02
	Short	0.40^***	0.03	0.01	0.27^**	0.03	0.01
	Medium	0.36^***	0.06^*	0.03	0.28^***	0.05^*	0.03
	Long	0.46^***	0.04	0.01	0.34^***	0.00	0.01

The percentage is reported for all connections, and separately for short-, medium-, and long-distance connections. The statistical significance is determined through permutation testing; significance level is indicated by asterisks.

***

p < 0.001, ^** p < 0.01, ^* p < 0.05.

Discussion

The current study evaluated the presence of the motion-correlated artifact in HCP rfMRI data and its removal by some common denoising strategies. Consistent with prior work (Power et al., 2014; 2012; Satterthwaite et al., 2013; Van Dijk et al., 2012), the current results suggest two separate classes of motion-correlated artifacts: global and spatially specific. The denoising strategies we tested differed in efficacy for these two classes. Consequently, a combination of methods—such as FIX denoising to reduce the spatially specific artifact and MGTR to reduce the global artifact—will be necessary to effectively address motion-corrected artifacts.

Denoising strategies had differential effects on motion-correlated artifacts

Motion regression, censoring, and FIX denoising seemed to have the largest effect on the spatially specific artifact. The grayordinate plots revealed that motion regression strongly affected some grayordinates, but did not affect others at all, and the sign of the effect varied across grayordinates. Censoring high-motion time points demonstrated a distance-dependent reduction in the ΔR plot for MPP rfMRI data. FIX reduced the distance-dependent artifact, as indexed by the comparison of FIX+MGTR versus MPP+MGTR data in the ΔR and QC-rsFC plots.

On the other hand, evidence of global artifact was only reduced substantially by MGTR. Global artifacts manifested as vertical bands in residual grayordinate plots, influences of high-motion time points across all distances in ΔR plots, and relationships between motion estimates and rsFC estimates across the entire brain in QC-rsFC plots. These global effects were significantly reduced for FIX+MGTR data compared with FIX data.

Motion-group differences were reduced by FIX and by censoring. However, these reductions were relatively modest. They did not appear more strongly for short-distance connections, as might be expected from a reduction in spatially specific artifacts. In contrast, MGTR reduced motion-group differences substantially across all distance bins. This pattern suggests that the bulk of motion-group differences resulted from globally distributed artifact, as opposed to spatially specific artifact.

FIX denoising seems to primarily address spatially specific artifacts while doing less to address global artifacts. Perhaps that should not be surprising. First, as with many other ICA-based denoising methods, FIX denoising uses spatial ICA algorithms that maximize the spatial independence of components. Consequently, global noise variance is less likely to be isolated into a separate component and removed by ICA denoising. Temporal ICA algorithms may be more likely to identify global components (Smith et al., 2012), which subsequently might be classified as noise. Second, the HCP FIX-ICA-denoising pipeline applies nonaggressive denoising, regressing only the portion of noise variance orthogonal to signal components. Nonaggressive denoising may not remove global noise that is shared across signal and noise components (Smith et al., 2013).

Censoring appears to reduce spatially specific artifacts. However, the significant slope of the QC-rsFC relationship for MPP and MPP+MGTR data suggests that censoring did not eliminate distance-dependent artifacts as reported in previous studies. We postulate that the reduced efficacy of the censoring procedure in HCP rfMRI data might result, in part, from increased noise variance in the FD and DVARS motion estimates, as discussed in the Supplementary Data.

Potential influences of physiological noise in HCP rfMRI data

The rfMRI time series from many HCP participants contain periodic fluctuations in signal intensity across the brain. These are visible as evenly spaced global bands in the grayordinate plots, resulting in periodic fluctuations in the mean grayordinate time series (cf., center of Fig. 1A). One intriguing question is whether some portion of this motion-correlated global artifact might be physiological in origin. Head motion and physiological artifact may be coupled if respiratory movements directly cause head motion (e.g., yawning, sneezing, deep breathing, or sighing), if ballistocardiographic forces directly cause head motion, or if effort exerted while moving the head and body in the scanner results in breath holding or change in heart rate.

Interestingly, the global artifact in HCP rfMRI data (Fig. 1A) and other rfMRI data (Power et al., 2014) seems to lag head motion by 10–20 sec. It has also been shown that the respiratory artifact has a response function that extends across a similar time frame (Birn et al., 2008; Chang and Glover, 2009). Respiration modulates BOLD signal across the entire brain, but that modulation is stronger in somatosensory, motor, and visual cortices than other brain regions (Birn et al., 2006; Wise et al., 2004). Compared with other areas, these regions show stronger correlations both with global signal (Fox et al., 2009) and with individual differences in head motion (Pujol et al., 2014; Yan et al., 2013a). These patterns could indicate that changes in respiration are frequently accompanied by head motion and that these respiration changes may ultimately be a core driver of global artifact in the BOLD signal. This hypothesis should be tested directly in future work.

Disadvantages of unmitigated global artifact

Concerns have been raised regarding global signal regression (e.g., Murphy et al., 2009; Saad et al., 2012; Schölvinck et al., 2010), which may also apply to MGTR. However, we argue that a major underaddressed concern in the literature is that unmitigated motion-correlated artifacts can be mistakenly attributed to meaningful individual or group differences. Researchers are rightly concerned about reducing sensitivity to individual and group differences by inadvertently discarding a signal. Nonetheless, we argue that it is critical to fully address motion-correlated artifacts when researching individual and group differences that covary with head motion.

If head motion is correlated with one's variable of interest, any motion-correlated artifact that is retained can be mistakenly attributed to that variable of interest. Denoising strategies that are too lenient could result in features similar to motion-correlated artifacts—such as distance-dependent correlations and globally increased correlations—being attributed to a wide variety of clinical, developmental, and psychological group differences. Unfortunately, the publication bias in scientific literature (Franco et al., 2014) may make it relatively easy to report significant correlations resulting from motion-correlated artifacts, but harder to refute such findings with correlations that are not significant after denoising. Therefore, we lean toward being more conservative with denoising to avoid contaminating the literature with artifactual findings.

In our opinion, time series denoising methods such as those investigated here are strongly preferred to leaving behind unmitigated motion artifacts in rfMRI data. In the current study, the denoising methods were chosen to target and remove specific aspects of time series variance that have been linked to artifacts previously in the literature (Birn et al., 2006; Friston et al., 1996; Wise et al., 2004) and empirically in the data (Power et al., 2014; Salimi-Khorshidi et al., 2014). Importantly, if rsFC estimates no longer correlate with individual and group differences after these time series denoising methods, it indicates that those variables only related to aspects of the time series that were removed by the denoising methods. If individual or group differences resulted entirely from aspects of the time series previously linked to artifacts in the literature, it leads to substantial doubt that those differences might be neural in origin.

To adequately address global artifacts while avoiding the drawbacks associated with global signal regression or MGTR, additional research is needed to identify denoising strategies that eliminate global artifacts related to head motion and physiological processes. For example, CompCor (Behzadi et al., 2007; Chai et al., 2012) uses white matter and cerebrospinal fluid signal as confound regressors to address motion and physiological artifacts without removing neural signal or inducing anticorrelations. However, Muschelli and colleagues (2014) (Supplementary Fig. S1) found increased correlations between FD and rsFC estimates across the brain after using CompCor, suggesting that the residual global artifact confounded with individual differences in head motion. Numerous other approaches (c.f., Yan et al., 2013c, for example) exist, where participant-level motion estimates are treated as covariates in group-level and individual difference analyses. These approaches may reduce the influence of global artifact on those analyses. However, because these approaches regress out variance from the rsFC estimates rather than the time series, they may be less capable of separating motion-correlated artifacts from meaningful individual differences confounded with motion.

After denoising HCP data with FIX+MGTR, we still observed a greater number of motion-group differences than expected by chance. It is possible that these motion-group differences reflect differences between the high- and low-motion participants in neural activity related to factors that might influence head motion, such as alertness, anxiety, and cognition. We believe that a comparison of within-participant motion effects versus between-participant motion effects may help to resolve this question. However, it is critical that future studies carefully equate the degree of head motion in within-participant with between-participant analyses to which they are being compared.

Conclusion

rfMRI data with high spatial and temporal resolution are available to the public from The HCP (http://humanconnectome.org). Although cutting-edge technological advances have led to many improvements, HCP rfMRI data still are affected by artifacts correlated with head motion and other physiological effects. An evaluation of artifactual changes in BOLD signal intensity suggests that the presence of spatially specific and global artifacts correlated with head motion. Several denoising techniques—including FIX-ICA denoising, motion regression, and censoring high-motion time points—primarily address spatially specific artifacts. However, these strategies leave substantial differences between rsFC estimates from high- and low-motion individuals. In contrast, MGTR primarily addresses global artifacts, and substantially reduces rsFC differences between high- and low-motion individuals. Consequently, a combination of denoising strategies that captures both spatially specific and global aspects of motion-correlated artifacts will be necessary for productive analysis of HCP rfMRI data.

Footnotes

Acknowledgments

Data were provided by the HCP, WU-Minn Consortium (Principal Investigators: David Van Essen and Kamil Ugurbil; 1 U54MH091657) funded by the 16 NIH Institutes and Centers that support the NIH Blueprint for Neuroscience Research; and by the McDonnell Center for Systems Neuroscience at Washington University. The authors would like to thank A.M. Winkler, T.E. Nichols, and M.J. Strube for suggestions regarding valid methods for statistical inference and Matt Glasser and David Van Essen for important perspectives and thoughts on this work. While the authors do not always agree on approaches and interpretations, their insights have provided valuable feedback on the work presented here. Finally, the authors thank the anonymous reviewers for their insightful comments and suggestions on an earlier version of the manuscript.

Author Disclosure Statement

No competing financial interests exist.

References

Behzadi

, Restom

, Liau

, Liu

. 2007. A component based noise correction method (CompCor) for BOLD and perfusion based fMRI. NeuroImage, 37:90–101.

Birn

, Diamond

, Smith

, Bandettini

. 2006. Separating respiratory-variation-related fluctuations from neuronal-activity-related fluctuations in fMRI. NeuroImage, 31:1536–1548.

Birn

, Smith

, Jones

, Bandettini

. 2008. The respiration response function: The temporal dynamics of fMRI signal fluctuations related to changes in respiration. NeuroImage, 40:644–654.

Carp

. 2013. Optimizing the order of operations for movement scrubbing: Comment on Power et al. NeuroImage, 76:436–438.

Chai

, Castañón

, Öngür

, Whitfield-Gabrieli

. 2012. Anticorrelations in resting state networks without global signal regression. NeuroImage, 59:1420–1428.

Chang

, Glover

. 2009. Relationship between respiration, end-tidal CO2, and BOLD signals in resting-state fMRI. NeuroImage, 47:1381–1393.

Epstein

, Casey

, Tonev

, Davidson

, Reiss

, Garrett

, et al. 2007. Assessment and prevention of head motion during imaging of patients with attention deficit hyperactivity disorder. Psychiatry Res, 155:75–82.

Feinberg

, Moeller

, Smith

, Auerbach

, Ramanna

, Glasser

, et al. 2010. Multiplexed echo planar imaging for sub-second whole brain FMRI and fast diffusion imaging. PLoS One, 5:e15710.

Fox

, Zhang

, Snyder

, Raichle

. 2009. The global signal and observed anticorrelated resting state brain networks. J Neurophysiol, 101:3270–3283.

10.

Franco

, Malhotra

, Simonovits

. 2014. Social science. Publication bias in the social sciences: Unlocking the file drawer. Science, 345:1502–1505.

11.

Friston

, Williams

, Howard

, Frackowiak

, Turner

. 1996. Movement-related effects in fMRI time-series. Magn Reson Med, 35:346–355.

12.

Glasser

, Sotiropoulos

, Wilson

, Coalson

, Fischl

, Andersson

, et al. 2013. The minimal preprocessing pipelines for the Human Connectome Project. NeuroImage, 80:105–124.

13.

Gordon

, Laumann

, Adeyemo

, Huckins

, Kelley

, Petersen

. 2014. Generation and evaluation of a cortical area parcellation from resting-state correlations. Cereb Cortex, 26:288–303.

14.

, Gotts

, Reynolds

, Bandettini

, Martin

, Cox

, Saad

. 2013. Effective preprocessing procedures virtually eliminate distance-dependent motion artifacts in resting state fMRI. J Appl Math, 2013:1–9.

15.

Kong

X-Z

, Zhen

, Li

, Lu

H-H

, Wang

, Liu

, et al. 2014. Individual differences in impulsivity predict head motion during magnetic resonance imaging. PLoS One, 9:e104989.

16.

Mamah

, Barch

, Repovš

. 2013. Resting state functional connectivity of five neural networks in bipolar disorder and schizophrenia. J Affect Disord, 150:601–609.

17.

Murphy

, Birn

, Handwerker

, Jones

, Bandettini

. 2009. The impact of global signal regression on resting state correlations: Are anti-correlated networks introduced?. NeuroImage, 44:893–905.

18.

Muschelli

, Nebel

, Caffo

, Barber

, Pekar

, Mostofsky

. 2014. Reduction of motion-related artifacts in resting state fMRI using aCompCor. NeuroImage, 96:22–35.

19.

Power

, Barnes

, Snyder

, Schlaggar

, Petersen

. 2012. Spurious but systematic correlations in functional connectivity MRI networks arise from subject motion. NeuroImage, 59:2142–2154.

20.

Power

, Mitra

, Laumann

, Snyder

, Schlaggar

, Petersen

. 2014. Methods to detect, characterize, and remove motion artifact in resting state fMRI. NeuroImage, 84:320–341.

21.

Power

, Schlaggar

, Petersen

. 2015. Recent progress and outstanding issues in motion correction in resting state fMRI. NeuroImage, 105C:536–551.

22.

Pujol

, Macià

, Blanco-Hinojo

, Martínez-Vilavella

, Sunyer

, de la Torre

, et al. 2014. Does motion-related brain functional connectivity reflect both artifacts and genuine neural activity?. NeuroImage, 101:87–95.

23.

Saad

, Gotts

, Murphy

, Chen

, Jo

, Martin

, Cox

. 2012. Trouble at rest: How correlation patterns and group differences become distorted after global signal regression. Brain Connect, 2:25–32.

24.

Salimi-Khorshidi

, Douaud

, Beckmann

, Glasser

, Griffanti

, Smith

. 2014. Automatic denoising of functional MRI data: Combining independent component analysis and hierarchical fusion of classifiers. NeuroImage, 90:449–468.

25.

Satterthwaite

, Elliott

, Gerraty

, Ruparel

, Loughead

, Calkins

, et al. 2013. An improved framework for confound regression and filtering for control of motion artifact in the preprocessing of resting-state functional connectivity data. NeuroImage, 64:240–256.

26.

Satterthwaite

, Wolf

, Loughead

, Ruparel

, Elliott

, Hakonarson

, et al. 2012. Impact of in-scanner head motion on multiple measures of functional connectivity: Relevance for studies of neurodevelopment in youth. NeuroImage, 60:623–632.

27.

Schölvinck

, Maier

, Ye

, Duyn

, Leopold

. 2010. Neural basis of global resting-state fMRI activity. Proc Natl Acad Sci U S A, 107:10238–10243.

28.

Smith

, Beckmann

, Andersson

, Auerbach

, Bijsterbosch

, Douaud

, et al. 2013. Resting-state fMRI in the Human Connectome Project. NeuroImage, 80:144–168.

29.

Smith

, Miller

, Moeller

, Xu

, Auerbach

, Woolrich

, et al. 2012. Temporally-independent functional modes of spontaneous brain activity. Proc Natl Acad Sci U S A, 109:3131–3136.

30.

Ugurbil

, Xu

, Auerbach

, Moeller

, Vu

, Duarte-Carvajalino

, et al. 2013. Pushing spatial and temporal resolution for functional and diffusion MRI in the Human Connectome Project. NeuroImage, 80:80–104.

31.

Van Dijk

KRA

, Sabuncu

, Buckner

. 2012. The influence of head motion on intrinsic functional connectivity MRI. NeuroImage, 59:431–438.

32.

Van Essen

, Smith

, Barch

, Behrens

TEJ

, Yacoub

, Ugurbil

, Consortium

W-MH

. 2013. The WU-Minn Human Connectome Project: An overview. NeuroImage, 80:62–79.

33.

Wise

, Ide

, Poulin

, Tracey

. 2004. Resting fluctuations in arterial carbon dioxide induce significant low frequency variations in BOLD signal. NeuroImage, 21:1652–1664.

34.

Yan

C-G

, Cheung

, Kelly

, Colcombe

, Craddock

, Di Martino

, et al. 2013a. A comprehensive assessment of regional variation in the impact of head micromovements on functional connectomics. NeuroImage, 76:183–201.

35.

Yan

C-G

, Craddock

, He

, Milham

. 2013b. Addressing head motion dependencies for small-world topologies in functional connectomics. Front Hum Neurosci, 7:910.

36.

Yan

C-G

, Craddock

, Zuo

X-N

, Zang

Y-F

, Milham

. 2013c. Standardizing the intrinsic brain: Towards robust measurement of inter-individual variation in 1000 functional connectomes. NeuroImage, 80:246–262.

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.73 MB

0.00 MB