Brain connectivity at rest predicts individual differences in normative activity during movie watching.

When exposed to the same sensory event, some individuals are bound to have less typical experiences than others. Previous research has investigated this phenomenon by showing that the typicality of one's sensory experience is associated with the typicality of their stimulus-evoked brain activity (as measured by intersubject correlation, or ISC). Individual differences in ISC have recently been attributed to variability in focal neural processing. However, the extent to which these differences reflect purely intra-regional variability versus variation in the brain's baseline ability to transmit information between regions has yet to be established. Here, we show that an individual's degree and spatial distribution of ISC are closely related to their brain's functional organization at rest. Using resting state and movie watching fMRI data from the Human Connectome Project, we reveal that resting state functional connectivity (RSFC) profiles can be used to predict cortex-wide ISC. Similar region-level analyses demonstrate that the levels of ISC exhibited by brain regions during movie watching are associated with their connectivity to other regions at rest, and that the nature of these connectivity-activity relationships varies as a function of regional roles in sensory information processing. Finally, we show that an individual's unique spatial distribution of ISC, independent of its magnitude, is also related to their RSFC profile. These findings contextualize reports of localized individual differences in ISC as potentially reflecting larger, network-level alterations in resting brain function and detail how the brain's ability to process complex sensory information is linked to its baseline functional organization.

Introduction

The ability to interpret information from the outside world in a manner similar to one’s conspecifics is critical to healthy human behavior. Recent work on the neural substrates of sensory information processing has accordingly focused on the processing of naturalistic stimuli that reflect the complex and continuous nature of everyday experiences. One key finding from this line of research is that rich, time-locked stimuli such as movies elicit strongly conserved brain responses across individuals (Hasson et al., 2004). The extent of this stimulus-evoked brain activity is most frequently indexed by intersubject correlation (ISC), a measure of how similar one’s activity time course in a given brain region is to that of another individual or group average (Hasson et al., 2004; Nastase et al., 2019).

Early research using naturalistic stimuli focused on aspects of stimulus-evoked activity that are shared across subjects, but more recent work has sought to understand individual differences in ISC, or why some individuals exhibit less typical brain responses to naturalistic stimuli than others (Campbell et al., 2015; Finn et al., 2018; Finn et al., 2020; Gruskin et al., 2020; Guo et al., 2015; Salmi et al., 2013). Across these studies, atypical brain responses have on several occasions been attributed to two general causes: failure to encode stimulus-specific information (e.g., a lapse in attention causes a subject’s gaze to deviate from the visual display) and/or idiosyncratic intra-regional processing of the stimulus (e.g., a lesioned brain region processes information in a unique manner). Importantly, recent studies using demanding cognitive tasks have shown that atypical task activations reflect aberrant inter-regional communication, rather than purely intra-regional deficits (Hearne et al., 2021; Mill et al., 2020). These findings suggest that atypical brain responses to naturalistic sensory information might similarly be linked to individual differences in intrinsic functional brain organization. The aims of the present manuscript are to test this hypothesis and identify general patterns in how the brain’s intrinsic functional organization is associated with normative stimulus-evoked activity.

Just as brain responses to naturalistic stimuli reflect individual-specific variations on a population-general theme, the brain’s resting functional architecture is organized into networks whose general structures are conserved from person to person but whose fine-grained topologies and strengths serve as stable identifiers of individuals (Finn et al., 2015; Gordon et al., 2017). These networks are classically derived from resting state functional connectivity (RSFC), which is defined as the correlation between activity time courses from any two regions measured while the brain is not performing an explicit task. A growing literature has shown that RSFC can be used to predict both normal and pathological variation in model-based measures of brain activity during challenging experimental tasks, providing compelling evidence for the idea that RSFC describes the routes along which task-relevant information travels (Cole et al., 2016; Hearne et al., 2021; Ito et al., 2017; Mill et al., 2020; van den Heuvel and Hulshoff Pol, 2010). However, the extent to which RSFC can predict individual differences in model-free neural responses to naturalistic stimuli and shed light on how inter-regional communication facilitates localized processing of such stimuli is still unclear.

To characterize how RSFC is associated with normative stimulus-evoked activity during movie watching, we analyzed fMRI data acquired across multiple days from healthy Human Connectome Project participants. First, we show that RSFC can be used to predict cortex-wide ISC in new subjects and using data from held-out stimuli. Next, we demonstrate that the resting state functional connections most associated with ISC in a specific region vary systematically across the brain and quantify how greater intrinsic connectivity to one region is associated with more or less typical activity in others. Finally, we show that individual-specific RSFC patterns are related to the spatial distribution of ISC across cortex. Taken together, these results provide important context for interpreting atypical responses to naturalistic stimuli and detail new relationships between intrinsic and task-driven brain function.

Materials and methods

Participants

Data used for this project come from the Human Connectome Project (HCP) Young Adult 7T release (Van Essen et al., 2013). Of the 184 subjects who underwent 7T fMRI scanning, eight participants did not complete every resting state and movie watching run. These subjects were excluded from all analyses for a sample size of n = 176 (106 females, 70 males). All participants were healthy individuals between the ages of 22 and 36 (mean age = 29.4 years, standard deviation = 3.3) and provided informed written consent as part of their participation in the study. Self-reported racial identity in this sample was 87.5% White, 7.4% Black or African American, 4% Asian/Native Hawaiian/Other Pacific Islander, and 1% unknown/not reported, and 1.7% of the sample identified as Hispanic/Latino.

fMRI data

FMRI data for the HCP 7T release were collected at the University of Minnesota on a 7T Siemens Magnetom scanner during four sessions spread across two or three days. Each day of data collection involved two resting state and two movie-watching runs across two sessions, with one rest run and both movie runs taking place during the same scan session and the second rest run taking place during another. The same echo-planar imaging sequence was used for all rest and movie scans and its key parameters are as follows: time of repetition (TR) = 1000 ms, echo time (TE) = 22.2 ms, number of slices = 85, flip angle = 45°, spatial resolution = 1.6 mm3.

All four rest runs had duration of 900 TRs/15:00 min, whereas the four movie runs had variable durations (921, 918, 915, and 901 TRs). Subjects were instructed to fixate on a projected bright crosshair on a dark background during the rest runs and passively viewed a series of video clips (with sound) during each movie run. Movie runs one and three each consisted of four unique clips from independent films, and movie runs two and four were each composed of three unique clips from major motion pictures. The same montage of brief (1–4 s) videos was included at the end of every run for test-retest purposes. Clip durations varied between 83 and 259 s. More information on these clips can be found in Finn and Bandettini (2021) and at https://db.humanconnectome.org. Each video clip was preceded by 20 s of rest (during which the word “REST” was projected against a dark background). Any TR that took place during one of these rest blocks was discarded from all analyses. The first 20 TRs after each rest block (corresponding to the start of each video clip) were also discarded to prevent onset transients from biasing our intersubject correlation (ISC) measurements (Dosenbach et al., 2006; Fox et al., 2005), leaving four movie runs of 697, 735, 691, and 718 TRs, respectively. Finally, rest and movie runs from the same day were normalized and concatenated, yielding one rest run (1800 TRs) and one movie run (1432/1409 TRs) for each of the two days of data collection.

Preprocessing and parcellation

ICA-FIX denoised CIFTI files (e.g., rfMRI_REST1_7T_PA_Atlas_hp2000_clean.dtseries.nii) at 2 mm resolution were downloaded from ConnectomeDB. Briefly, these data were pre-processed using motion and distortion correction, high-pass temporal filtering, and MNI alignment, followed by regression of 24 motion parameters as well as a set of independent component analysis-derived confound time courses (Glasser et al., 2013; Griffanti et al., 2014). Because head motion-related and other artifacts may persist in fMRI data even after ICA-FIX, additional denoising was performed. Following recent work with HCP data (Finn and Bandettini, 2021; Li et al., 2019), the average time course of all grayordinates and its temporal derivative was regressed from each scan. Additionally, high-motion frames were censored from all rest scans. To identify these frames, a 0.2 Hz (12 breaths per minute) low pass filter was first applied to each scan’s framewise displacement (FD) trace to account for respiratory artifacts found in fMRI data (Fair et al., 2020; Gratton et al., 2020). Any volume exceeding 0.2 mm FD post-filtering was flagged, as were all runs of fewer than five contiguous volumes. Finally, delta functions corresponding to each of the censored volumes were included along with the two global signal time courses for each rest scan’s denoising design matrix. Frame censoring was not performed on the movie scans to ensure that ISC would be calculated with the same TRs across individuals.

Following denoising, functional data were parcellated into 360 cortical regions of interest using the Glasser HCP-MMP parcellation (Glasser et al., 2016). The Glasser parcellation was chosen for its anatomically specific labels, association with the HCP dataset and CIFTI format, and compatibility with the Cole-Anticevic Brain-wide Network Partition (Ji et al., 2019), whose granular and interpretable network definitions (e.g., separate language and auditory networks) are especially suitable for the analysis of BOLD signals evoked by complex audiovisual stimuli.

Resting state functional connectivity

Using data from the two concatenated rest scans, time courses from all possible pairs of the 360 parcels were (Pearson) correlated to create two symmetric 360 × 360 resting state functional connectivity (RSFC) matrices for each subject, one for each day of data collection. Frames flagged as having high motion as per Section 2.3 were excluded from the correlation calculations.

Intersubject correlation

ISC analysis was used to quantify the typicality of each individual’s BOLD responses to the Day 1 and Day 2 movie stimuli. For each of the 360 cortical parcels, each participant’s BOLD signal time course was normalized and (Pearson) correlated with that parcel’s average BOLD signal time course across all other participants to yield a 360 parcels × 176 subjects intersubject correlation matrix. This was repeated separately for both concatenated movie scans such that every participant had two independent ISC values for each of the 360 parcels, each reflecting the typicality of that individual’s BOLD responses to that day’s movie stimuli in that parcel.

Ridge connectome-based predictive modeling

Overview:

Ridge connectome-based predictive modeling (rCPM) was used to relate RSFC to global ISC (gISC), defined as the average of each participant’s 360 parcel-level ISC values. Our use of this measure as a cortex-wide indicator of brain activity typicality is motivated by its simplicity as well as by previous work showing global relationships between ISC and behavior (Gruskin et al., 2020). CPM is an established technique for predicting behavior (e.g., fluid intelligence, symptom severity) from RSFC (Finn et al., 2015; Rosenberg et al., 2016; Shen et al., 2017). Here, we use a variant of CPM, rCPM (Gao et al., 2019), to predict gISC from RSFC. Although the prediction of gISC is conceptually similar to the prediction of any other measure, gISC’s derivation from BOLD data leaves it susceptible to motion and other fMRI acquisition artifacts. More importantly, individual differences in these artifacts are likely consistent across rest and movie scans such that relationships between RSFC and gISC may be driven by artifacts rather than neural signals of interest. Therefore, in addition to the conservative preprocessing approach out-lined in Section 2.3, all correlations between RSFC and (g)ISC performed in this paper included head motion (mean FD; calculated as the average of the filtered and uncensored FD traces) and temporal signal-to-noise ratio (tSNR; calculated as the whole-brain average of the means of all grayordinate time series divided by their standard deviations for a given subject/scan) measures as covariates, as both of these variables were found to be correlated with gISC (mean FD: Day 1 ρ = −0.24, P = 1.2 × 10−3; Day 2 ρ = −0.20, P = 8.8 × 10−3; tSNR: Day 1 ρ = 0.27, P = 3.3 × 10−4; Day 2 ρ = 0.28, P = 2.1 × 10−4). Although detailed descriptions of CPM can be found elsewhere (Shen et al., 2017), provided below is a summary of the specific approach used in this paper.

Generating the cross-validated models:

Because the HCP 7T dataset is composed of data from individuals of varying degrees of genetic relatedness (monozygotic and dizygotic twins, non-twin siblings, and un-related individuals; 93 unique families), all individuals from the same family were randomly assigned to one of two groups of 88 (i.e., split-half cross-validation), with one group being used to train a model that would then be tested on the other (and vice versa). The following approach was then applied to 100 of these random splits of the data to assess the performance of rCPM across different training/testing sets and to build a bagged model that is more robust to overfitting (O’Connor et al., 2021).

Calculate leave-one-out ISC according to the method described in Section 2.4 separately for each group of 88 subjects. Re-calculation of ISC is necessary within each randomly assigned group because calculating an individual’s ISC using data from the 175 subjects across both groups would compromise the independence of the test set by using training subject data to calculate test subjects’ ISC values (Scheinost et al., 2019). We note that ISC estimates have been shown to stabilize at samples of around ~30 subjects (Pajula and Tohka, 2016), so the magnitudes of the ISC values calculated using either the average of data from 87 or 175 other subjects would not be expected to differ significantly. Next, gISC is calculated for each individual by taking the average of their Fisher z-transformed parcel-level ISC values.

Identify FC edges whose connectivity strength is most associated with gISC in the 88 subject training group. Spearman partial correlations are performed between each of the 64,620 unique FC edges and gISC, controlling for tSNR and mean FD in both the concatenated rest and movie day 1 scans. Spearman rank correlations are used here and throughout the rest of the paper because the distribution of gISC values exhibited a significant left skew (as did ISC values in general; gISC skewness Day 1 = −1.2, Day 2 = −1.7). All edges whose strengths were correlated with gISC at P < .01 were retained, although this feature selection step is not strictly necessary for rCPM and exists largely to reduce computational demands (Gao et al., 2019).

Fit a regularized linear model using the features (RSFC edges and mean FD measures) identified in the previous step as predictors and gISC as the response. Hyperparameters for this model include alpha (the ridge coefficient), whose optimal value has been shown to be near-zero for true ridge regression, and lambda, which is calculated in an inner fold using the method and code of (Gao et al., 2019).

Predict gISC for each participant in the test set by multiplying that participant’s RSFC edge strengths for the edges that passed the feature selection step by each edge’s ridge coefficient and adding the intercept. The ridge coefficients, intercept, and optimal lambda value for the present linear model are saved to allow for future bootstrap aggregation (described in the next section).

Flip the training and testing groups and repeat steps 1–4.

Evaluate model performance by calculating the Spearman partial correlation between the predicted and observed ISC values across participants in the test set, controlling for mean FD and tSNR from the rest and movie scans.

Significance testing the cross-validated models:

A permutation scheme was used to assess the statistical significance of the prediction coefficients generated in the 100 split-half iterations. First, the order of the gISC and movie watching tSNR and mean FD matrices was shuffled such that one participant’s resting state FC, tSNR, and mean FD values were associated with the movie BOLD time courses, tSNR values, and mean FD values from a random participant. Steps 1 through 6 from the previous section were then performed and the whole process was repeated 10,000 times. The resulting 10,000 correlation coefficients serve as a null distribution with which the following permutation p-value was calculated (Finn and Bandettini, 2021): p = sum(rnull > x̃) + 1/10,001, where x̃ is the median of the 100 correlation coefficients obtained through the true (unshuffled) models.

Building and testing the bagged model:

The cross-validation paradigm described in the previous two sections was used to ask whether a model trained on data from one group of subjects can predict gISC from RSFC in a novel group of subjects. Importantly, the gISC values used here all come from the same set of stimuli (the Day 1 movie clips). It could then be the case that the CV models rely on stimulus-specific signals and may fail to predict gISC during viewing of a different movie. A bootstrap aggregating, or “bagging,” approach was used to test whether the 200 linear models trained on Day 1 movie watching and resting state data could predict Day 2 gISC (derived from a different set of stimuli) from Day 2 RSFC, as previous work has shown bagged CPM models to be more accurate and more generalizable than their non-bagged counterparts (O’Connor et al., 2021). To construct the bagged model, RSFC edges that passed the P < .01 feature selection step in at least 10% (20/200, reflecting the 100 iterations of split-half cross-validation) of iterations were identified, yielding 1437 edges total. Next, for each of these edges, the ridge coefficients for a given edge across all iterations in which that edge was selected were averaged. The average intercept across the 200 bootstraps was then combined with the average ridge coefficients to yield a singular composite linear model. Each subject’s Day 2 RSFC matrix was then submitted to this model (by multiplying the correlation coefficients of selected RSFC edges by the ridge coefficients and adding the intercept) to generate 176 predicted gISC values, which were then Spearman correlated with the “true” gISC values (calculated through leave-one-out ISC using the full sample), controlling for mean FD and tSNR in the Day 2 rest and movie scans. To characterize network-level contributions of edges to the bagged model, the 1437 ridge coefficients were averaged according to their corresponding edge’s network memberships separately for positive and negative weights.

To further evaluate possible relationships between nuisance variables and gISC prediction, we divided our subjects into low- and high-motion groups (n = 88 each) using a median split of mean FD over all four scans (median mean FD = 0.038) and repeated the above rCPM analysis independently in each group (still controlling for mean FD and tSNR). Finally, global intrasubject correlation (gIntraSC), measured by correlating one subject’s BOLD time courses from two viewings of the same stimulus and averaging correlation coefficients across parcels, could also be considered as a measure of data quality. Taking advantage of the 83 s test-retest clip, we calculated gIntraSC values for each subject and each day of data collection and repeated the above analysis with this new covariate. Code for these rCPM analyses was adapted from https://github.com/YaleMRRC/CPM (Greene et al., 2020).

Inbound and outbound analyses

To investigate relationships between RSFC and ISC in each of the 360 Glasser parcels, we used two related partial correlation analyses. First, we (Spearman) correlated ISC values for a given parcel with the RSFC coefficients for each edge involving that parcel across participants, controlling for mean FD and tSNR in the rest and movie scans. This “inbound” analysis describes how ISC in a given parcel is related to that parcel’s resting state functional connections to every other parcel. Next, the “outbound” analysis involves finding the Spearman partial correlation between the RSFC coefficients for each connection between a given parcel and every other parcel and ISC in those other parcels, again controlling for rest and movie scan mean FD and tSNR. This analysis shows how a given parcel’s resting state functional connections to other parcels are associated with ISC in those parcels. Both the inbound and outbound analyses yield a 360×1 map (with one undefined value for the reference parcel) across the brain for every parcel and for each day of data collection, allowing for the visualization of results on the cortical surface.

To identify parcels in which ISC was most positively/negatively associated with RSFC across the brain, the inbound and outbound maps were Fisher z-transformed, averaged, and converted back to units of Spearman’s ρ. Next, these inbound/outbound averages were grouped by resting state network (RSN) membership according to the Cole-Anticevic Brain-wide Network Partition (Ji et al., 2019) to facilitate visualization of network-level trends. Finally, mean outbound values from both days were Spearman-correlated with rankings along a sensorimotor-association axis, where low ranks correspond to low-level sensory regions (e.g., V1, A1) and high ranks correspond to high-level association areas (e.g., dmPFC, PCC; developed by Sydnor et al. (2021) and downloaded from https://pennlinc.github.io/S-A_ArchetypalAxis/). Only left hemisphere parcels were used for this analysis, as the sensorimotor-hierarchy ranking created by Sydnor et al., did not include the Glasser parcellation’s right hemisphere.

Because the number of parcels belonging to each Cole-Anticevic RSN is an arbitrary value based on the resolution of the Glasser parcellation, we used bootstrap confidence intervals (CIs) to test whether the averages of the mean inbound and outbound maps within RSN were significantly different from zero. Non-parametric bootstrap CIs were constructed by randomly sampling 176 subjects with replacement, repeating the inbound and outbound analyses, and recalculating the RSN averages 1000 times. To account for multiple comparisons across the twelve RSNs, 95% CIs were Bonferroni corrected using the formula 1 – (alpha/number of RSNs), leading to an effective CI of 99.6% for each RSN average on each day.

The inbound analysis introduced earlier in this section addresses how connectivity from parcel A to parcels B-Z is associated with ISC in parcel A. Although this analysis neatly parallels the dimensionality of the outbound analysis, it ignores how resting state functional connections that don’t involve parcel A (i.e., connectivity between parcels B-Z) might be associated with ISC in that parcel. To address this limitation, we used a “full” inbound analysis that correlates all edges of the RSFC matrix with ISC in the reference parcel.

Spatial permutation tests

Similarity between the Day 1 and Day 2 inbound/outbound maps was quantified using Spearman correlation. Because (1) the number of samples in these correlations is determined by parcellation resolution and (2) spatial autocorrelation is present among adjacent parcels, a spatial permutation test was used to assess the significance of these correlations (and all other spatial correlations performed here) following the method and code of Váša et al. (2018). Briefly, this test works by inflating one of the spatial maps from each correlation to a sphere, randomly rotating this projection, and calculating the Spearman correlation between the empirical map and the randomly rotated projection. This procedure is repeated until a desired number (here, 10,000) of permutations have been performed, and the final permutation P value reflects the number of null permutations for which the resulting Spearman correlation is greater than the observed correlation divided by the total number of permutations. Although recent work has identified that this approach relies upon often unrealistic statistical assumptions (Weinstein et al., 2021), we note that no alternative method is compatible with our effect maps which are undefined at the individual subject level. As such, the present spin method can be seen as the best available significance test for the day-to-day consistency of our observed effects.

Principal component analysis

Principal component analysis (PCA) was performed separately for outbound and (full) inbound matrices from both days of data collection to simplify the matrices into a small number of orthogonal factors. Only the first two PCs of the outbound PCA are visualized in Fig. 4 as these components were sufficient to account for > 50% of the variance in the outbound maps on both days.

RSFC and ISC topological similarity

One pairwise RSFC similarity matrix was generated for each day of data collection by (Spearman) correlating RSFC profiles (i.e., each subject’s set of 64,620 unique FC edge strengths) from all possible pairs of subjects, leaving two 176 subjects x 176 subjects matrices in which higher values correspond to greater RSFC profile similarity for those two individuals. The same procedure was repeated using each subject’s ISC topologies, again yielding two 176 subjects x 176 subjects matrices. Here, we noted that four subjects, two different subjects on each day, had markedly different ISC topologies from the rest of the sample. These individuals also had the lowest gISC values for their respective movie scans, suggesting minimal engagement with the stimulus. Data from these outlier scans were discarded from the rest of the analyses described in this section for a sample size of 174 subjects, although we note that including these subjects does not change the interpretation of our results. These subjects were not excluded from the previous analyses as they were not identified as outliers in the Figs. 2D and 3 scatter plots.

Pairwise motion, tSNR, and gISC similarity matrices were also generated by calculating the negative absolute value of the difference between every pair of subjects’ mean framewise displacement for each day’s rest and movie scans as well as their (Fisher z-transformed) gISC values. Finally, a Spearman partial correlation was performed between each day’s RSFC and ISC similarity matrices, controlling for pairwise similarity in motion, tSNR, and gISC. A permutation test was used to evaluate the significance of the RSFC-ISC topological similarity relationships. Specifically, the subject-level ISC, gISC, and movie watching FD and tSNR matrices were shuffled using the same randomly generated order across subjects and pairwise ISC/gISC/tSNR/FD similarities were recalculated. The Spearman partial correlation described above was then performed using the shuffled movie vectors and the original resting state vectors such that each pair’s mean movie watching FD and tSNR, gISC, and ISC pattern similarity values were associated with the mean resting state FD and tSNR and RSFC pattern similarity values from a random pair. This process was repeated 10,000 times to generate the null distributions shown in Fig. 6B.

Replication with Schaefer 400 subject-specific parcellations

To probe the sensitivity of our results to parcellation choice and individual differences in parcel boundaries, we repeated several of our analyses after re-parcellating our data with individual-specific versions of the Schaefer 400 parcellation generated by Kong et al. (2021) (downloaded from https://github.com/ThomasYeoLab/CBIG/tree/master/stable_projects/brain_parcellation/Kong2). Parcellations were not available for eight of our 176 subjects, leaving us with a sample size of 168 for these analyses (and 167/166 for the topological similarity analysis due to the excluded subjects).

Results

Functional connectivity at rest predicts a global measure of normative movie-evoked activity

Ridge connectome-based predictive modeling (rCPM) was used to predict each participant’s average ISC across all parcels (Fig. 1). This cortex-wide measure, which we refer to as global ISC (gISC), serves as a low-dimensional marker of brain synchronization, or how similar an individual’s temporal patterns of neural activity during movie watching are to the rest of the group’s. We used rCPM to ask two related questions: to what extent can a model trained on data from one set of subjects watching one set of movie clips predict gISC in (1) a different set of subjects watching the same clips, and (2) the same subjects watching a different set of clips.

Fig. 1.

Ridge CPM pipeline. (A) Predicting gISC from RSFC in held-out subjects: RSFC matrices were created by (Pearson) correlating activity time courses from all pairs of parcels for each individual and day of data collection. Subjects were divided into two commensurate groups, and ISC values for every parcel were calculated for each subject as the Pearson correlations between their activity time course in a given parcel and that of the group average. These 360 parcel-wise values were averaged to yield 1 gISC value per subject per day of data collection. A model was then trained to predict gISC values from RSFC data in one group and tested on data from the held-out group, with model performance being evaluated as the Spearman partial correlation between predicted and actual Day 1 gISC values, controlling for head motion and tSNR. This procedure was repeated 99 more times using different train/test group splits, yielding 200 linear models and 100 correlation coefficients reflecting prediction accuracies. (B) Predicting gISC from RSFC in held-out stimuli: The 200 linear models were aggregated into one composite model which was tested on Day 2 gISC and RSFC data from all 176 subjects. Once again, model accuracy was evaluated using the Spearman partial correlation between predicted and actual Day 2 gISC values, controlling for Day 2 head motion and tSNR.

Using only the Day 1 rest and movie data, we performed 100 iterations of rCPM with split-half cross-validation (CV) to predict gISC from RSFC in unseen subjects. Parcel-level ISC (averaged across individuals; Fig. 2A) and individual-level gISC (averaged across parcels; Fig. 2B) were consistent across the two days, as reflected by Spearman correlations of ρ = 0.97 (Pspin < 1 × 10−4 ) and ρ = 0.68 (P = 2.2 × 10−24 ) respectively. We note that the ρ = 0.68 test-retest reliability value serves as a sort of upper bound on prediction performance as it represents how well gISC can predict itself across the different stimuli. Because the reliability of gISC was calculated while controlling for tSNR and head motion (quantified by mean framewise displacement; mean FD) from both movie scans, rank residuals were used to create the scatterplots shown in Fig. 2 and elsewhere in the paper when appropriate (although plots showing raw values can be found in Supplementary Fig. S1).

Fig. 2.

RSFC predicts gISC during movie watching in novel stimuli and in novel subjects. (A) ISC analysis reveals a non-uniform cortical distribution of shared stimulus-evoked activity that is consistent across Day 1 and Day 2 stimuli. (B) Individual differences in gISC were consistent across the Day 1 and Day 2 scans. In this and subsequent scatter plots, each dot represents an individual subject unless noted otherwise. As these plots represent Spearman partial correlations (controlling for head motion and tSNR), the x- and y-axes are in units of rank residuals so that the slope of the best fit line represents the corresponding Spearman correlation coefficient. (C) A bagged rCPM model trained on Day 1 data predicts Day 2 gISC from RSFC with considerable accuracy. (D-E) Heatmaps illustrating that the bagged rCPM model’s negative and positive features were broadly distributed across functional networks. (F) Models trained on Day 1 data from one set of subjects predict Day 1 gISC in held-out subjects at above-chance accuracy. Red dots represent prediction accuracies from the actual 100 CV iterations, whereas gray dots reflect null prediction accuracies from re-running the rCPM analysis using permuted gISC and movie FD values 10,000 times.

The bagged model was able to predict Day 2 gISC from Day 2 RSFC with significant accuracy (ρ = 0.51, P = 1.4 × 10−12; Fig. 2C) indicating that an individual’s ability to exhibit normative stimulus-evoked activity is closely related to their brain’s intrinsic functional architecture. To visualize the functional distribution of RSFC edges that contributed to the model’s predictive performance, we averaged positive and negative edges separately according to resting state network (RSN) membership as defined by the Cole-Anticevic Brain-wide Network Partition (Ji et al., 2019). Predictive edges were widely distributed across RSNs according to no obvious pattern, although positive relationships between RSFC and gISC were more pronounced than negative relationships (as seen in the scales of Figs. 2D–E).

Importantly, individual differences in RSFC are known to depend on factors that could complicate the interpretation of our results, including head motion, parcellation selection (Bryce et al., 2021), and individual differences in areal topologies (Gratton et al., 2018; Gordon et al., 2017). To control for head motion more conservatively, we assigned subjects to a low or high motion group based on a median split over mean FD and repeated this CPM analysis in both groups independently, finding that Day 2 gISC could be predicted in both groups (low motion ρ = 0.55, P = 8 × 10−8; high motion ρ = 0.53, P = 1.7 × 10−7 ; Supplementary Fig. S2). RSFC could also be used to predict Day 2 gISC when controlling for individual differences in global intrasubject correlation, a potential measure of data quality (Supplementary Fig. S3). Finally, we repeated the rCPM analysis using individualized versions of the Schaefer 400 parcellation, again finding that Day 2 gISC could be predicted from RSFC (Supplementary Fig. S4).

To examine how models trained on data from one group of subjects can predict gISC in held-out individuals watching the same stimulus, we plotted the prediction accuracies obtained from the 100 CV iterations (Fig. 2F, red dots) against a distribution of null accuracies (Fig. 2F, gray dots) generated by permuting gISC values across subjects and re-running the rCPM pipeline 10,000 times. As shown in Fig. 2F, the median CV accuracy was significantly greater than would be expected by chance (median ρ = 0.15, permutation P = .037).

Inbound and outbound analyses characterize region-specific RSFC-ISC relationships

In the previous section, we demonstrated that an individual’s RSFC profile is predictive of their cortex-wide ISC during movie watching. Next, to investigate region-specific RSFC-ISC relationships, we introduce the related “inbound” and “outbound” analyses. The inbound analysis, demonstrated in Fig. 3A using Day 1 data, describes how ISC in a given parcel is associated with that parcel’s resting connectivity to every other parcel (raw values are shown in Supplementary Fig. S5). The intuition behind and results of this analysis are first exemplified using the left frontal eye field (FEF), a parcel chosen for its circumscribed role in the well-studied saccade pathway (Felleman and Van Essen, 1991; Paus, 1996). To create the brain maps shown in the top row of Fig. 3A, left FEF ISC was correlated with RSFC between left FEF and the other 359 cortical parcels while controlling for head motion and tSNR from the corresponding rest and movie scans. This analysis revealed that individuals whose left FEF was more connected to visual and parietal areas at rest exhibited more typical FEF activity during movie watching, consistent with the FEF’s known inputs from visual and parietal cortex (Felleman and Van Essen, 1991; Paus, 1996).

Having demonstrated that the inbound analysis reveals expected RSFC-ISC relationships in a relatively unimodal and well-studied area, we next sought to apply the same analysis to a parcel in the temporo-parieto-occipital junction (TPOJ2), a multimodal area in the temporoparietal junction/posterior superior temporal sulcus region whose functioning is less understood but has been implicated in sensory integration during narrative processing (Lerner et al., 2011; Patel et al., 2019). Reflecting TPOJ2’s multimodal function, greater ISC in this parcel was associated with increased connectivity to unimodal sensory (e.g., visual, auditory, somatomotor) cortices and decreased connectivity to higher order areas (e.g., angular gyrus, dorsolateral prefrontal cortex; dlPFC) at rest (Fig. 3A, lower row). Although only Day 1 results are visualized in Fig. 3, Spearman spatial correlations confirmed that the inbound maps for both left FEF and TPOJ2 were consistent across both days of data collection (spatial ρ = 0.90, P < 1 × 10−4 and spatial ρ = 0.98, P < 1 × 10−4 respectively).

Fig. 3.

Inbound and outbound relationships between Day 1 RSFC and ISC in FEF and TPOJ2. (A) Inbound analysis: This analysis describes how RSFC between the reference parcel (yellow) and some other parcel (red) is associated with ISC in the reference parcel. Greater resting connectivity from FEF (or TPOJ2) to red-shaded parcels is associated with greater FEF (or TPOJ2) ISC during movie watching. Example inbound maps are shown for parcels covering left FEF (top row) and TPOJ2 (bottom row). Corresponding scatter plots illustrate how the inbound correlation coefficients for the bordered parcels are calculated. (B) Outbound analysis: This analysis describes how RSFC between the reference parcel (yellow) and some other parcel (red) is associated with ISC in the other parcel. Greater RSFC from FEF (or TPOJ2) to red-shaded parcels is associated with greater ISC in those parcels during movie watching. Corresponding scatter plots illustrate how the outbound correlation coefficients for the bordered parcels were calculated.

While the inbound analysis describes how a brain region’s functioning (ISC) is related to its intrinsic connectivity (RSFC) to other regions, a distinct but similarly informative relationship is how a region’s RSFC to other areas is associated with ISC in those areas. To quantify this “outbound” relationship, we calculated the Spearman partial correlation between the RSFC edge strengths for all connections involving a given parcel and ISC in every other parcel. The FEF outbound map therefore illustrates that greater RSFC between FEF and dlPFC is associated with greater ISC in dlPFC (Fig. 3B), again replicating known hierarchical relationships between these two areas (Arkin et al., 2020; Corbetta et al., 2008). The TPOJ2 outbound map shows a diffuse pattern that is notably attenuated in somatosensory, auditory, and visual regions compared to the corresponding inbound map, illustrating that RSFC between sensory regions and TPOJ2 is more relevant to TPOJ2 ISC than to sensory cortex ISC. Both the FEF and TPOJ2 outbound maps were consistent across the different days of data collection (spatial ρ = 0.84, P < 1 × 10−4 and spatial ρ = 0.85, P < 1 × 10−4, respectively). In addition to our FEF and TPOJ2 maps, we have also included inbound and outbound maps for the unimodal sensory areas A1 and V1 in Supplementary Fig. S6 for interested readers.

Network-level inbound/outbound relationships and means

The human cortex is organized into large-scale networks that sub-serve different sensory and cognitive functions during movie-watching. To characterize how movie-evoked activity in these networks is associated with intra- and inter-network RSFC, we averaged inbound and outbound values according to network membership to create the heatmaps shown in Fig. 4A. Greater intra-network RSFC was consistently associated with greater ISC, as indicated by the red diagonals of both heatmaps. On the other hand, inter-network RSFC showed more variable relationships with ISC, but these relationships tended to be consistently positive or negative within networks and inbound/outbound analyses. For example, sensory networks like visual 2, somatosensory, and auditory networks tended to exhibit more positive outbound relationships with other networks. Meanwhile, the posterior and ventral multimodal networks (PMN/VMN) exhibited positive inbound relationships with all networks except the frontoparietal and default mode networks (FPN/DMN), indicating that these areas encode more stimulus-related information when they are more broadly connected to sensory areas at rest.

Fig. 4.

Network-level trends in mean inbound and outbound maps. (A) Heatmaps illustrate inbound (shown in columns) and outbound (shown in rows) relationships averaged within functional networks. (B) Each parcel’s values on these surfaces reflect the average of the 359 values comprising its inbound map (e.g., the surfaces in Fig. 3A) on each day. (C) Same as B for the outbound maps. (D) Boxplots group the parcels from the surfaces in (B) according to their network affiliation, such that the y-coordinates in these plots are the same as the values on the above surfaces. The left and right boxes for each network reflect Day 1 and Day 2 results, respectively. The posterior and ventral multimodal networks had mean inbound values that were significantly greater than zero on both days. Bonferroni-corrected confidence intervals for all other networks contained zero on one or both days. (E) Same as D for the outbound analysis; only visual 2 and frontoparietal networks had average outbound values that were significantly different from zero on both days. (F) Relationships between mean outbound values and sensorimotor-association hierarchy ranks for all left hemisphere parcels (Day 1 on left, Day 2 on right). Because Spearman correlations were used for this analysis, the least-squares lines visualized here do not directly correspond to the accompanying correlation coefficients.

Classical theories of sensory information processing place cortical areas on a spectrum ranging from specialized to integrative (Mesulam, 1998). To identify integrative hubs whose normal functioning during movie watching is most associated with greater intrinsic connectivity to the rest of the brain, we visualized the average of each parcel’s inbound map and grouped these average values by RSN affiliation to illustrate network-level trends (Fig. 4B). At the network level, the PMN and VMN displayed the highest average inbound values (Fig. 4D). These were also the only RSNs whose average values were found to be significantly different from zero on both Day 1 (PMN mean ρ = 0.048, 99.6% bootstrap CI [0.020, 0.077]; VMN mean ρ = 0.045, 99.6% bootstrap CI [0.012, 0.078]) and Day 2 (PMN mean ρ = 0.050, 99.6% bootstrap CI [0.041, 0.084]; VMN mean ρ = 0.035, 99.6% bootstrap CI [0.0056, 0.069]).

We next averaged the outbound maps to identify source-like parcels to which greater connectivity at rest was associated with greater ISC during movie watching (Figs. 4C and E). The visual 2 and auditory networks exhibited the highest and second-highest average outbound values across Days 1 (mean visual 2 ρ = 0.044, 99.6% bootstrap CI [0.017, 0.068]; mean auditory ρ = 0.027, 99.6% bootstrap CI [−0.001, 0.056]) and 2 (mean visual 2 ρ = 0.047, 99.6% bootstrap CI [0.023, 0.060], mean auditory ρ = 0.025, 99.6% bootstrap CI [0.024, 0.076]). On the other hand, the FPN and DMN exhibited negative average outbound relationships on Day 1 (FPN mean ρ = −0.039, 99.6% bootstrap CI [−0.054, −0.022]; DMN mean ρ = −0.030, 99.6% bootstrap CI [−0.050, −0.0096]) and Day 2 (FPN mean ρ = −0.024, 99.6% bootstrap CI [−0.054, −0.023]; DMN mean ρ = −0.017, 99.6% bootstrap CI [−0.021, 0.021]). Because averaging positive and negative values together can potentially obscure interesting effects, we have included mean inbound/outbound values using of positive- and negative-only values in Supplementary Fig. S7.

Although the Day 1 auditory and Day 2 DMN CIs contained zero, this overall pattern of results indicates that for the average parcel, greater connectivity at rest to unimodal sensory areas (i.e., auditory and visual networks) was most associated with greater ISC during movie watching, while the opposite was true regarding greater connectivity to higher-order networks (i.e., FPN and DMN). This suggests that a parcel’s mean outbound value reflects its position along a sensory-to-cognitive gradient. To test this hypothesis explicitly, we correlated outbound scores for each parcel with sensorimotor-association hierarchy rankings from (Sydnor et al., 2021). We found strong negative relationships between these two variables on both days (Fig. 4F; Day 1 spatial ρ = −0.71, P < 1 × 10−4; Day 2 spatial ρ = −0.70, P < 1 × 10−4 ), confirming that greater connectivity to more sensorimotor parcels at rest is associated with greater ISC during movie-watching (Fig. 4F). Finally, to determine the sensitivity of the mean inbound/outbound topologies identified here to parcellation choice, we generated the mean inbound/outbound maps using the Schaefer 400 individualized parcellations and found qualitatively similar maps (Supplementary Fig. S8).

Principal component analysis reveals low-dimensional embedding of inbound and outbound maps

Although the network-level analysis shown in Fig. 4A serves as a convenient way to reduce the dimensionality of the inbound and outbound maps, it is unable to reveal connectivity-activity relationships that cut across network boundaries. To explore the underlying structure of the 360 inbound and outbound maps without making a priori assumptions about network organization, we used principal component analysis (PCA), a dimension reduction technique that preserves more of the information in the original feature space than the averaging shown in Section 3.3. PCA of the outbound maps revealed that two PCs (Fig. 5A) were sufficient to account for a majority of the outbound map variance on both days, describing ~30% and ~20% of variance, respectively (Fig. 5B). Here, parcels with higher (redder) scores in the upper rows have outbound maps that look like the PCs shown in the lower rows. Relatedly, given the complementary nature of the outbound and inbound matrices, parcels with stronger loadings on these PCs have inbound maps that look more like the cortical maps shown in the upper rows.

Fig. 5.

Principal component analysis of the outbound maps. (A) Cortical surfaces display the scores and loadings for the first two principal components (PCs) of the outbound maps. Parcels in the lower row of surfaces that appear redder have inbound maps that look more like the surfaces in the upper row. Due to the symmetry of the inbound/outbound maps, redder parcels in the upper row of surfaces have outbound maps that look more like the surfaces in the lower row. (B) Scree plots show that the majority of outbound map variance can be accounted for by two PCs, with the percentage of variance explained tapering off after the third PC on both days.

The first PC was most positively expressed in visual, cingulo-opercular, and dorsal attention parcels, and was negatively expressed in DMN parcels. This PC loaded most strongly onto left TPOJ, but DMN parcels loaded either negatively or only weakly positively on this component. In other words, DMN parcels exhibited greater ISC when they were more connected to other DMN (and less connected to visual) parcels in this PC. The second PC, most positively represented in the superior temporal sulci (STS) and somatomotor cortices and negatively represented in FPN parcels, loaded most heavily onto left dlPFC. Although TPOJ and dlPFC had among the highest average inbound values in Fig. 4, the former’s loading onto PC1 and the latter’s loading onto PC2 suggests that greater ISC in these two integrative areas is most associated with connectivity to separate sets of regions. This difference is especially apparent when considering connectivity to STS, which was positively associated with dlPFC ISC (in PC2) but negatively associated with ISC in TPOJ (in PC1). Similarly, although DMN and FPN parcels both had similarly low outbound values (shown in Fig. 4A/C), their differential expression in PC 1 vs. PC2 suggests heterogeneous connectivity-activity relationships between these sets of parcels and the rest of the brain.

While the inbound analysis introduced in Section 3.2 neatly parallels the dimensionality of the outbound analysis, it only captures how ISC in parcel A is associated with its RSFC to parcels B-Z, ignoring RSFC between parcels B-Z. To fill this gap, results from a version of the inbound PCA that relates RSFC across all edges with ISC are visualized in Supplementary Fig. S9, which shows that ISC in different functional systems is most associated with notably different full RSFC patterns.

Subjects with more similar RSFC fingerprints share more similar spatial patterns of ISC

Our analyses thus far have focused on relationships between the magnitudes of RSFC and ISC across participants. However, the spatial pattern of ISC across a person’s brain likely reflects additional individual-specific aspects of brain function. An emphasis on topological patterns is apparent in functional connectome fingerprinting studies, which have shown that individuals can be identified by their unique set of (relative) RSFC edge strengths (Finn et al., 2015). However, the relevance of these whole-brain RSFC profiles to spatial patterns of ISC is still unclear.

To investigate whether RSFC profiles are related to how ISC, independent of its overall magnitude, is distributed across the brain, we used a pairwise similarity analysis to determine if individuals with more similar functional connectomes also shared more similar spatial patterns of ISC. First, we correlated RSFC and ISC profiles separately across participants to generate the similarity matrices shown in Fig. 6A. As expected, patterns of RSFC and ISC were largely conserved across participants (Day 1/Day 2 RSFC mean pairwise similarity ρ = 0.52/0.53, std. = 0.015/0.012; ISC mean pairwise similarity ρ = 0.74/0.74, std. = 0.017/0.021). On each day, the ISC topologies for two individuals were markedly different compared to the rest of the sample (seen as dark stripes in the ISC similarity heatmaps in 6A, different individuals on each day), so these four outlier movie scans were discarded from the following analysis. Next, we correlated the pairwise ISC and RSFC similarity values, controlling for similarity in gISC, and rest and movie scan head motion and tSNR. Across both days, we found that participant pairs with more similar RSFC profiles also shared more similar spatial patterns of ISC (Day 1 ρ = 0.16, permutation P = 1.0 × 10−3; Day 2 ρ = 0.18, permutation P = 5.0 × 10−4; Fig. 6B, upper row; raw values shown in Supplementary Fig. S10). Null distributions generated to test this relationship are shown in Fig. 6B (lower row) and illustrate that the topological similarity relationships are significantly stronger than would be expected by chance. Because individual differences in areal topologies could complicate the interpretation of this result, we repeated this analysis using the individualized Schaefer 400 parcellations, again finding very similar effects (Supplementary Fig. S11).

Fig. 6.

Individuals with more similar resting state connectomes have more similar spatial distributions of ISC. (A) Heatmaps show that RSFC and ISC profiles are generally conserved across individuals. (B) Scatterplots show the positive relationships between RSFC and ISC profile similarity, where each dot reflects one unique pair of subjects. To preserve visual clarity, the x- and y-axes reflect unranked ISC/RSFC similarity values after partialling out the effects of head motion, tSNR, and gISC similarity. The null distributions in the lower row visualize the results of matching one pair’s RSFC similarity with the ISC similarity of a random pair and repeating the Spearman partial correlation 10,000 times. Plotting the strength of the true RSFC/ISC profile similarity correlation coefficients on these distributions reveals that the observed effects are much stronger than would be expected by chance. We note that because unranked values were graphed on the scatterplots, the slopes of the best fit lines only indirectly correspond to the vertical lines on the null distributions.

Discussion

Here, we extend research into the neural processing of naturalistic stimuli by showing that an individual’s brain connectivity at rest is associated with their degree and distribution of normative movie-evoked brain activity. First, we demonstrated that distributed networks of RSFC edges can be used to predict global ISC across subjects and scanning sessions. Next, we explored regional RSFC-ISC relationships, finding that greater ISC in multimodal parcels is related to higher RSFC to unimodal sensory cortex and lower RSFC to frontoparietal and default mode networks. Finally, we showed that individuals with more similar RSFC fingerprints share more similar spatial patterns of ISC. Taken together, these results reveal new relationships between intrinsic brain connectivity and task-driven function while opening up additional lines of inquiry for the study of idiosyncratic responses to naturalistic stimuli.

Individual differences in ISC magnitude have previously been associated with behavioral traits such as top-down attentional control (Campbell et al., 2015). Building upon this literature, we used rCPM to identify a network of RSFC edges (i.e., the 1437 RSFC edges included in the bagged model) whose combined strengths accounted for a sizable proportion of between-subject variance in a global measure of ISC. Extending the work of Campbell and colleagues, this RSFC network could reflect a general signature of attentional ability, similar to networks identified by (Rosenberg et al., 2016; Rosenberg et al., 2020). In this case, individuals who express this RSFC profile more so than others would be expected to attend more to movie clips, leading to increased levels of ISC. At the parcel level, we would expect ISC in different parcels to be associated with the same 1437 RSFC edges that comprise the predictive network, as greater attention to the movie clips should lead to distributed increases in ISC (Ki et al., 2016; Regev et al., 2019; Song et al., 2021). Instead, in our full inbound analysis, we observed that relationships between parcel-level ISC and whole-brain RSFC patterns differed considerably across parcels. Specifically, 10+ principal components were needed to account for a majority of the full inbound matrix variance shown in Fig. S9. Given this result, the RSFC network might instead be more representative of the functional architecture of bottom-up information transfer, independent of top-down influences. Future studies that simultaneously evaluate relationships between RSFC, ISC, and behavioral traits like attentional control will provide valuable insight into the cognitive implications of the present findings.

In addition to identifying RSFC edges that predict global ISC, we also characterized region-specific RSFC-ISC relationships through our inbound and outbound analyses. First, our inbound analysis showed that the resting state connections most associated with ISC in a given parcel are not simply a function of its network affiliation. For example, although the FEF parcel studied here was shown to be more connected to cingulo-opercular areas than to visual cortex at rest (Ji et al., 2019), our inbound analysis showed that typical activity in this parcel during movie-watching was more strongly associated with its connectivity to visual and dorsal attention networks rather than cingulo-opercular cortex. Furthermore, although RSFC derived from bivariate correlation is inherently undirected, comparing the inbound and outbound maps for a given parcel allows for the generation of directional hypotheses. That visual-FEF connectivity is more associated with FEF ISC than visual cortex ISC suggests that, consistent with established models of the visual processing hierarchy (Felleman and Van Essen, 1991), stimulus-related information is more likely to flow from visual cortex to FEF. This expected finding, as well as the stability of these results across different rest and movie scans, motivates the use of inbound/outbound analyses to investigate how connectivity relates to function in more arcane regions and functional systems, as we demonstrated with the TPOJ.

After showing how the inbound and outbound analyses reveal stable RSFC-ISC relationships in specific regions, we next investigated network-level trends in these relationships. By averaging across the inbound maps, we found that parcels in the ventral and posterior multi-modal networks (located in temporoparietal cortex) exhibited the highest average inbound map values. In other words, activity in these parcels was more typical when they were more positively connected to the rest of the brain. This finding is consistent with previous work showing that these posterior multimodal regions consolidate audiovisual information (Baldassano et al., 2017; Ji et al., 2019; Lerner et al., 2011) and suggests that the integrative function of these regions is related to their distributed intrinsic connectivity. Complementing this result, we found that parcels in visual and auditory networks exhibited some of the greatest average outbound values, indicating that greater connectivity to these areas at rest was associated with greater ISC. This likely reflects the fact that auditory and visual areas are among the most synchronized during movie-watching and can be thought of as sources of stimulus-driven signal. Alternatively, connectivity to frontoparietal and default mode parcels, which tend to be among the least synchronized by movies, was on average inversely related to ISC. Subsequent PCA of the outbound maps clarified this result by showing that connectivity to default mode parcels was negatively associated with TPOJ ISC but positively associated with ISC in default mode parcels themselves. That these parcels exhibited more typical activity when they were more connected to functionally similar (and less connected to dissimilar) areas at rest is consistent with the idea that more modular cortical organization is associated with more efficient information flow (Wig, 2017). Lastly, the loading of left dlPFC and TPOJ onto different PCs indicates that ISC in both areas is closely related to their RSFC to different sets of parcels, suggesting that these areas may play important but distinct roles in the processing of naturalistic stimuli.

Beyond informing theories of cortical information processing, our results have implications for the understanding of circumscribed ISC deficits reported in patient populations. Our inbound/outbound findings provide evidence for the hypothesis that low ISC in a particular area may be driven by pathological alterations in inter-regional connectivity rather than dysfunction in the region itself. Future work outlining the causal direction of these RSFC-ISC relationships will inform the clinical applications of naturalistic fMRI (Eickhoff et al., 2020), as whether atypical activity in a region is self-limited or a consequence of its connectivity to other areas is an important consideration in developing targeted therapies (Fox et al., 2012).

Although naturalistic imaging studies have tended to focus on individual differences in the magnitude of ISC, it is important to note that this measure is susceptible to session effects that may obscure meaningful trait-related variance. For example, an attentive subject scanned at the end of a long day may exhibit less ISC than a more distractible individual scanned when they are most alert. The low test-retest reliability of activation magnitudes is a general issue in task fMRI research (Elliott et al., 2020). However, by asking “where” task activations occur instead of “how much” activation is present, researchers have been able to more reliably elicit individual differences and excise state-dependent noise from their analyses (Kragel et al., 2021). Moreover, individuals’ activation maps derived from modeled responses to difficult tasks have been accurately predicted from their RSFC profiles, suggesting that unique spatial activation patterns are intimately related to variation in the brain’s intrinsic wiring (Cole et al., 2016; Tavor et al., 2016). Extending this line of work to model-free naturalistic paradigms, we found that subjects with more similar resting state connectomes shared more similar spatial patterns of ISC during movie watching, independent of their similarity in ISC magnitude. This finding indicates that spatial distributions of ISC capture subject-specific information that ISC magnitudes alone do not, echoing recent reports that have established multivariate patterns as superior detectors of individual differences (Kragel et al., 2021).

Several limitations of this study constrain the interpretation and generalizability of our results. First, the short (3–5 min) clips used here may have failed to adequately engage cortical areas with longer temporal receptive windows (Baldassano et al., 2017). Future studies that utilize continuous stimuli may be better able to characterize RSFC-ISC relationships in these areas. Second, functional brain topologies have been shown to differ meaningfully across individuals and task states (Gratton et al., 2018; Salehi et al., 2020). Subjects identified as having lower RSFC and/or ISC in our analyses may instead have functional topologies that are less represented by the Glasser parcellation, complicating the interpretation of our inbound/outbound results. Our use of Spearman partial correlations to link RSFC and ISC presents its own set of limitations. Some RSFC-ISC relationships may be non-monotonic, such that both hypo- and hyper-connectivity between a pair of regions is associated with lower ISC. The correlation method used here also prevents us from assessing the unique contribution of single RSFC edges to parcel-level ISC, as would be possible with multiple regression. Finally, although we controlled for motion and tSNR while performing our correlations, other imaging artifacts and/or non-monotonic effects could still have contributed to our RSFC-ISC relationships.

Despite these limitations, this study constitutes an important step towards characterizing the complex relationships between normative stimulus-evoked activity and the brain’s intrinsic functional architecture. By linking resting state connectivity with task-driven function at the parcel, network, and cortex-wide levels, these findings further develop our understanding of sensory information processing and anchor individual differences in ISC in the resting state connectome. Additional research into these connectivity-activity relationships will be instrumental in exploring the cognitive and perceptual relevance of resting state fMRI and contextualizing new discoveries from the growing field of naturalistic imaging.

Citation diversity statement

Recent work in several fields of science has identified a bias in citation practices such that papers from women and other minority scholars are under-cited relative to the number of such papers in the field (Caplar et al., 2017; Dion et al., 2018; Dworkin et al., 2020; Maliniak et al., 2013; Mitchell et al., 2013). Here we sought to proactively consider choosing references that reflect the diversity of the field in thought, form of contribution, gender, race, ethnicity, and other factors. First, we obtained the predicted gender of the first and last author of each reference by using databases that store the probability of a first name being carried by a woman (Dworkin et al., 2020; Zhou et al., 2020). By this measure (and excluding self-citations to the first and last authors of our current paper), our references contain 5.56% woman(first)/woman(last), 3.7% man/woman, 27.78% woman/man, and 62.96% man/man. This method is limited in that a) names, pronouns, and social media profiles used to construct the databases may not, in every case, be indicative of gender identity and b) it cannot account for intersex, non-binary, or transgender people. Second, we obtained predicted racial/ethnic categories of the first and last author of each reference by databases that store the probability of a first and last name being carried by an author of color (Ambekar et al., 2009; Sood and Laohaprapanon, 2018). By this measure (and excluding self-citations), our references contain 8.59% author of color (first)/author of color(last), 16.32% white author/author of color, 21.51% author of color/white author, and 53.58% white author/white author. This method is limited in that a) names and Florida Voter Data to make the predictions may not be indicative of racial/ethnic identity, and b) it cannot account for Indigenous and mixed-race authors, or those who may face differential biases due to the ambiguous racialization or ethnicization of their names. We look forward to future work that could help us to better understand how to support equitable practices in science.

Data Availability

The raw HCP data used for this project can be downloaded from ConnectomeDB (db.humanconnectome.org). Code for all analyses can be found in the following Github repository: https://github.com/davidgruskin/hcp_rsfc_isc. All analyses were performed in MATLAB (R2020a). All cortical surface visualizations were performed with Connectome Workbench (Marcus et al., 2011).