Impact of methodological variables on functional connectivity findings in autism spectrum disorders

Participants

The study included 54 ASD and 51 TD participants. Of these 27 ASD and 24 TD participants were recruited for the VS study. After exclusion for excessive head motion (>1.5 mm), 19 ASD and 19 TD participants, matched for age, handedness, sex, IQ, and head motion, were retained for the VS study. For the RSVP and RS studies, 27 ASD and 27 TD participants were scanned. For the RS study, five ASD and two TD participants were excluded for excessive motion, reducing the sample for this dataset to 21 ASD and 25 age, handedness, and IQ‐matched TD participants. An additional two participants, one from each group, were excluded from the RSVP sample due to red–green color‐blindness affecting task performance, reducing the final sample for this study to 20 ASD and 24 age, handedness, and IQ‐matched TD participants (Table 1). No significant differences or interactions were detected between groups for errors or response time on the VS task. For the RSVP task, no significant differences between groups were found for response time, but the ASD group made significantly more errors than the TD group (Table 1).

For all three studies, clinical diagnoses were confirmed using the Autism Diagnostic Interview–Revised (ADI‐R) [Rutter et al., 2003], the Autism Diagnostic Observation Schedule (ADOS) [Lord et al., 2001], and expert clinical judgment according to DSM‐IV criteria. Children with ASD‐related medical conditions (e.g., Fragile‐X syndrome, tuberous sclerosis), and other neurological conditions (e.g., epilepsy, Tourette's Syndrome) were excluded. Participants in the TD group had neither reported personal or family history of ASD nor any reported history of other neurological or psychiatric conditions. Informed assent and consent was obtained from all participants and their caregivers in accordance with the University of California, San Diego and San Diego State University Institutional Review Boards.

Data Acquisition

All MRI data were acquired on a GE 3T MR750 scanner with an eight‐channel head coil. High‐resolution structural images were acquired with a standard FSPGR T1‐weighted sequence (TR: 11.08 ms; TE: 4.3 ms; flip angle: 45°; FOV: 256 mm; 256 × 256 matrix; 180 slices; 1 mm³ resolution).

Visual search (VS) data were acquired in four functional runs, each consisting of 128 whole‐brain volumes acquired in 39 or 40 interleaved slices using a single‐shot, gradient‐recalled, echo‐planar pulse sequence (TR: 2500 ms; TE: 30 ms; flip angle: 90°; matrix: 64 × 64; slice thickness: 3.2 mm; in‐plane resolution: 3.4 mm²). Target stimulus was an upright “T” and distractors were Ts rotated 90°, 180°, and 270°. A total of 240 search trials, and 232 null trials were presented [for details see Supporting Information Figure 1 and Keehn et al., 2013b].

For rapid serial visual presentation (RSVP), four 6 min functional runs were obtained using a single‐shot gradient‐recalled, echo‐planar pulse sequence consisting of 180 whole‐brain volumes (TR: 2000 ms; TE: 30 ms; slice thickness: 3.4 mm; in‐plane resolution: 3.4 mm²). Stimuli included three streams of varying numbers presented in four iterations per trial (120 ms per iteration, 480 ms per trial; for details see Supporting Information Figure 2).

RS data were acquired in one 6:10 min scan using the same protocol as for RSVP. For this scan, participants were asked to keep their eyes fixated on a white cross‐hair in the center of a black rear‐projection screen, to relax, and to try not to fall asleep.

Preprocessing and Activation Analysis

All data were preprocessed and analyzed using the Analysis of Functional Neuroimages suite [AFNI; Cox, 1996]. The first four time points for each run were discarded to remove effects of signal instability, and slice‐time correction was performed. Functional data were coregistered to Talairach space, resampled to isotropic 3 mm voxels, and spatially smoothed to a global full‐width at half‐maximum (FWHM) of 5 mm, using AFNI's 3dBlurToFWHM.

VS and RSVP fMRI data were then scaled to a grand signal mean of 100, and concatenated to create a single time‐series. The hemodynamic impulse response function (IRF) for each stimulus type was estimated using a general linear model. Variable‐shape IRFs for each stimulus type were estimated using piecewise linear B‐spline (tent) basis functions. For RSVP, within‐ and between‐group t‐tests were used to assess activation for target‐present versus target‐absent trials. For VS, within‐ and between‐group t‐tests were used to assess activation for all search trials versus null trials. All statistical maps were corrected for multiple comparisons to a cluster corrected threshold of p < 0.05.

Seed Selection and Regions of Interest

The selection of dataset‐specific ROIs served two goals. First, it allowed us to determine whether group bias in seed ROI selection had any effect on functional connectivity findings. In simple terms, the question was whether a seed that was based on activation patterns detected in one group (rather than the other group or both groups combined) would affect the overall pattern of results. Second, selection of ROIs was required for the examination of the field of view variable, to compare effects within a given activation‐defined network (e.g., visual search network) to effects detected across the entire brain. As described above, this variable had been found crucial in our precursor study [Müller et al., 2011], as ASD underconnectivity studies were mostly those limited to effects within task‐specific networks, whereas studies reporting mixed or overconnectivity findings mostly examined effects across the whole brain.

Five ROIs each were identified for VS and RSVP based on the strongest activation clusters each from three separate one‐sample t‐tests (for TD group, ASD group, and both groups combined), generating three sets of ROIs: TD‐biased, ASD‐biased, and combined activation ROIs. These three sets of five ROIs each constituted the restricted FOV for the VS and RSVP datasets. For whole‐brain FOV analysis, one seed, corresponding to the cluster of strongest activation, was selected from each set (see Supporting Information Tables I and III for ROI coordinates). For the resting‐state dataset, the restricted FOV was the default mode or task negative network as defined by Fox and colleagues [2005] (see Supporting Information Table III for ROI coordinates), and the primary seed for whole brain FOV analysis was in posterior cingulate cortex. All ROIs had approximately equal volume of c.1,000 µl (30–40 voxels). For each ROI, a 6 mm radius sphere was drawn, centered at the reported coordinates, and voxels falling outside the group grey matter mask were removed. If this removal resulted in an ROI with <30 voxels, the radius was increased until an ROI volume in the target range was identified.

An additional set of ROIs was also used for direct comparison across all three datasets. These ROIs were adopted from the landmark study by Just et al. [2004], which was the first to report functional underconnectivity in ASD (see Fig. 1 for ROIs). Freesurfer parcels [Desikan et al., 2006; Fischl et al., 2004] corresponding to the regions, for which underconnectivity was reported by Just and colleagues, were trimmed to the size of 30–40 voxels each, to maintain consistency of ROI volumes across all analyses. A seed in LIFG was used as primary seed for whole brain FOV analysis, based on the main activation effect reported in this study for the TD group and the finding of underconnectivity for ROI pairs involving IFG.

FCMRI Data Analysis: Co‐activation Versus Intrinsic Pipelines

Preprocessed data were subject to two analysis pipelines. The first pipeline was a co‐activation analysis (as described above), with no low‐pass filtering (NLP) or task regression. The second pipeline was an intrinsic fcMRI analysis, including low band‐pass (LBP) filtering at 0.008 < f < 0.08 Hz (applied to seed time series and all nuisance regressors) and task regression for the two task‐related datasets (VS, RSVP) to isolate intrinsic network‐specific BOLD fluctuations [Cordes et al., 2001; Lowe et al., 2000]. Orthogonal task regressors for each trial type were used (14 for VS, 9 for RSVP). Task regressors were convolved with a variable hemodynamic response, using AFNI's TENT function. Two different FOVs were applied in each pipeline: network‐specific (limited to ROIs) and whole brain. The ROIs adopted from Just et al. [2004] were only implemented for the intrinsic fcMRI pipeline, since no task‐activated data from the original study were available.

For each dataset, the average BOLD time course was extracted from each ROI in each set (as described above). In view of the known impact of head motion on BOLD correlations [Power et al., 2012; Van Dijk et al., 2012], several steps beyond conventional motion correction were taken in the preprocessing and analysis of each dataset to reduce the effects of head motion. Six rigid‐body motion parameters estimated from realignment of functional volumes were modeled as nuisance variables and removed with regression. Furthermore, time points with excessive motion (>1.5 mm displacement from one time point to the next) and their two immediate neighbors were censored from statistical analysis. A univariate ANOVA testing for group differences in motion (using the root mean square of displacement) across datasets indicated no main effects for dataset (p = 0.69) or group (p = 0.40), and no interaction effect between group and dataset (p = 0.89). For each dataset, time series from white matter and ventricular masks were also extracted for each participant using Freesurfer [Desikan et al., 2006; Fischl et al., 2004], trimmed to avoid overlap with grey matter. These time series were entered as nuisance variables in the general linear model testing for BOLD time series correlations between each seed ROI and other voxels in the brain.

Field of View

Analyses were performed either for all voxels across the entire brain or restricted to regions of interest. For ROI analyses (restricted FOV), data were analyzed in two ways: (1) For correlation matrices (presented in Fig. 1A,B), averaged time courses were extracted from each ROI and correlated against each other in a pairwise fashion. This served as a first‐pass approach for summary visualization of differences between co‐activation and intrinsic fcMRI pipelines (Fig. 1A,B). (2) For voxelwise ROI analyses (presented in Table 2), the averaged time course from each seed was tested for correlation with each voxel within the ROIs of the given network. This voxelwise approach served to compare the effects of restricted versus whole brain FOV, as part of an overall test of methodological variables presented in Table 2.

Table 2.

Overview of functional connectivity findings for different datasets and methodological variables [Color table can be viewed in the online issue, which is available at http://wileyonlinelibrary.com.]

Open in a new tab

Font color (web version only) indicates two main processing pipelines (red: co‐activation [retaining task effects without low‐pass filtering]; blue: intrinsic FC [isolating intrinsic fluctuation through low bandpass filtering and task regression]). NLP, no low‐pass filter; LBP, low bandpass filter. Seeds for each analysis are stated in the main text. Note that the LIFG seed corresponds to main ROI from Just et al. [2004], based on their activation and FC findings. Fill color in right column indicates main patterns of between group findings (UR ≥ 0.25: underconnected; UR ≤ −0.25: overconnected; all else: mixed). [Color Table can be viewed in the online issue, which is available at http://wileyonlinelibrary.com.]

Main effect p < 0.05;

(U − O)/(U + O)U, O = number of underconnected, overconnected voxels (p < 0.05, uncorrected).

Underconnectivity Ratio

For this test of methodological variables, a univariate ANOVA was run using type of dataset, seed selection, type of pipeline (co‐activation vs. intrinsic), and field of view as independent variables and an underconnectivity ratio (UR) for each condition as the dependent variable (Table 2). The UR was calculated as: UR = $\frac{U-O}{U+O}$, where U and O equal the number of underconnected and overconnected voxels (p < 0.05, uncorrected), respectively, derived from direct between‐group comparisons. For detailed justification of the selected significance threshold, see Supporting Information (p. 2).

Underconnectivity Versus Overconnectivity Results

Using cutoffs as stated in Table 2 for qualitative labeling of the overall pattern of findings for each analysis (listed in the right most column), we found that 9 of 24 analyses resulted in underconnectivity findings. All nine of these were generated by the co‐activation pipeline (no task regressor or low bandpass filter). Among these nine, six were generated by analyses with FOV restricted to activation ROIs. Two of the nine resulted from the VS, six from the RSVP, and one from the RS dataset. Conversely, 10 of the 24 analyses resulted in overconnectivity findings. Of these 10, six were generated by the intrinsic FC pipeline (with task regressors and low bandpass filter), and four by the co‐activation pipeline. Of the latter, three were generated by analyses with whole‐brain FOV. Eight of the 10 resulted from the VS, one from the RSVP, and one from the RS dataset.

An additional 5 of 24 analyses resulted in mixed connectivity findings. Four of these five were generated by the task‐regressed pipeline with low‐pass filter. All five were generated by analyses with whole‐brain FOV. Three of the five resulted from the RSVP dataset, and two from the RS dataset.

Individual Methodological Variables

Patterns of results differed for type of dataset, with two under‐ and eight overconnectivity findings for VS dataset, contrasted by six under‐ and only one overconnectivity finding for RSVP dataset. Results for RS dataset were overall mixed.

For the ROI selection variable, three of six analyses using seeds based on TD activation yielded underconnectivity, two overconnectivity, and one mixed connectivity findings. For the six analyses using seeds based on ASD activation, four yielded overconnectivity findings, and two underconnectivity findings. For the six analyses using seeds based on combined activation of both groups, three yielded underconnectivity, two overconnectivity, and one mixed connectivity findings.

For the variable co‐activation versus intrinsic FC analysis (i.e., task effect and temporal filtering), 6 of 10 analyses for task‐related datasets (VS and RSVP) yielded overconnectivity findings for intrinsic FC analysis, while the other four produced mixed connectivity effects. Of the 14 co‐activation analyses using no low band‐pass filter and task regressors, nine resulted in underconnectivity, four in overconnectivity, and one in mixed connectivity findings. Correlation matrices for all ROI pairings contrasting results from co‐activation versus intrinsic FC pipelines, and examples of contrasting findings are shown in Figure 1.

For the FOV variable, six of seven analyses with FOV restricted to activation ROIs yielded underconnectivity findings, with one yielding overconnectivity findings. Of the 17 analyses using whole‐brain FOV, three resulted in underconnectivity findings, nine in overconnectivity finding, and four in mixed connectivity findings.

For the seed in LIFG selected based on Just et al. [2004], only the intrinsic FC pipeline with whole brain FOV was applied for all three datasets (as task‐activated data for sentence comprehension were not available). Two of these analyses yielded overconnectivity, and one mixed findings.

Underconnectivity Ratio

Results from the univariate ANOVA for the UR indicated three main effects (Fig. 2). We detected a main effect for type of pipeline (temporal filtering/task regression; F[1,22] = 8.77, p = 0.007), with the co‐activation pipeline (task‐activated, no low‐pass filter) obtaining a higher mean UR (M = 0.38) than the intrinsic FC pipeline (task‐regressed, low band‐pass filter; M = −0.35). A main effect was also found for FOV (F[1,22] = 7.87, p = 0.01), with a higher mean UR for ROI‐restricted FOV (M = 0.62) than for whole‐brain FOV (M = −0.15). An additional main effect was found for type of dataset (F[1,21] = 4.38, p = 0.04), with higher mean UR for the RSVP (M = 0.35) than the RS dataset (M = −0.10), and lowest mean UR for the VS dataset (M = −0.41). No significant interaction effects were detected.

Additional Analyses

Following exclusion for head motion, samples were no longer sex‐matched for RSVP and RS datasets and age‐matched for the RS dataset. Since the main comparisons in this study were between data analysis pipelines for identical datasets, this was unlikely to affect the pattern of results. Nonetheless, we performed secondary analyses for fully sex‐ and age‐matched subsamples. The overall pattern of findings was similar for these subsamples (see Supporting Information Table IV).

Intrinsic Versus Co‐activation FC

As hypothesized, differential data processing along co‐activation versus intrinsic FC analysis pipelines had significant effects on the overall pattern of between‐group findings, with co‐activation yielding mostly underconnectivity and intrinsic FC analysis mostly overconnectivity effects in the ASD group. This highlights the importance of a strict distinction between co‐activation and intrinsic FC analyses. Insufficient awareness of this distinction can be considered the major reason for the often stark inconsistencies in the fcMRI literature on ASD and the resulting misperceptions [Hughes, 2007]. As described in the Introduction, intrinsic fcMRI findings are more likely to reflect task‐independent network organization. Spontaneous neuronal activity fluctuations are observed at many temporal scales, from the very low frequencies detected in BOLD fMRI to high gamma [>80 Hz, Canolty et al., 2007], which can only be detected through intracranial recording and possibly in magnetoencephalography [Khan et al., 2013]. Oscillations at different frequencies have been found to be linked in a nesting relationship of cross‐frequency coupling, with the phase of lower‐frequency oscillations driving amplitude at higher frequencies [Canolty and Knight, 2010]. Specifically relevant to fcMRI is the observation that low‐frequency oscillations regulate distributed neuronal networks and are susceptible to environmental entrainment [Canolty and Knight, 2010]. Learning induced changes in network‐specific synchronized low‐frequency fluctuations have indeed been detected in fcMRI studies [Lewis et al., 2009; Schultz et al., 2012; Wang et al., 2012]. Intrinsic fcMRI, which isolates such distributed and synchronized very‐low frequency oscillations, is therefore considered to reflect large‐scale Hebbian synaptic changes, providing a window onto fundamental mechanisms of functional network organization [Johnson, 2011]. The finding of partial overconnectivity in ASD from intrinsic fcMRI analyses will be interpreted with regard to its possible neurofunctional relevance further below in the context of field of view.

The mentioned fundamental neurobiological relevance of intrinsic fcMRI does not imply that co‐activation analyses, which tend to yield underconnectivity findings—both based on the current study and the broader literature [Müller et al., 2011]—are not informative. However, one needs to be aware that findings from the co‐activation technique will predominantly reflect the effects of online cognitive and sensorimotor processing. If this processing occurs in a domain of known or suspected cognitive‐behavioral impairment, the finding of reduced BOLD correlations is expected. As discussed further below, the inverse is true for domains of relative strength in ASD [Keehn et al., 2013b].

Field of View

The detection of a significant main effect of FOV has both pragmatic and conceptual implications, although this finding requires some caution, given the small number of voxels entered into the computation of the underconnectivity ratio for several pipelines. Focus on ROIs selected based on their known or expected relevance for a functional domain of interest (e.g., face perception or theory of mind) has been rather common in ASD fcMRI studies because it may be considered to protect from detection of unexpected effects elsewhere in the brain that are hard to interpret. Selective presentation of findings was obvious, for example, in two studies [Agam et al., 2010; Ebisch et al., 2011], which also showed partial overconnectivity findings in Figures (5 and 2–3, respectively; confirmed by personal communication with authors), but exclusively discussed underconnectivity effects occurring within ROIs. Neglect of unexpected overconnectivity findings is equivalent to the exclusion of empirical “anomalies” in the context of normal science [Kuhn, 1962]—here: the theory of general underconnectivity in ASD—and therefore not surprising. However, such neglect runs the risk of missing crucial features of aberrant network organization in ASD.

A few recent studies have shown some indication of the underlying network characteristics. These studies suggest that network differentiation and specialization, which are emerging features of development in the typical brain [Johnson, 2011], may be impaired in ASD. Two of these [Rudie et al., 2012; Shih et al., 2011] reported reduced network differentiation in ASD, using seed‐based fcMRI analyses. Concordant findings for both functional and anatomical network organization from graph theory were more recently reported by Rudie et al. [2013]. While conceptually related, findings from the two groups differ in important respects. Whereas the studies by Rudie et al. [2012, 2013] suggest that overconnectivity findings may actually reflect reduced negative connectivity observed in TD groups, which could be interpreted as reduced inhibitory connectivity, this explanation does not account for findings by Shih et al. [2011] or similar results from other studies [Keehn et al., 2013b; Shen et al., 2012], all of which detected overconnectivity in ASD groups between brain regions that showed positive BOLD correlations in TD control groups. The latter was also true for analyses yielding overconnectivity effects in the present study.

Divergent findings may relate to the ongoing debate about global signal regression [GSR; Fox et al., 2009; Murphy et al., 2009]. Jones et al. [2010] first suggested that overconnectivity in ASD might be exclusively observed after GSR. In partial support, supplementary analyses in Rudie et al. [2012] show expansion of overconnectivity effects in ASD after GSR. However, while the shift in the distribution of correlations through GSR may surely have an effect on between‐group comparisons, it is unlikely that overconnectivity findings can be explained in this way. Several studies reporting (partial) overconnectivity in ASD did not employ GSR [Keehn et al., 2013b; Monk et al., 2009; Noonan et al., 2009; Shen et al., 2012; Shih et al., 2010, 2011; Welchew et al., 2005]. Since our precursor study [Müller et al., 2011] did not suggest any consistent effects of GSR, we did not implement GSR as a variable in the current study. Overconnectivity was found along numerous pipelines in three different datasets (Table 2), always in the absence of GSR. In addition to these empirical findings, methodological considerations also appear to speak against the use of GSR. With the implementation of nuisance regressors extracted from white matter and CSF, it is likely that some global signal fluctuations caused by signal instability or other global noise components can be statistically removed without GSR. Since the global signal probably includes true BOLD fluctuations across widespread brain networks [Schölvinck et al., 2010], additional use of GSR may therefore do more harm than good [Gotts et al., 2013; Saad et al., 2012].

Type of Dataset

We also found a significant main effect for type of dataset. Whereas the VS dataset yielded predominantly overconnectivity effects even when task effects were regressed out as nuisance variables, the RSVP dataset yielded mostly underconnectivity effects. The most likely explanation relates to the domain of processing the two different tasks tapped into. Visual search has been referred to as an “island of sparing,” with normal or even superior function in ASD [Keehn et al., 2013b; O'Riordan et al., 2001]. Conversely, behavioral performance on the RSVP task was reduced in ASD (compared to TD participants), suggesting that this task tapped into a subdomain of attention that is impaired in ASD [Allen and Courchesne, 2001; Keehn et al., 2013b] Our findings for VS are consistent with an often overlooked prediction of the underconnectivity theory by Just et al. [2004] that functional connectivity may actually be increased for domains of strength in ASD.

While performance differences between VS and RSVP would obviously affect co‐activation results, as indeed observed, note that even intrinsic FC analyses were impacted by task domain since seeds were selected based on activation for the respective tasks. However, findings for a seed in LIFG suggest that dataset‐specific differences were further affected by factors other than seed selection. Although predominance of overconnectivity (reflected in a negative UR) was found for all three datasets processed along the intrinsic whole brain FOV pipeline, its degree varied from robust (UR ≤ −0.75) for VS and RS to subtle (UR = −0.13) for RSVP. In particular, differential results for VS and RSVP may suggest that even thorough measures to remove effects of task‐related online processing (through multiple task regressors for each stimulus type and low bandpass filtering) were not entirely successful in isolating intrinsic non‐task driven BOLD fluctuations [cf. Fair et al., 2007]. However, note that cohorts for the two studies differed (e.g., VS included ASD participants under 12 years of age), which may have also contributed to divergent results.

The possibility of residual task effects confounding intrinsic fcMRI results may favor the use of RS, but caution is warranted. It is likely that individuals with ASD do not enter a typical self‐referential “default mode” as readily during instructed “resting” as TD controls [Kennedy et al., 2006]. In fact, little is known about what cognitive states and types of mind wandering [Mason et al., 2007] people with ASD go through in the absence of a clearly defined task. Since the RS does not permit monitoring of cognitive processing through measures of task performance, potential confounds due to systematic differences in response between TD and ASD cohorts are hard to rule out definitively. As a general conclusion, the recommendation must be to implement both RS and task‐regressed paradigms to identify consistency of intrinsic FC effects across datasets. The implementation of simple, but continuously engaging tasks, the effects of which may be more cleanly regressed out, could be advantageous.

Seed Selection

No significant effect was found for seed selection. However, this null finding may be attributed to the inclusion of only two task datasets in the study, which limited the number of pipeline comparisons. As discussed above, the two task datasets also differed heavily in the distribution of underconnectivity ratios (being shifted towards overconnectivity for VS, but underconnectivity for RSVP), which may have further reduced our ability to detect effects of seed selection. Note that inspection of individual analyses mostly showed differences in mean URs in the expected direction (higher for seeds based on TD activation than for those based on ASD activation), when co‐activation analyses were considered. Our null finding therefore by no means indicates that careful avoidance of group bias in the selection of seeds would be irrelevant. On the contrary, such selection is of great importance specifically in co‐activation analyses. However, compared to other methodological issues discussed above, seed selection bias has been avoided in many fcMRI studies of ASD and therefore does not constitute a major topic of debate.

Limitations

While we tested for between‐group effects along a large number of analysis pipelines, numerous other methodological, demographic, and clinical variables might affect fcMRI findings in ASD. We focused on contrasting findings from co‐activation and intrinsic fcMRI analyses. For pragmatic reasons (avoiding an unmanageable plethora of findings from even more pipelines), we combined task regression and low bandpass filtering, both of which serve the isolation of spontaneous low‐frequency fluctuations, and were therefore unable to test for any differential effects of these variables. Our study resorted to datasets with sample sizes of c.20 per group, which were not available for other task domains. In view of the stark differences between VS and RSVP datasets, it is likely that inclusion of data from different task domains might have affected some of our findings. Note, however, that the direction of differences in findings between pipelines was overall consistent across datasets.

Our cohorts were limited to children and adolescents, leaving the possibility that fcMRI methods may affect adult ASD in different ways. However, this appears unlikely given that our literature meta‐analysis [Müller et al., 2011], which included many adult studies, yielded highly consistent results. Our ASD cohort was also limited to high‐functioning participants, given the need for compliance with motion restraints. Very little is known about functional brain network organization at the low‐functioning end of the autism spectrum.