From correlation to communication: Disentangling hidden factors from functional connectivity changes

The brain is a dynamic system that integrates information in order to fulfill its functions. This system is characterized by both short- and long-range connections between regions that are organized into large-scale brain networks (Fox et al., 2005; Power et al., 2011; Yeo et al., 2011). The coupling in activity between different brain regions is referred to as functional connectivity (FC). FC is typically measured by indices of statistical dependence, such as correlation coefficients, coherence, or transfer entropy (Baccalá & Sameshima, 2001; Patel et al., 2006; Vicente et al., 2011). In human functional magnetic resonance imaging (fMRI), the Pearson correlation derived from blood oxygen level–dependent (BOLD) time series is the dominant means of quantifying FC and drawing inferences on neural coupling (Biswal et al., 1995; David et al., 2004; Fiecas et al., 2013; Mahadevan et al., 2021).

While functional connectivity patterns are largely stable over time within an individual (Gordon et al., 2017; Gratton et al., 2018; Laumann et al., 2015), there is evidence that task context induces subtle changes in FC (Arbabshirani et al., 2013; Cole et al., 2014; Gonzalez-Castillo & Bandettini, 2018; Gratton et al., 2016, 2018; Krienen et al., 2014; Wu et al., 2021). Changes in FC between rest and task can help us understand how functional networks support human cognition (Allen et al., 2014; Chang & Glover, 2010; Kietzmann et al., 2019), but its interpretation remains controversial (Hutchison et al., 2013; Laumann & Snyder, 2021). Some authors have raised concerns about interpreting change in correlation as change in neural coupling (Behseta et al., 2009; Friston, 2011; Laumann et al., 2017). For example, Duff et al. (2018) suggested that many changes in BOLD correlation can be attributable to “additive signal change,” where a separate signal is added to one or both neighboring nodes. Additive signal changes include changes in noise levels or the amplitude of a common signal driving correlation between two nodes. By making an inference from the relationship between changes in variances and covariances, the authors distinguished additive signal changes from changes involving more complex combination of increases and decreases in the strength of existing signal components. A conjunction analysis incorporating variances and covariances was also proposed by Cole et al. (2016) to investigate the shared and unshared signals in BOLD data.

In the current project, we distinguish two types of factors driving the change in BOLD FC (Figure 1). The first type involves a change in the local coupling between two nodes. Local coupling can be altered by increased synaptic strength (e.g., Yao et al., 2007; Zucker & Regehr, 2002) in the cortical circuit. Local coupling can also be instantiated when the oscillation in a neuronal group is entrained by the oscillation from another group (Fries, 2005). Local coupling is closely related to the causal influence¹ (Friston, 2011) between two regions, regardless of anatomical distance.

The second type of factor involves input from the rest of the network that modulates the FC of an edge even in the absence of causal interaction between two neighbors. We differentiate two subtypes here. The first subtype, shared input, is comprised of neural signals that simultaneously affect both neighbors. For example, increases in the correlation between visual regions during visual stimulation could be explained by shared activity induced by the same stimulus. There is also considerable evidence showing that the ascending arousal system plays a role in modulating functional network topography globally (Guedj et al., 2017; Shafiei et al., 2019; Shine, 2019; Shine et al., 2016, 2018; Turchi et al., 2018). Projections from critical neuromodulatory sources like the dorsal raphe, locus coeruleus, ventral tegmental area, pedunculopontine nucleus, and the cholinergic basal forebrain are widespread, thus making these regions potential sources of shared input (Aston-Jones & Cohen, 2005; Mena-Segovia & Bolam, 2017; Morales & Margolis, 2017; Selden et al., 1998; Shine, 2019). The second subtype, unshared input, involves neural signals that are idiosyncratic to each node of a pair. For example, signals from different modalities (e.g., visual and auditory) can project to separate, but coupled regions in the higher order cortex (Sepulcre et al., 2012). Unshared input can also include local neural noises (Daunizeau et al., 2012) or measurement errors. Regardless of the mechanism, the two subtypes of input are exogenous to the causal interaction between two neighbors, and we collectively refer to them as (nonlocal) network input.

While neither network input, nor local coupling can be observed directly with fMRI, we aim to estimate their contributions to task-induced FC change by incorporating information regarding the temporal variance of the BOLD signal. Specifically, we focus on the FC changes under “task state” that occur over a block (multiple minute period). Like correlation, BOLD signal variance can change over contexts. For example, resting-state BOLD data has been reported to have a higher variance than task BOLD data (Ito et al., 2020). However, for a pair of neighbors, the shared and unshared input variance can have opposing effects. Shared signals are expected to increase the correlation between the pair while unshared signals to attenuate the correlation. We aim to quantify this relationship in both simulation and real data by estimating the network input variance from BOLD signals.

Our simulation not only highlights the multiple factors modulating FC but also provides insight on a new measure, communication change, that may be used to track local coupling more closely than FC. Communication change utilizes a learned relationship between FC and input variances during a baseline (or resting state) to mitigate the impact of input fluctuation on task-induced FC change. However, unlike correlation or partial correlation (see below) that can be used as a “standalone” metric for functional connectivity, communication change is a method to disentangle factors underlying FC changes across task contexts.²

Previous work has used partial correlation (Marrelec et al., 2006) to disentangle the impact of local coupling from the mutual dependencies on other brain regions. Partial correlation calculates the correlation between a pair of time series after the portion of their variance explained by all other observed time series is removed. Thus, partial correlation can capture a similar concept as communication change, but with a different modeling and computational approach. To evaluate how partial correlation is affected by the changes in local coupling and network input, we analyzed it in the same simulation framework.

When examining empirical data, we found that network input account for a significant but moderate amount of task-induced FC change. To estimate FC and input variances during task performance, we used BOLD signal residuals after removing evoked task effects following the “background connectivity” approach (Al-Aidroos et al., 2012; Fair et al., 2007). Depending on the task context, communication change sometimes suggests local coupling changes in the opposite direction as FC change would indicate. To gain insights into the informational content of communication change, we investigated whether it can reliably discriminate task type.

We simulated two interacting neural regions with a two-dimensional Wilson–Cowan dynamic model (Wilson & Cowan, 1972.) to illustrate that (1) FC change is sensitive to local coupling as well as to network input and (2) the impact of network input is dampened in the communication change measure where the contribution from the rest of the network, proxied by shared and unshared input variance (see section Communication Change in Methods), is removed from FC changes. BOLD signals were computed from simulated neural activity via the Balloon–Windkessel model (Friston et al., 2003). FC and communication changes were derived from the simulated BOLD series. Throughout, we refer to the temporal correlation between nodes in this model as “FC,” in equivalence to BOLD correlation.³

With all parameters (including the local coupling weight and shared/unshared neural signal variances) fixed, we computed baseline FC, baseline input variance for each node, and regression coefficients of FC with respect to the input variance. These values will be used to compute communication change in the following scenarios: (1) after varying local coupling strength, (2) after varying shared neural signal variance, and (3) after varying unshared neural signal variance.

As expected, simulated FC is sensitive to both changes in local coupling and changes in input signal variance, confounding the interpretation of FC measures. Both stronger local coupling (Figure 2A) and amplified shared signal (Figure 2B) lead to higher FC, while amplified unshared signals (Figure 2C) lowered the FC. In comparison, communication change shows a dampened relationship to shared and unshared signal variance while remaining sensitive to the local coupling (Figure 2A–C). For example, with the current simulation setup, both FC and communication were equally sensitive to local coupling. However, FC was 5 times more sensitive to shared signal (the slope in Figure 2B) and 2.5 times more sensitive to unshared signal (the slope in Figure 2C) relative to communication.

Given that partial correlation is often used to account for the mutual dependence on other brain regions, we computed partial correlation in the same three scenarios (see Methods). We asked: Does partial correlation have less dependence on input signal variance than FC? In the current simulation framework, partial correlation had similar sensitivity as FC with respect to changes in local coupling and input signals (Figure 2D–F). Therefore, it cannot discriminate among those factors.

This analysis suggests that communication, compared to correlation and partial correlation, more closely tracks changes in local coupling relative to changes in input variance. Communication diminishes the impact of network input by using the statistical relationship between FC and variance in a baseline state. Notably, this is an imperfect proxy because input signals are unobservable in real data. Therefore, communication is not guaranteed to have zero sensitivity to shared and unshared input variances, but their impact is considerably dampened as observed from the simulation. This property allows us to infer the change in local coupling from real data when neither the coupling nor input signals are directly observable.

To compare communication and FC in empirical data, we used Midnight Scan Club (MSC) dataset (Gordon et al., 2017), a precision fMRI dataset with 10 sessions of fMRI rest and task data from nine participants. See Methods for BOLD data acquisition, preprocessing, FC, and variance calculation for resting-state and task data. We computed BOLD FC, shared input variance, and unshared input variance for each of the 10 resting-state scans. We fit a linear model for each edge with the FC as the dependent variable and shared/unshared input variance as independent variables. As expected, FC increased with shared input variance and decreased with unshared input variance in the MSC data. Figure 3 shows the distribution of the regression coefficients of all participants. Network input fluctuations jointly explained 24.7 ± 0.1% (adj. R²) of the FC variability across resting-state sessions. Results from individual participants all show a similar pattern and are included in Supporting Information Figure S1, indicating robustness of the results.

Thus, the empirical relationship between FC and network input is consistent with what we observed in the simulation: FC increases with the shared input variance and decreases with unshared variance idiosyncratic to one of the nodes. We next use the regression coefficients learned from resting state to estimate the contribution of input variance in FC during task-induced contexts.

Since the previous simulations suggest that communication tracks local coupling while dampening the impact of input fluctuation, we next directly compare FC and communication measures in real data. We compared the task-induced changes in each measure, grouping edges according to a widely used large-scale network designation (Gordon et al., 2016). Communication and FC changes during tasks differed in many networks and these differences were modulated by the task type (Figure 5). Sometimes communication was higher magnitude than FC, and other times it was lower; on occasion it even showed relationships in the reverse direction. For example, when participants performed the categorization task⁴ (Figure 5A), FC decreased on average within the ventral attention network, while communication (calculated from FC after removing the contribution from input variance, see Methods, Equation 2) increased. This highlights the danger of interpreting an FC decrease directly as reduced interregional interaction.

To probe the informational content in the communication measure, we ask to what extent it can discriminate among tasks with distinct task demands. To this end, we computed the intraclass correlation coefficient (ICC; see Methods) between task types for each edge across all scan sessions (4 task types including the resting state with 10 sessions each). ICC measures the percentage of variance accounted for by task types among all changes observed for the same edge (Bartko, 1976; McGraw & Wong, 1996). If communication change contains task information, ICC⁵ is expected to be closer to 1 (rather than to 0).

Across all edges and all participants, the average communication ICC = 0.699 ± 0.001, with 71% of the edges showing a a p value less than 0.05 (F test), indicating that most edges contain reliable information about task type (Figure 6). For comparison, we computed the ICC for FC changes using the same methodology. Similar to communication, FC changes reliably represent task information (ICC = 0.462 ± 0.001; Figure 6). However, the ICC based on the communication measure is significantly higher than the ICC of FC (z = 140, p < .00001), suggesting that communication change increases the discrimination of task type.

The same pattern was observed consistently for all participants (Figure S4). The difference between FC and communication in terms of task discrimination is not solely driven by the distinction between the resting state and any task state. ICC among the three task states exhibits a similar pattern (although the difference is narrower, z = 77, p < .00001; see Figure S5) indicating that communication change contains more task-specific information relative to FC.

In summary, communication changes during tasks exhibit significant differences from FC in various networks of the brain. Importantly, communication contains more informational content to discriminate specific task types, compared to FC.

While BOLD correlation is widely used to capture both static functional connectivity and its changes across contexts, the interpretation of this measure is often ambiguous. The entanglement of multiple latent factors including local coupling of two neighbors and network inputs from the rest of the network limits the scope of the conclusions that can be drawn from correlation measures alone. Here, we developed a method to compute the shared and unshared input variance to an edge and quantify their impact on task-induced FC changes. We further proposed a metric, termed “communication change,” that dampens the impact of input signals and thereby provides a promising candidate for more closely tracking local coupling changes.

We first analyzed local and network impact on FC by conducting a simulation based on a well-established model of neural dynamics (Wilson & Cowan, 1972). The model suggests that FC can be driven by multiple factors. Increased local coupling as well as increased variance of shared inputs result in higher FC, while increases in the variance of unshared (idiosyncratic) input has the opposite effect. Although the conclusion seems intuitive, explicit modeling of the three latent factors leads to a novel metric, communication, that dampens the impact of both shared and unshared inputs using only observable, statistical relationships. Our simulation shows that communication change is more likely to track local coupling.

Interestingly, partial correlation, a metric often used to disentangle the impact of local coupling from the mutual dependencies on other brain regions (Marrelec et al., 2006; Smith, 2012), shows similar sensitivity as FC change in the current simulation framework. Partial correlation removes the contemporaneous signal dependence by using a linear residualization. However, as the BOLD signal exhibits temporal nonlinearity, the transformation from neural activity to BOLD signal (via Ballon–Windkessel model) is state and path dependent (Buxton et al., 2004; Friston et al., 2000). It is possible, therefore, that the correlation between residualized signals is still affected by the network input through their “history.” We are not aware of other simulative studies illustrating the (lack of) efficacy of partial correlation, except for Cole et al. (2016), where the authors also reported that partial correlation retains limitations of the FC measure. In comparison, to compute communication change, we first aggregated inputs from the rest of the network into three terms, shared variance and two unshared variances idiosyncratic to each node. Communication then removes the input contribution via a “learned” statistical relationship between FC and variance at baseline (here operationalized as the average resting-state signal). Therefore, the method makes few assumptions about the underlying structure of the input signals.

With empirical data, we computed the sensitivity of FC with respect to shared and unshared input during resting-state scans. As expected, FC in real BOLD data increased with the variance of shared input and decreased with unshared input. We then analyzed task-induced FC changes using the coefficients learned from the resting state. Input fluctuation explained a significant, but moderate, amount of FC changes during task performance. After removing the input related contributions, the remaining FC change, that is, communication change, often had the opposite sign of the original FC change in various tasks and brain networks. This highlights the ambiguity of interpreting task-induced FC change by itself. If the correlation of an edge increases, for example, it can be driven by an increase in shared input from the rest of the network, a decrease in an idiosyncratic input affecting only one of the nodes, an enhanced neural coupling between the nodes, or any combination of the above factors.

Although we cannot directly observe how communication differentiates mechanisms behind FC changes in empirical data, evidence suggests that it shows better discrimination across the task types relative to the standard FC measure. Both simulation and empirical analysis suggest that communication may be an improved measure (from FC) for tracking changes in causal interactions between regions across task contexts. Such changes can occur at the synaptic level due to short-term plasticity (Krystal et al., 2003; Yao et al., 2007; Zucker & Regehr, 2002), such as an elevated dopamine tone in the same cortical circuit (Vijayraghavan et al., 2007) and via increased coherence or phase locking (Fries, 2005).

While FC is strongly related to shared and unshared network inputs, communication shows a diminished relationship to these latent signals. Shared input can arise from signals with widespread influence stemming from subcortical projections (Raut et al., 2021; Shine, 2019), while unshared (idiosyncratic) input can come from, for example, unique inputs from different modalities (visual or somatosensory) that only project to selected regions. The communication metric dampens the impact of network input changes using a statistical relationship and remains agnostic to the underlying biophysiological mechanism. This allows contextual factors to guide the interpretation of the difference between FC and communication across a variety of states.

Other researchers have also explored the idea of incorporating variance (or covariance) to disambiguate the interpretation of FC (Duff et al., 2018; Cole et al., 2016). For example, Duff et al. (2018) combined correlation and variance to classify FC changes into cases attributable to “additive signal change” (a type of “shared or unshared input”) and changes arising from modifying existing signals. While the current analysis echoes the importance of augmenting FC with variance/covariance, we propose a framework to quantitatively decompose FC changes into contributions from different underlying mechanisms.

Communication shows higher intraclass correlation across tasks compared to FC, suggesting that the proposed measure better discriminates task types. This suggests that local coupling is modulated by specific task demands, while overall signal variance contains less task-specific information. Previous work has shown a simultaneous task-induced reduction in FC and variance (Ito et al., 2020). It is possible that task performance, relative to the resting state, reduces synchronized, spontaneous activity that affects a large region of the brain and leads to an FC reduction in a pattern that is similar across task types. By extracting the “excess” FC change after accounting for the variance-related effects, we can gain more insight into how specific task demand modulates the brain network.

Furthermore, the current framework provides a way to decompose exogenous network inputs into shared and unshared variance. One powerful approach may be to use several of these metrics in combination. For example, an increase in communication implies stronger local coupling. If a decrease in FC is simultaneously observed, we can ask whether the discrepancy is due to a reduction of shared signal or an increase of idiosyncratic inputs by examining changes in shared and unshared variance simultaneously. Therefore, examining communication together with FC may allow researchers to draw inferences regarding how a particular region adjusts its interactions with its neighbors across a range of different contexts, providing an improved handle on whether local or shared inputs are at work in any given scenario.

One pitfall of our approach is that the interpretation of our simulation results is dependent on the validity of our computational model. The communication change metric was inspired by a simplified neural dynamic simulation. Although Wilson–Cowan dynamics have been widely used and have had a major impact on computational neuroscience, it is a simplification of the underlying biological process (Chow & Karimipanah, 2020). Many extensions and generalizations have been proposed since its initial introduction in 1972. To highlight the basic interaction between several statistical measures of interest, we adopted the simplest form of the model in a two interacting nodes setup. By doing this, we imposed the assumption that any effect from the rest of the system could be represented by Gaussian signals. Furthermore, we fixed the mean of the input signals and focus on the impact of their variance/amplitude. Because the neural dynamics involve sigmoid activation functions, the shift in the average level of the input signal would have a nonlinear impact on the resulting fMRI BOLD signals (Ito et al., 2020). Future work will explore how the average activity level interacts with the relationship between FC and input variances discussed in this article.

When it comes to empirical data, a key premise of our analysis is based on evidence that fMRI signals are indirectly related to underlying neural population activity. Such a view is widely supported in the literature (Logothetis et al., 2001; Ma et al., 2016), but certain MRI artifacts can ambiguate the results (Friston et al., 1996). We employed multiple artifact mitigation measures (scrubbing, motion regression, global signal regression, white matter signal regression, and ventricle signal regression) to minimize the impact of artifactual factors.

Moreover, the BOLD signal analyses utilized nine highly sampled subjects from the MSC dataset (Gordon et al., 2016). This is a small sample, but there is evidence that results from MSC data extend to other independent datasets (Gordon et al., 2018). Importantly, the single subject results presented here are consistently observed in all individuals, demonstrating a certain degree of robustness. However, communication change should be examined in other datasets in future work, especially with task conditions that impose varying degrees of cognitive demand and a larger number of subjects.

There are unanswered questions that are beyond the scope of the current investigation that pertains to communication change. Task-induced communication increases in some regions (e.g., visual network, dorsal attention network), and the magnitude of the changes seems to be modulated by the type of task. In addition, the spatial pattern of communication exhibits significant individual differences (Figure S3), consistent with findings from previous work (Gratton et al., 2018; Porter et al., 2023). An intriguing question for future research is how differential communication is related to task performance.

Finally, task discrimination is only one possible venue to explore the informational content in communication and its advantage over FC. Research directly modifying potentially widespread signals such as arousal (shared input), signals that might be expected to be more local such as selective stimulus-response features (unshared input), or manipulation at the local synaptic level (coupling) may help to validate the interpretation of communication change.

FC, computed as a correlation between paired BOLD signals, is not only determined by local coupling between pairs of regions but also by inputs from the rest of the network, among other factors. Utilizing the statistical relationship between FC and input variance, we developed a method to extract a component (communication change) from FC change with a dampened influence of network input. Examining changes in communication, along with FC, allows us to draw inferences regarding how the interaction between a particular region with its neighbors changes across different contexts.

We used a combination of numerical simulations and empirical analysis to examine the effects of local interactions and input fluctuations on BOLD correlation and variance. These results led us to develop a new measure, communication change, that more closely tracks the change in coupling between two neural regions compared to the standard FC measure (BOLD correlation) across task contexts. Below, we first describe our modeling approach, then the derivation of the communication metric, and finally we include a description of a dataset and analyses used to empirically examine differences in communication change and correlation.

We developed a novel metric to better capture the interregional coupling shift in task-induced FC change, which we termed communication change. To compute the communication change of an edge, we first need to estimate the shared and unshared network input variance.

Shared/unshared input variance is computed from the variance of the BOLD signals from common and idiosyncratic neighbors. Specifically, assume node1 and node2 are the nodes of the edge under consideration. N_s is the set of neighboring nodes that connect to both node1 and node2. N₁ is the set of nodes that only connect to node1, and N₂ is the set of nodes that only connect to node2. That is, Ns, N₁, and N₂ are mutually exclusive. Shared input variance is the total variance of the BOLD signals from N_s, computed by the trace of the covariance matrix. We residualized all the signals in N₁ against all signals in N_s. Unshared input varaince₁ is the total variance of the residualized signals in N₁. Likewise, unshared input varaince₂ is the total variance of the residualized signals in N₂ against N_s.

All edges and neighbors in empirical data (including rest and task) are defined by 5% sparsity based on the average resting-state correlation. We reevaluate the key findings using 10% sparsity and reached similar conclusions (Figure S6).

The matrix W = $w11w12w21w22$ describes the information exchange within and between regions. In particular, the off-diagonal terms, w₂₁, w₁₂, determine the strength of interregion coupling. For the baseline scenario, the matrix was set to W = $2112$⁠. The diagonal entries were fixed throughout the simulation. The parameter values (and those in the following) are adapted from Ito et al. (2020) and adjusted to obtain a realistic baseline correlation.

Shared input s(t) is a Gaussian signal with shared signal variance V_sh, which was set to V_sh = 0.25 in the baseline scenario. V_sh indexes the magnitude of the common exogenous input from the rest of the brain network to both regions.

I₁(t), I₂(t) are two independent Gaussian signals representing idiosyncratic signals with unshared signal variance parameter V_un, which was set to 0.35 in the baseline. V_un indexes the magnitude of the idiosyncratic input to one of the neighbors.

f is the sigmoid input-output activation defined as f(x) = 1/(1 + e^(−k * x)), where k is fixed at 0.5 throughout the simulation.

Neural activities simulated from the above equations were transformed to BOLD signals using the Balloon–Windkessel model (Friston et al., 2003). The transformation assumes normalized deoxyhemoglobin content, normalized blood inflow, resting oxygen extraction fraction, and normalized blood volume. State equations and parameters were taken from previous work (Friston et al., 2003). Correlation (FC), partial correlation, and communication change were then computed on the simulated BOLD signals as specified below. Each simulation (1 run) generated a pair of BOLD time series with 2,000 steps each.

For the baseline simulation, we ran the scenario with fixed parameters 500 times (500 runs). We obtained the average FC, average shared/unshared input variance, and the regression coefficients as in Equation 3. The baseline simulation mimics the resting state as in real data.

Next, we computed the sensitivity of correlation, communication, and partial correlation with respect to the input signal variance and local coupling parameters. To this end, we carried out the same simulation as above, but now independently varied the shared signal variance V_sh, unshared signal variance V_un, and the coupling parameter w₁₂ (and symmetrically w₂₁). We used inputs for V_sh ranging from 0.15 to 0.35, V_un ranging from 0.25 to 0.45, and w₁₂ from 0.4 to 1.6 (centered around baseline 1). The range was chosen such that the parameters had comparable impact on FC. For each scenario, the model is run 500 times to obtain the mean and standard deviation of all quantities of interest.

Correlation (FC) was then computed between a pair of simulated BOLD timeseries.

Communication change was computed as described above. All the variances required in the formula were computed from the simulated BOLD signals. To simulate the input BOLD signals, we applied the Balloon–Windkessel model to transform the shared signal s(t) and the idiosyncratic signals I₁(t), I₂(t).

Partial correlation was computed as the correlation between a pair of residualized series, after removing the common variance from the rest of the network. In the current simulation framework, inputs from the rest of the network are represented by the shared/unshared signals. Therefore, we residualized the BOLD signals derived from x₁(t), x₁(t), against the input BOLD signals derived from I₁(t), I₂(t), and s(t), and computed the correlation between the residualized signals.

In the second portion of our investigation, we examined the relationship between FC, variance, and communication in an empirical dataset. We used the Midnight Scan Club (MSC), a precision fMRI dataset (Gordon et al., 2017) that contains large quantities of both task and resting-state data was selected, as large amounts of data are required to achieve highly reliable FC measures (Gordon et al., 2017; Laumann et al., 2015). The dataset includes data from 10 participants who each participated in 10 separate fMRI sessions. During each session, each participant completed three different tasks (a motor task, a categorization task, and a mixed dot coherence/verbal discrimination task) and rest. This dataset is described in detail elsewhere (see Gordon et al., 2017), but we will review parameters relevant to this study below. In addition, a dataset of 120 healthy adults (WashU-120) was used as a reference set for the network definition and the calculation of graph theoretic metrics (Power et al., 2011; Power et al., 2013). This dataset has been described in detail in Power et al. (2013).

Each resting-state session was 30 min, during which participants were asked to lie still while fixating on a white cross presented against a black background. Each motor task session included two runs (7.8 min total) of a blocked motor task adapted from the Human Connectome Project (Barch et al., 2013). In each block (15.4 s in duration), participants were cued to move either their left or right hand, left or right foot, or tongue. Each run included two blocks of each type of movement, as well as three fixation blocks (15.4 s). Each mixed task session included two runs (14.2 min total) that each included four blocks, two of the semantic tasks and two of the coherence tasks. Both tasks had a mixed block/event-related design, modeled after the tasks in Dubis et al. (2016). During the semantic task, 30 individual trials were presented, consisting of words presented for 0.5 s with jittered 1.7–8.3 s intervals. Participants were asked to respond whether the words were nouns or verbs (50% nouns and 50% verbs were included). Forty-four second fixation periods separated blocks. During the coherence task (also 30 trials), individual trials consisted of arrays of Glass patterns (Glass, 1969) that were varied in how concentrically they were arranged (50% or 0% coherence to a concentric arrangement, displayed with equal frequency). Participants were asked to respond whether dots were arranged concentrically or randomly. Each categorization task session included three runs of an event-related categorization task with implicit repetitions, with a separate run per stimulus type (face, scenes, and words). Within each run, participants viewed 24 images, repeated three times and made binary decisions.

High-resolution T1-weighted (224 slices, 0.8 mm³ isotropic resolution, TE = 3.74 ms, TR = 2,400 ms, TI = 1,000 ms, flip angle = 8 degrees), T2-weighted (224 slices, 0.8 mm³ isotropic resolution, TE = 479 ms, TR = 3,200 ms) both with 0.8 isotropic resolution, and resting-state BOLD data were collected on a Siemens 3T Magnetom Tim Trio with a 12-channel head coil (Gordon et al., 2017). Functional scans were collected with a gradient-echo EPI sequence, isotropic 4-mm³ voxels, TE of 27 ms, and TR of 2,200 ms (Gordon et al., 2017). The MSC dataset is considered a deep (precision) dataset with 5 hours of resting-state data per subject in 30-min blocks over 10 separate sessions (Gordon et al., 2017).

The T1-weighted images were processed via automatic segmentation of the gray matter, white matter, and ventricles in Freesurfer 5.3 (Fischl et al., 2002). The default recon-all command in Freesurfer was then applied to produce the anatomical surface for each subject (Dale et al., 1999). These surfaces were manually edited to improve the quality of the registration and the surfaces were registered to the fs_LR_32k surface space via the procedure outlined in Glasser et al. (2013). Functional data preprocessing began with field inhomogeneity distortion correction using the mean field map, which was applied to all sessions. Slice timing correction was implemented using sinc interpolation to account for temporal misalignment in slice acquisition time. This step was followed by motion correction that was conducted within and across BOLD runs (aligned to the first frame of the first run) through a rigid body transformation. Then whole-brain intensity values across each BOLD run were normalized to achieve a mode value of 1,000 (Miezin et al., 2000). Functional BOLD data was first registered to a T2-weighted image and then to the T1. An affine transformation was employed for registration. The T1-weighted image was aligned to a template atlas (Lancaster et al., 1995) conforming to Talairach stereotactic atlas space (Talairach & Tournoux, 1988) using an affine transformation. All computed transformations and resampling to 3-mm isotropic voxels were simultaneously applied at the end of these steps.

Steps were taken to attenuate the impact of artifacts on BOLD time series. The effect of nuisance signals was mitigated through regression of average signal from the white matter, ventricles, whole brain (global signal), motion parameters, and derivatives and expansion terms of the motion parameters (Friston, 1998; Power et al., 2014). Removal of frames with FD > 0.2 mm, in addition to sequences containing less than five contiguous low-motion frames, the first 30-s of each run, and runs with < 50 low-motion frames was performed as an additional safeguard against the impact of motion (Power et al., 2014). A special filtering procedure was applied to two MSC subjects (MSC03 and MSC10) with respiratory contamination in their motion parameters. This procedure involved applying a low-pass filtered at 0.1 Hz to reduce the effects of respiratory artifacts on motion estimates stemming from the short-TR multiband acquisition, all before censoring high-motion frames (Fair et al., 2020; Gratton et al., 2018; Laumann et al., 2017). Then a filtered FD threshold of 0.1 mm was applied to censor frames. In all cases, flagged head motion frames were removed and the time points were replaced with interpolated data using a power-spectral matched approach (Power et al., 2014), and this was followed by the application of a bandpass filter (0.009–0.08 Hz). Previous research has shown that ∼30 min of low-motion data is necessary to achieve high reliability of FC (Gordon et al., 2017). We therefore removed a participant (MSC08) that exhibited high degree of motion (Gordon et al., 2017). The excluded participant also reported a high level of drowsiness (Gordon et al., 2017; Laumann et al., 2017; Tagliazucchi & Laufs, 2014).

The processed BOLD data were mapped to each individual’s native midthickness surface via the ribbon-constrained sampling procedure (Marcus et al., 2013) Then, the mapped data were registered to the fsaverage surface in one step using the deformation map generated from the ribbon-constrained sampling procedure described in Glasser et al. (2013). This was followed by spatial smoothing conducted via a geodesic Gaussian smoothing kernel to the surface registered data (FWHM = 6 mm, sigma = 2.55) (Gordon et al., 2016; Marcus et al., 2011). Next, temporally interpolated frames were removed (prior to functional connectivity analysis).

We examined a previously published set of 333 brain parcels associated with 12 networks somatomotor dorsal, somatomotor lateral, visual, auditory, cinguloopercular, fronto-parietal, dorsal attention, ventral attention, salience, default, parietal memory, and retrosplenial (Gordon et al., 2016).

In the case of the task scans, task related activation was captured with a general linear model (GLM) using in-house software written in IDL (Research Systems, Inc.), as described in Gratton et al. (2018), using a finite impulse response approach for trial-level regressors (Ollinger et al., 2001). Functional connectivity was then measured for the task scans following the “background connectivity” method (Al-Aidroos et al., 2012; Fair et al., 2007) with the residuals from the GLM analysis being used for BOLD time series correlations. The residuals were subjected to the aforementioned connectivity processing pipeline prior to functional connectivity calculation. There are multiple methods to measure task FC, that here we focus on background connectivity because of its simplicity and ability to track task changes at the resolution we are interested in, and that it is well suited for comparison with a resting baseline.

We used the Python package Pingouin (Vallat, 2018) to evaluate how much variability in FC and communication across scans is attributable to task type. Specifically, we used the function intraclass_corr, where the “targets” are the task types (rest, mixed, categorization, motor), raters are the session index (1–10). We used “ICC1k” output indicating that the target (task type) is rated by different raters (the session) and the raters are selected at random. This is equivalent to the one-way ANOVA fixed effects model (Shrout & Fleiss, 1979).

Derek M. Smith would like to thank The Therapeutic Cognitive Neuroscience Fund.

Supporting information for this article is available at https://doi.org/10.1162/netn_a_00290. Midnight Scan Club data is publicly available (https://openneuro.org/datasets/ds000224). Code related to the analysis in this paper is located at https://github.com/DijunQuant/communication_public/.

Yuhua Yu: Conceptualization; Formal analysis; Investigation; Methodology; Visualization; Writing – original draft; Writing – review & editing. Caterina Gratton: Funding acquisition; Supervision; Writing – review & editing. Derek M. Smith: Conceptualization; Data curation; Funding acquisition; Methodology; Writing – original draft; Writing – review & editing.

Derek M. Smith, National Institutes of Health (https://dx.doi.org/10.13039/100000002), Award ID: T32NS047987. Caterina Gratton, National Science Foundation (https://dx.doi.org/10.13039/100000001), Award ID: CAREER2048066. Caterina Gratton, National Institutes of Health (https://dx.doi.org/10.13039/100000002), Award ID: R01MH118370.

Local coupling:

The causal interaction between two neural regions that can be modulated by influences such as the synaptic strength or entrained oscillation.

Shared input:

Neural signals that simultaneously affect two coupled network regions.

Unshared input:

Neural signals that are idiosyncratic to one of the coupled network regions.

Communication change:

The part of task-induced functional connectivity change that cannot be attributed to the change of shared and unshared input.

Partial correlation:

The correlation between a pair of time series after the portion of their variance explained by all other observed time series is removed.

Intraclass correlation coefficient:

The percentage of variance accounted for by the group difference out of the overall population variance.