Altered Global Signal Topography in Schizophrenia

Participants

We studied two independent clinical samples that were combined for the purposes of the present analysis (Table S1). One sample comprised 90 chronic SCZ patients and 90 demographically matched HCS recruited at the Olin Neuropsychiatry Research Center through outpatient clinics and community mental health facilities in the Hartford. The second sample, comprised of 71 SCZ and 74 demographically matched HCS, was obtained from a publicly distributed dataset provided by the Center for Biomedical Research Excellence (COBRE) (http://fcon_1000.projects.nitrc.org/indi/retro/cobre.html). All participants provided informed consent approved by the institutional review board of relevant institutions. Across samples, all subjects met identical neuroimaging exclusion criteria, underwent identical preprocessing and analyses. SCZ symptom severity, across both samples, was determined using the Positive and Negative Syndrome Scale (PANSS) (Kay et al. 1987) (Table S1).

Inclusion & Exclusion Criteria

Patients were included based on the following criteria: (1) SCZ diagnosis as determined by the Structured Clinical Interview (SCID) for the Diagnostic and Statistical Manual of Mental Disorders-IV (First et al. 2002), administered by experienced MA or PhD-level research clinicians; (2) no major medical or neurological conditions (e.g. epilepsy, migraine, head trauma with loss of consciousness); and (3) premorbid IQ > 70 as assessed by National Adult Reading Test, Wide Range Achievement Test or Wechsler Test of Adult Reading depending on the study protocol (Spreen and Strauss 1998). As in prior studies, these measures were normed and converted to IQ equivalents for each subject. If more than one premorbid achievement measure was available per subject the scaled scores were averaged per standard practice (Lezak 1995). To increase ecological validity of the patient sample (Krystal et al. 2006), comorbid Axis I anxiety disorders and/or history of substance abuse (fully remitted >6 months before the study) were allowed. HCS participants were excluded if they had: (1) history of medical or neurological disorders, (2) history of mental retardation, (3) history of severe head trauma with >5 min loss of consciousness; (4) history of substance abuse or dependence within the last 12 months, (5) current or lifetime Axis I psychiatric disorder as assessed by SCID-NP; (6) history of psychotic disorders in first-degree relatives (reported by detailed family history). Across samples we accomplished matching on a number of relevant demographic variables. Patients and HCS did not significantly differ on any of the variables (see Table S1), apart from educational attainment and premorbid intellectual functioning (IQ), which were lower for SCZ. These differences are impacted by the illness course (Glahn et al. 2006) (43) and thus were not included as a covariate. All medication was converted to chlorpromazine equivalents (Andreasen et al. 2010) and verified by trained raters. Notably, smoking status and medication dose did not alter reported effects (see Figure S4 for comprehensive confound analyses).

Neuroimaging Data Acquisition

At the Olin Neuropsychiatry Research Center 90 SCZ and 90 HCS were scanned  using a Siemens-Allegra 3 T scanner with identical acquisition parameters (Anticevic, Brumbaugh et al. 2012). Images sensitive to BOLD signal were acquired with axial slices parallel to the anterior-posterior commissure (AC-PC) using a T2*-weighted gradient-echo, echo-planar sequence (TR/TE = 1500/27 ms, flip angle = 60°, field of view = 24 × 24 cm, acquisition matrix = 64 × 64, voxel size = 3.43 × 3.43 × 4 mm). The acquisition lasted 5.25 min and produced 210 volumetric images per subject (29 slices/volume, inter-slice gap = 1 mm). Subjects were instructed to lay awake in the scanner and keep their eyes open. Subjects were monitored on a video camera to ensure that they stayed awake and were removed from the analyses if they fell asleep during the scan or if their head movement exceeded 1 mm along any axis. Structural images were acquired using a T1-weighted, 3D magnetization-prepared rapid gradient-echo (MPRAGE) sequence (TR/TE/TI = 2200/4.13/766 ms, flip angle = 13°, voxel size [isotropic] = 0.8 mm, image size = 240 × 320 × 208 voxels), with axial slices parallel to the AC-PC line.

As noted, 71 additional SCZ patients and 74 matched HCS underwent data collection at COBRE using a Siemens Tim-Trio 3 T scanner. Full acquisition details are reported previously (Stephan, Penny et al. 2009; Hanlon et al. 2011; Mayer et al. 2012). Briefly, BOLD signal was collected with 32 axial slices parallel to the AC-PC using a T2*-weighted gradient-echo, echo-planar sequence (TR/TE = 2000/29 ms, flip angle = 75°, acquisition matrix = 64 × 64, voxel size = 3 × 3 × 4 mm). The acquisition lasted 5 min and produced 150 volumetric images per subject. Structural images were acquired using a 6-minute T1-weighted, 3D MPRAGE sequence (TR/TE/TI = 2530/[1.64, 3.5, 5.36, 7.22, 9.08]/900, flip angle = 7°, voxel size [isotropic] = 1 mm, image size = 256 × 256 × 176 voxels), with axial slices parallel to the AC-PC line. Described parameters were obtained via the website (http://fcon_1000.projects.nitrc.org/indi/retro/cobre.html).

Neuroimaging Preprocessing

Preprocessing followed prior work (Cole et al. 2011; Anticevic, Brumbaugh et al. 2012); complete fMRI acquisition and preprocessing details are presented in Supplement together with some limitations to these legacy acquisition and preprocessing methods. Notably, as part of standard preprocessing, we used multiple regression to model the time-varying BOLD signal in each voxel, including a nuisance regressor for the mean signal across gray matter aka the GS (as well as other typical nuisance regressors, e.g. for white matter, ventricular signal, and rigid body motion). This step was performed after the data were scrubbed for high motion (see Supplementary Material). From this multiple regression, we obtained beta weights for each regressor variable. Presented analyses concern the voxel-wise beta weights obtained for the nce regressor for the mean signal across GS, assessed in SCZ and HCS.

Global Gray Matter Signal Beta Map Calculation

To obtain GS beta values, we first performed GSR using standard widely adopted procedures (Cole et al. 2011; Anticevic, Brumbaugh et al. 2012). The GS timeseries for each subject was obtained by calculating mean raw BOLD signal averaged over all gray matter voxels for each time point, explicitly excluding ventricles and white matter signal (which are defined as separate nuisance regressors). This GS timeseries was used as nuisance predictor term within a multiple linear regression model along with other nuisance predictor terms (ventricular signal, white matter signal, movement parameters, and the first derivatives of each of these with respect to time). Critically, we did not use the first derivative of the GS as an additional nuisance regressor, as this would have produced 2 GS-related beta weights rather than one. We computed the multiple regression with a single GS-related variable to yield a more parsimonious interpretation of the obtained voxel-wise GS beta weight. More formally, we used the following multiple regression analysis:

where $\mathrm{BOL}{\mathrm{D}}_k^{\mathrm{raw}}\left(t\right)$ represents the raw BOLD signal in voxel k as a function of time, t. b₀ is the intercept, containing trend parameters for the scan, X_i represents the ith nuisance regressor (e.g. ventricles or GS), b_i is the corresponding beta weight computed for regressor X_i. The last term is the residual signal that is not accounted for by the nuisance regressors. In other words, the residual represents the preprocessed BOLD signal at voxel k. As noted, the n nuisance regressors include GS(t), ventricular signal V(t), white matter signal WM(t), 6 movement parameters, and first derivatives of the latter 8 with respect to time. The GS beta weights reported are represented by the b_GS values obtained from this multiple regression. GS(t) is the spatial average of time-varying BOLD signal across all gray matter voxels:

The “Mean GS beta weight” computation in Figures 1F, 2D–E, 3C–D, and 4 is done by simply fitting a generalized linear model (GLM) to each voxel's BOLD time series to obtain the GS beta weight ($b_{\text{GS}}$) as shown above, and then averaging across each voxel belonging to the region of interest. Importantly, the voxel-wise whole-brain map of GS beta weights is more akin to a task-based GLM analysis than to a functional connectivity measure. In other words, GS beta weights are not functional connectivity values and should not be interpreted as such.

Functional Network Definition

In addition to the voxel-wise regression analyses, we used an a priori functional network parcellation to explicitly quantify the presence (or absence) of GS beta spatial shifts in SCZ. For this purpose, we employed a preexisting functional parcellation computed by Power et al. (2011). We defined boundaries for the following higher-order association networks: default mode network (DMN), frontoparietal control network (FPCN), and ventral attention network (VAN). Similarly, we isolated the visual, auditory, and somatosensory networks. These network boundaries were used to provide independent validation of voxel-wise effects (Power et al. 2011).

Defining Spatial Localization of the GS in Healthy Individuals

We first examined the spatial pattern of GS representation in healthy subjects. Of note, instead of computing voxel-wise correlations to the GS—as was done in earlier work (Aguirre et al. 1997; Zarahn et al. 1997; Fox et al. 2009), here we utilize beta coefficients obtained using multiple regression (see Exp. Procedures available in Supplementary Material). This measure of GS representation has some advantages when aiming to account for GS along with multiple nuisance regressors simultaneously. Critically, the GS component in a given voxel's BOLD signal is expressed as $b_{\mathrm{GS}}\mathrm{GS}\left(t\right),$ where the regressor variable GS(t) is a function of time, and the beta coefficient b_GS is spatially-varying across voxels (see Exp. Procedures available in Supplementary Material). We hypothesized that this proportion of GS (b_GS) at each voxel may show a distinct spatial pattern across the brain in healthy subjects. We quantified the beta coefficients, b_GS, across all gray matter voxels in a sample of 164 healthy subjects, producing a group-level T-map of voxel-wise GS beta values (Fig. 1A,B, see Table S2 for region coordinates and statistics surviving whole-brain type I error correction). We observed a strikingly nonuniform representation of GS across functional brain networks: b_GS values were significantly greater in sensory than in association network regions (Fig. 1F) [t(163) = 9, paired samples test, P < 4 × 10⁻¹⁶]. To clarify, this analysis compared 2 measures (mean association and mean sensory region GS beta weight) within each subject and not between subjects: The voxel-wise GS beta weights were averaged across association and sensory regions within a given single subject first—iterating this process produced a pair of subject-specific values for each subject, using only data from that subject. These subject-specific values for mean association and mean sensory GS beta weights were then compared for all healthy subjects in a paired samples t-test. Further, mean b_GS in sensory versus association networks exhibited a significant anti-correlation across healthy subjects—those with larger sensory b_GS had smaller association b_GS (Fig. 5C) [r = −0.37, P < 10⁻⁶].

GS Representation is Spatially Altered in SCZ

Prior work showed that GS may be altered in SCZ relative to HCS (Yang et al. 2014). We hypothesized that the spatial pattern of GS representation may be altered in SCZ. To test this, we mapped b_GS values in SCZ, exactly as was done for the healthy subjects (Fig. 1 & S1). Then we computed a T-map comparing SCZ versus HCS groups to identify areas in SCZ that might be altered in their statistical relationship to the GS. Results revealed significant group differences in whole-brain voxel-wise b_GS values (Fig. 2), with distributed clusters of significant elevations and reductions in voxel-wise b_GS for SCZ (Fig. 2A, see Table S3 for region coordinates and statistics surviving whole-brain type I error correction). In addition, the unthresholded T-map (Fig. 2B) and Cohen's D map (Fig. 2C) show a pattern of b_GS increases and reductions that qualitatively co-localize to association versus sensory networks respectively (see quantitative analysis in Fig. 4). The whole-brain T-maps of GS beta weights for SCZ and HCS groups are separately shown in Figure S1.

Notably, the overall sample was collected across 2 sites (Table S1). We observed a similar qualitative pattern of GS beta changes when examining unthresholded between-group maps separately for each sample (Fig. 3A,B). For the full sample, effects extracted from the type I error corrected ROIs quantitatively confirmed increases [t(322) = 13, P < 3 × 10⁻¹⁶, Cohen's d = 1.5] and reductions [t(267) = 11, P < 3 × 10⁻¹⁶, Cohen's d = 1.2] in SCZ relative to HCS (Fig. 2D,E). Finally, this bi-directional effect was quantified via a 2-way ANOVA analysis examining a “Group × ROI” interaction [F(1,323) = 192.2, P < 2 × 10⁻¹⁶]. In addition, these “Group” effects replicated when examining samples separately by site (Fig. 3C,D) for both the mean elevated GS beta ROI [t(160) = 9, P < 3 × 10⁻¹⁶, Cohen's d = 1.4 for Sample 1; t(128) = 9, P < 3 × 10⁻¹⁶, Cohen's d = 1.6 for Sample 2], and the mean reduced GS beta ROI [t(168) = 10, P < 3 × 10⁻¹⁶, Cohen's d = 1.6 for Sample 1; t(108) = 7, p < 2 × 10⁻¹⁰, Cohen's d = 1.2 for Sample 2]. Again, these effects were quantified via a 2-way ANOVA analysis examining a Group x ROI interaction [F(1,178) = 121.8, P < 2 × 10⁻¹⁶ for Sample 1; F(1,143) = 87.3, P < 2 × 10⁻¹⁶ for Sample 2].

In turn, we quantified the spatial similarity between the obtained T-maps for samples 1 & 2, using a measure of voxel-wise similarity (eta²) (Cohen et al. 2008). The aim was to quantify the spatial reproducibility of the between-group differences in GS topography shifts (see SI Appendix for details). Figure 3 shows the quantified similarity index (eta² = 0.994 for full range of T-values), revealing a high degree of spatial reproducibility for the between-group GS topography differences across samples.

Next, we examined if altered GS localization increases and reductions in SCZ were quantitatively related. Using the T-map results shown in Figure 2A, we defined specific regions of significantly elevated b_GS in SCZ (the mean elevated GS beta ROI). Similarly, we defined regions of significant b_GS reductions in SCZ (the mean reduced GS beta ROI). We averaged voxel-wise b_GS values in each large ROI. Results indicated that b_GS values for the 2 sets of ROIs were negatively correlated for all groups (SCZ-only: r = −0.68, P < 10⁻¹⁵; HCS only: r = −0.37, P < 10⁻⁵; and all subjects: r = −0.69, P < 10⁻¹⁵, Fig. 2D).

Of note, we identified a significant relationship between head motion and b_GS values in the 2 sets of ROIs. We observed that higher head motion in patients correlated with higher [r = 0.26, P < 0.001] and lower b_GS values [r = −0.4, P < 4 × 10⁻⁷] in ROIs that were associated with mean elevated and reduced GS localization, respectively. To account for the effect of head motion on our measures, we computed an ANCOVA with head motion as a continuous covariate (measured in % of frames flagged for movement, see Exp. Procedures available in Supplementary Material). Group main effects remained highly significant even after accounting for motion for both ROIs [F(1,321) = 155.1, P < 2 × 10⁻¹⁶ for the elevated, and F(1,321) = 101.9, P < 2 × 10⁻¹⁶ for the reduced GS localization ROI]. Furthermore, after repeating analyses for a subset of explicitly movement-matched subjects (130 SCZ & 130 HCS) all effects remained highly significant (Figures S2 & S3). Finally, we verified the association versus sensory network-preferential effects via an a priori network analysis presented below.

Characterizing Spatial Patterns of GS Beta Changes in SCZ

Voxel-wise results revealed a spatial dissociation between regions of b_GS increases and decreases in SCZ, which qualitatively matched association and sensory areas, respectively. To formally quantify the spatial mapping of these shifts, we computed the probability of overlap for these effects with a priori defined large-scale functional networks. Here, we defined association and sensory networks using a functional cortical parcellation (Power et al. 2011; Cole et al. 2013). For the “association network” ROI, we combined FPCN, DMN, and the VAN (Fig. 4A). Similarly, we defined the sensory network ROI by combining the somatosensory, auditory, and visual networks (Fig. 4A). “Non-association” and “non-sensory” network ROIs were defined as all cortical gray matter voxels not part of the relevant large-scale network (Fig. 4B,C). Consistent with qualitative observations, voxels showing significantly elevated b_GS in SCZ preferentially co-localized to association networks [binomial test for proportions, P < 0.001, Fig. 4B,D]. In turn, voxels showing significantly reduced b_GS preferentially co-localized to sensory networks [binomial test for proportions, P < 0.001, Fig. 4C,D]. A complementary analysis revealed significantly different proportions with respect to spatial coverage of association versus non-association networks by the “mean elevated GS beta” ROI voxels [P < 0.001, Fig. 4D]. Similarly, the “mean reduced GS beta” ROI showed significantly different spatial coverage of sensory versus non-sensory networks [P < 0.001, Fig. 4D]. Collectively, the overlap analyses indicate that the mean elevated GS beta ROI was predominantly co-localized with association cortices, whereas the mean reduced GS beta ROI was predominantly co-localized with sensory cortices.

Independently Examining GS Representation Shifts in Association versus Sensory Regions

Presented results, in line with recent empirical reports (Yang et al. 2016), suggest network-specific b_GS shifts in SCZ. Overlap analyses indicate preferential elevations in b_GS in association networks, contrasting with preferential b_GS reductions in sensory networks in SCZ. Next, we sought to confirm these effects via an independent a priori analysis examining association versus sensory network ROIs. We explicitly examined b_GS value differences for SCZ versus HCS across a priori functionally defined large-scale networks. Specifically, we averaged voxel-wise b_GS across 3 large-scale association networks: FPCN, DMN, and VAN (Fig. 5A), extending our initial focused network characterization (Fig. 4). We found significantly elevated voxel-wise b_GS across association regions in SCZ compared with HCS [t(322) = 3.4, P < 0.001, Cohen's d = 0.37] (Fig. 5B,C). Next, we averaged voxel-wise b_GS across 3 large-scale sensory networks (Fig. 5A), which revealed significantly reduced voxel-wise b_GS across sensory regions in SCZ compared with HCS [t(322) = 3.8, P < 0.0002, Cohen's d = 0.42] (Fig. 5B,C). This effect was confirmed via a Group x Network interaction (2-way ANOVA for SCZ vs. HCS, association vs. sensory networks [F(1,323) = 19.7, P < 2 × 10⁻⁵]). In addition, whereas b_GS weights are significantly higher in sensory compared with association regions in HCS [t(163) = 9, P < 4 × 10⁻¹⁶], the effect size for this network-level difference is reduced in SCZ [t(160) = 3, P = 8 × 10⁻⁴] (Fig. 5B). We further verified that these independent effects were not driven by amount of head motion or medication dose (Fig. S4A–D). For subjects with available smoking status information, we verified that there were no significant differences between smoker and nonsmoker SCZ patients with respect to reported effects (Fig. S4E,F). Collectively, these empirical effects independently suggest that association and sensory regions show dissociable GS localization changes in SCZ.

Finally, we again verified that the association and sensory b_GS effects were quantitatively related (Fig. 5C). As noted above, this analysis revealed a significant negative correlation between mean association and sensory region b_GS values across all HCS [r = −0.37, P < 10⁻⁶]. This negative relationship remained significant in SCZ [r = −0.16, P = 0.05], but the effect size was significantly attenuated [test for difference between correlation coefficients: z = −2.03, P = 0.042, 2-tailed]. Collectively, these results suggest that the pattern of b_GS changes in SCZ may relate to an overall loss of “differentiation” in association versus sensory networks with respect to GS representation.

Of note, some studies (Mueller et al. 2013) have reported greater inter-subject variability in association regions (compared with non-association regions) with respect to functional connectivity. We were thus motivated to examine whether a nonfunctional connectivity measure, such as voxel-wise GS beta weight, might also show greater variability in the association regions. We formally tested this via homogeneity of variance tests, using Bartlett's test and the Fligner–Killeen test. The null hypothesis in this case may be stated as: Let A = the mean voxel-wise GS beta value in the association cortex for a subject; B = the mean voxel-wise GS beta value in the sensory cortex for a subject. The null hypothesis is that the distributions of A and B across subjects have the same variance. This hypothesis was rejected in the healthy subjects group [Bartlett's K-squared = 8.8, P = 0.003, df = 1; Fligner–Killeen chi-squared = 6.4, P = 0.01, df = 1], with mean association GS beta weight per subject being more variable than the mean sensory GS beta weight per subject, as one might expect if there is overall more variability between individuals in the association regions. Interestingly, the null hypothesis was not rejected for the SCZ group [Bartlett's K-squared = 2.3, P = 0.13, df = 1; Fligner–Killeen chi-squared = 2.4, P = 0.12, df = 1], again consistent with the idea of less “differentiation” between the association and sensory networks in SCZ.

Characterizing GS Topography in Healthy Individuals

While the principal objective of this study was to characterize GS topography shifts in SCZ, the issue of GS topography in healthy individuals has received little quantitative attention in the literature. Seminal work on this topic indeed showed a generally positive correlation between the GS and voxel-wise BOLD signals in the brain's gray matter (Aguirre et al. 1997; Zarahn et al. 1997; Fox et al. 2009). Yet, an important knowledge gap remained regarding potential network-by-network differences in the distribution of the GS in healthy adults. As noted, we predicted a distinct pattern of GS topography across sensory and association networks given their distinct computational roles in processing time-locked incoming external stimuli versus internal maintenance and representation of parallel information, respectively. As an initial demonstration that voxel-wise GS contributions are indeed nonuniform, we formally quantified GS topography using beta coefficients obtained from a multiple regression. This approach accounts for contributions of GS to the voxel-wise signal along with contributions of multiple nuisance regressors simultaneously, with the GS beta coefficient yielding a quantitative measure of voxel-wise GS contribution. As hypothesized, we found a strong difference in GS representations across a priori defined association versus sensory networks (Power et al. 2011), potentially highlighting an intriguing mathematical consequence of the differing computational tasks these networks are specialized to accomplish. However, it remains possible that these dissociable contributions reflect a spatially preferential yet pervasive noise component (e.g. breathing) (Golestani et al. 2015), rather than a true underlying distribution of neuronal signals. Future experimental work that explicitly dissociates physiological and neuronal artifact will be needed to address the underlying mechanisms of this nonuniform GS topography. Nevertheless, it is vital to characterize if this observed GS topography may be altered in severe mental illnesses such as SCZ—a central objective of this study. In turn, such alterations may inform strategies for development of pharmacotherapies, guided by insights from neuroimaging.

Localization of Widespread BOLD Signal Disruption in SCZ

Accumulating evidence from SCZ studies implicates deficits in brain regions responsible for higher-order cognition (Pearlson et al. 1996; Ross and Pearlson 1996; Cannon et al. 2002), with well-established disruptions in executive and working memory functions performed by the frontoparietal networks (Goldman-Rakic 1991; Weinberger et al. 1991; Carter et al. 1998; Barch and Ceaser 2012; Cole et al. 2014). In parallel, different (and perhaps less severe) abnormalities are observed in sensory regions (Woodward et al. 2012; Anticevic et al. 2014; Yang et al. 2014; Yang et al. 2016). This study provides a novel insight into findings implicating sensory versus association regions in SCZ. Put differently, the same pathological process can potentially produce different effects across brain regions that exhibit distinct neurophysiological properties. In healthy subjects, we observed that sensory region BOLD signals showed higher GS resemblance, whereas association areas showed lower GS resemblance. In SCZ, this distinction was shifted whereby association areas resembled the GS more strongly, but sensory areas resembled the GS less, relative to HCS. This bi-directional shift in GS representations could occur through a pathological process that “blurs” the distinction between association and sensory areas, effectively “homogenizing” cortical computations. This would, on average, yield GS elevations in association regions and reductions in sensory regions, as observed in the current study. In the next section, we consider potential neuropathological mechanisms in SCZ that could disrupt patients’ ability to (1) broadly integrate sensory information—leading to reduced GS in sensory-motor areas, and (2) maintain independent information representations in association areas—leading to impairments in higher-order cognitive control functions, which may be reflected in increased contribution to GS from association areas.

Understanding Mechanisms of Network-Level Cortical Disruptions in SCZ

As noted, identified cortical GS representation shifts in SCZ follow a functional dissociation between association and sensory cortices. These effects are qualitatively in line with findings observed in corresponding thalamo-cortical circuits (Marenco et al. 2012; Woodward et al. 2012; Anticevic et al. 2014). An open question for future investigation is the directionality of the disturbance in GS representation. Can the pattern of GS shifts be linked to a network-specific pathology, or does it involve a dynamical disturbance across inter-connected systems that makes directionality less central (Loh et al. 2007)?

Several studies suggest association network alterations identified via neuroimaging in SCZ may arise from a localized dysfunction (Lewis 2000; Hashimoto et al. 2003; Kolluri et al. 2005). Notably, the PFC is involved in gating information flow within the brain, and can exert inhibitory top-down control over thalamic nuclei through projections via the basal ganglia (Haber and McFarland 2001). In contrast, sensory regions are considered a source of drive onto the thalamus. Within this framework, we could hypothesize a “primary” association network disruption impacting GS representations as well corresponding thalamo-cortical circuits. In turn, this could lead to “secondary” disruptions, altering the balance of activity amongst thalamic nuclei, leading to concurrent (but not primary) global effects in sensory networks. Such a pattern of thalamo-cortical disruptions is consistent with prior reports (Marenco et al. 2012; Woodward et al. 2012; Anticevic et al. 2014).

Conversely, other SCZ neuropathology frameworks favor a more spatially distributed primary pathology across cortical circuits. Growing evidence from preclinical, pharmacological, and postmortem clinical studies implicates disrupted E/I balance across cortical microcircuits in SCZ (Marin 2012). This could reflect changes in a number of candidate neurotransmitter systems, involving glutamate (Krystal et al. 2003; Macdonald and Chafee 2006), GABA (Lewis et al. 2005), and dopamine pathways (Laruelle et al. 2003). Such changes potentially impact neural development, causing global impairments in integration/segregation of neural processing, as is suggested by studies of SCZ and their siblings (Repovs et al. 2011; Repovs and Barch 2012). How can we reconcile these hypothesized cellular-level disruptions with observed network-level BOLD abnormalities? Our prior computational modeling work identified altered E/I balance in cortical circuits as a key mechanism for disrupted long-range interactions between cortical areas (Anticevic, Gancsos et al. 2012), with preferential impact upon association networks at rest (Yang et al. 2016). This work offers a potential mechanism for the observed “blurring” between association and sensory networks with respect to GS localization. Here we show that the tightly anti-correlated relationship between these networks becomes attenuated in SCZ relative to HCS. This reduction in anti-correlated patterns between association and sensory GS representations may relate to the same disinhibition-mediated mechanism that reduces anti-correlations between task-positive and task-negative networks in SCZ (Anticevic, Gancsos et al. 2012). Furthermore, the observed opposing effects in association versus sensory networks may relate to preexisting intrinsic differences between these networks that alter the impact of E/I imbalance on these networks (Yang et al. 2016). In summary, these cortical models of global E/I imbalance suggest that disrupted long-range anti-correlations (Anticevic, Gancsos et al. 2012) and network-dependent pathological changes (Yang et al. 2016) might de-stabilize cortical information flow in ways observed presently: manifesting as altered GS representations for sensory versus association networks in SCZ. Future pharmacological and animal studies will be needed close the gaps between these emerging cellular hypotheses and system-level observations in SCZ.

Implications for Disrupted Cortical Information Flow in SCZ

Disruptions in E/I balance in SCZ may profoundly affect large-scale cortical function (Anticevic et al. 2013; Uhlhaas 2013). Indeed, prior computational studies show that E/I imbalance can alter functional connectivity and variance of BOLD signals across cortex (Yang et al. 2014, 2016). However, clinical studies of SCZ and HCS show that between-group differences across both measures are sensitive to GSR (Yang et al. 2014). This indicates that group differences may need to be considered in context of potential group-level changes in the GS representation itself, even if they partially reflect nuisance signal. Here, we quantified group-level GS shifts using 2 well-powered independent samples. Across samples, SCZ patients showed significantly reduced GS representation in sensory networks. Concurrently, we observed increased GS representation in association networks. Overall, network-level differences in GS representation were significantly reduced in SCZ. In contrast, HCS showed significantly greater GS localization in sensory compared with association networks. As noted, this complex spatial shift in cortical GS representation provides whole-brain evidence that, in SCZ, GS localizations may not differentiate well between the association and sensory cortices, possibly reflective of a “blurring” between sensory-association signals flowing through the brain.

Of note, reported effects show that these hyper- and hypo-localizations of GS are highly related: HCS with lower voxel-wise GS beta values in association networks exhibited increased GS beta values in sensory networks. In SCZ, this anti-correlated relationship was significantly reduced. This is an important consideration for treatments designed to restore neural function. For instance, interventions with regional selectivity, like low-frequency transcranial magnetic stimulation, might treat some symptoms (Hoffman et al. 1999; Fox et al. 2012; Demirtas-Tatlidede et al. 2013). However, such interventions may disrupt tightly anti-correlated processes occurring in different large-scale neural networks by further de-coupling such systems. Alternatively, a careful characterization of internetwork dynamics may inform the refinement of future interventions toward re-establishing a healthy functional regime across networks. These are speculative possibilities, but prospective investigations incorporating detailed computational models of large-scale cortical network interactions may deepen our intuition for such complex dynamics (Yang et al. 2016).

Relating Findings to Diagnosis, Treatment, and Active Symptoms

As noted above, SCZ subjects differed from HCS in 3 key measures with respect to GS beta weights computed from their rs-fMRI data (also see Fig. 5): (1) SCZ showed significantly increased mean association region b_GS values relative to HCS, (2) SCZ showed significantly reduced mean sensory region b_GS values relative to HCS, and (3) network-level differences (association vs. sensory) were substantially reduced in SCZ. We cannot rule out the possibility that these diagnosis-related group-level differences could relate to treatment effects, since all SCZ patients were receiving medication, however none of the effects were significantly related to medication dosage (Figure S4). We did not discover any significant relationships between our findings and active symptoms obtained from patient PANSS scores.

Considering Inter-Subject Variability in Association Regions

Prior studies in healthy subjects (Mueller et al. 2013) have observed greater inter-subject variability in association regions (compared with non-association regions) with respect to functional connectivity. We report a similar finding of increased inter-subject variability in association regions in healthy subjects, with respect to a different measure: network-averaged voxel-wise GS beta weights, which impact functional connectivity measures but are not themselves measures of pair-wise connectivity between regions in the sense that a Pearson's correlation in a connectivity matrix would be. This finding motivates several questions for potential future study: what mechanism, if any, could plausibly link increased association region variability in functional connectivity to the same increased variability in network-averaged GS beta weights? Why is this elevation in between-subject variability in mean association b_GS values (compared with mean sensory b_GS values) seemingly absent in the SCZ cohort? Finally, if the GS topography is nonuniform in healthy subjects and differently distributed in SCZ, could GSR during preprocessing affect signals in some brain regions more than others?

Limitations

Some important limitations need to be considered. The combined SCZ sample, although large and well characterized, was associated with comorbid history of drug/alcohol abuse/dependence to provide more generalizability to patients typically encountered in the population. In addition, it is important to establish if treatments can attenuate observed effects. We could not address this question in our study. While present control analyses argue against medication confounds, medication remains a vital consideration for future studies. Notably, patients were medicated. Although medication dose did not alter effects statistically, long-term medication might impact patterns of network-level cortical activity independently of the illness. Thus, it will be critical to replicate findings across un-medicated, prodromal, at-risk, or first-degree relatives of patients. Because of correlational measures, it is unclear whether changes reflect the cause or the consequence of the illness, which may be associated with dynamical circuit alterations over time. To disambiguate these causal possibilities, it will be important to determine if observed disturbances relate to illness duration, number of psychotic episodes, and occur in at-risk populations. We observed no significant relationships between head motion, medication levels, smoking and GS beta weights in association and sensory networks (Figure S4).

A final set of limitations relates to use of single-band (SB) BOLD acquisition and complementary volume-based registration and processing methods that limit spatial and temporal resolution (Glasser et al. 2013; Robinson et al. 2014). This leads to extensive mixing across tissue types and brain areas and makes the spatial patterns interpretable only at large scales (e.g. the broad sensory vs. cognitive network distinctions). Thus, these findings should be replicated and characterized via high spatial resolution functional multi-band (MB) BOLD data processed using cutting-edge surface-based registration (Glasser et al. 2013; Robinson et al. 2014). Also, future studies using faster temporal resolution MB data will permit more sophisticated de-noising strategies (e.g. data-driven approaches like ICA+FIX), which will facilitate better detection and elimination of motion and physiological artifact. Relatedly, this study could not separate sources of GS gray matter variance. For instance, several physiological parameters explain a substantial fraction of the GS. These include differences in end tidal pCO2 concentration, breathing rate/depth, and heart rate variations (Golestani et al. 2015). Put differently, due to the relatively low temporal resolution of our data, it is likely that heart beat effects are aliased into the BOLD frequency range. This is an issue commonly affecting fMRI studies with slower TRs (Cordes et al. 2014). As noted, this may be addressable via faster TR sampling available using multiband acquisition that permits sub-second TRs. Therefore, our present data cannot be used to address the question of whether the observed GS topography differences arise from neurological versus physiological sources. This is an important consideration, as we do not definitively know to what extent differences in GS magnitude, or spatial distribution, reflect physiological or neural differences between SCZ and HCS. Addressing this issue will critically inform further refinement of reported effects and their viability as a disease marker.

Despite these caveats, this study involved a well-powered large sample collected via “legacy” SB approaches, which permits a broad neuroanatomical and functional analysis of the GS pattern in SCZ. Furthermore, our study replicated effects across 2 samples from different sites, making it less likely that effects were driven by systematically consistent physiological artifacts across both samples. Another promising aspect of our findings is that they are present bilaterally and dissociate across association and sensory regions, in contrast to findings in studies evaluating physiological artifacts in BOLD signals across the brain (Cordes et al. 2014; Jennings et al. 2016). These studies largely seem to find effects in both association and sensory regions without showing dissociation between them (Cordes et al. 2014). Nevertheless, a more careful spatial and temporal characterization of this effect, along with simultaneous monitoring of heart beat and respiration, would be vital for future work.