Abstract
Networks hubs represent points of convergence for the integration of information across many different nodes and systems. Although a great deal is known on the topology of hub regions in the human brain, little is known about their temporal dynamics. Here, we examine the static and dynamic centrality of hub regions when measured in the absence of a task (rest) or during the observation of natural or synthetic visual stimuli. We used Magnetoencephalography (MEG) in humans (both sexes) to measure static and transient regional and network-level interaction in α- and β-band limited power (BLP) in three conditions: visual fixation (rest), viewing of movie clips (natural vision), and time-scrambled versions of the same clips (scrambled vision). Compared with rest, we observed in both movie conditions a robust decrement of α-BLP connectivity. Moreover, both movie conditions caused a significant reorganization of connections in the α band, especially between networks. In contrast, β-BLP connectivity was remarkably similar between rest and natural vision. Not only the topology did not change, but the joint dynamics of hubs in a core network during natural vision was predicted by similar fluctuations in the resting state. We interpret these findings by suggesting that slow-varying fluctuations of integration occurring in higher-order regions in the β band may be a mechanism to anticipate and predict slow-varying temporal patterns of the visual environment.
SIGNIFICANCE STATEMENT A fundamental question in neuroscience concerns the function of spontaneous brain connectivity. Here, we tested the hypothesis that topology of intrinsic brain connectivity and its dynamics might predict those observed during natural vision. Using MEG, we tracked the static and time-varying brain functional connectivity when observers were either fixating or watching different movie clips. The spatial distribution of connections and the dynamics of centrality of a set of regions were similar during rest and movie in the β band, but not in the α band. These results support the hypothesis that the intrinsic β-rhythm integration occurs with a similar temporal structure during natural vision, possibly providing advanced information about incoming stimuli.
- brain rhythms
- centrality
- functional connectivity
- magnetoencephalography
- natural stimulation
- resting-state networks
Introduction
Brain hubs are a set of topologically central regions considered critical for connecting different functional brain networks and integrating information among separate functional domains (Hagmann et al., 2008; van den Heuvel and Sporns, 2013; de Pasquale et al., 2017). Damage to cortical hubs has been shown to produce multiple cognitive deficits (Buckner et al., 2009; Crossley et al., 2014; Warren et al., 2014). Although the topology of hub regions has been described at different spatial scales (van den Heuvel and Sporns, 2011; Harriger et al., 2012; de Reus and van den Heuvel, 2013; Nigam et al., 2016), much less is known about their dynamics in time and frequency. This is largely due to the low temporal resolution of functional magnetic resonance imaging (fMRI), the method that has been predominantly used to investigate brain networks.
Cortical hubs flexibly reconfigure in response to a wide range of cognitive tasks (Cole et al., 2013; Crossley et al., 2013), but their resting topography remains stable during task performance (Buckner et al., 2009; Jin et al., 2014). This supports the idea that behavior requires the modulation of an otherwise fairly stable intrinsic network structure (Cole et al., 2014; Krienen et al., 2014; Gratton et al., 2016; Schultz and Cole, 2016). One leading idea is that this intrinsic functional architecture, sculpted into its adult state by development, individual experience, and learning (Lewis et al., 2009; Harmelech and Malach, 2013; Harmelech et al., 2013), constitutes a spatial prior of task networks (Raichle, 2011). An alternative hypothesis is that intrinsic connectivity reflects “idling” cortical circuits that are suppressed to facilitate task-specific synchronization. For instance, recent fMRI work indicates that behavioral tasks cause robust reorganization of the spatial patterns of functional connections (Krienen et al., 2014; Spadone et al., 2015; Kim et al., 2017), including hubs (Gratton et al., 2016). In addition, the observation of movie scenes reduces MEG band-limited power (BLP) connectivity in the α band between regions of the same network while increasing β-, θ-, and γ-band correlation between regions of different networks (Betti et al., 2013).
Here, we examined these two hypotheses by measuring, at high temporal resolution with MEG in healthy observers, functional connectivity during visual fixation (rest) and audiovisual stimulation. Specifically, we addressed two questions. First, are the static overall network topology and transient centrality of cortical hubs modified during the observation of movie clips compared with rest? In a previous work (de Pasquale et al., 2016), we showed that different brain networks interact transiently at rest through a stable set of cortical regions and that these regions show joint fluctuations of centrality over time (dynamic core network). Here, we investigated whether the overall network topology and dynamics of this core network are maintained (prior) or modified (idling) when going from rest to movie observation.
Second, we compared two viewing conditions: one in which the sequence of scenes was played normally (natural vision) and another in which the scenes were temporally scrambled (scrambled vision). In the visual cortex, the statistics of intrinsic brain activity are more similar to natural than artificial or synthetic stimuli (Fiser et al., 2004; Berkes et al., 2011) and the slow-varying timescale of natural stimuli resembles that of neural fluctuations in associative cortex (Hasson et al., 2008; Honey et al., 2012; Boly et al., 2015). Accordingly, the prior hypothesis predicts that the topology and hub dynamics would be more similar between rest and natural vision than between rest and scrambled vision. In contrast, the idling hypothesis predicts significant departures in topology and dynamics from rest in both movie conditions. Finally, in light of previous findings that movie observation induces frequency-specific alterations of BLP connections both within and between networks (Betti et al., 2013), we measured BLP connectivity in both the α and β band, the most consistent MEG correlates of fMRI resting-state networks (RSNs) (de Pasquale et al., 2010; Brookes et al., 2011; Hipp et al., 2012).
Materials and Methods
Participants
Twelve healthy participants (mean age 24.7, range 21–31, five females; same sample as in Betti et al., 2013) without significant psychiatric or neurologic diseases and with normal or corrected-to-normal visual acuity participated in the study. All subjects were right-handed as judged by the Edinburgh Handedness Inventory (Oldfield, 1971). The study was approved by the Ethics Committee of G. d'Annunzio University, Chieti (Italy) in accordance with the standards of the 1964 Declaration of Helsinki. Before the experiment, participants gave their written informed consent.
Audiovisual stimuli
Audiovisual stimuli were taken from the Italian version of the film “The Good, the Bad, and the Ugly” directed by Sergio Leone. Three segments lasting ∼5 min each were obtained from 16 min 48 s to 31 min 48 s (movie). To generate scrambled movie sequences, each movie clip was cut in shorter fragments 250, 500, 750, or 1000 ms long. The cuts were randomly applied but made to maintain a similar frequency distribution of the four kinds of fragments. Finally, the fragments were temporally permuted in pseudorandom order. The end result was a video clip of the same duration (∼5 min) in which, the temporal structure of the movie narrative was disrupted (Fig. 1A). This has been shown to affect activity in higher-level visual regions that have longer temporal integration and much less activity in early visual cortex (Hasson et al., 2008). This procedure was repeated for each movie clip to generate three different scrambled movie clips. The audiovisual stimuli were presented via an LCD projector onto a translucent screen positioned in front of the subject (23 × 11.5 visual degree, 1024 × 576 pixels). The timing and frame rate of the audiovisual stimulation was controlled across subjects (on a subset of subjects, the difference between frame rates was smaller than 1.2 × 10−2 Hz). The screen was placed inside the MEG magnetically shielded room. A set of MEG-compatible headphones was used to listen to the audiovisual stimuli.
Schematic overview of visual stimuli and analysis strategy. A, Visual stimuli used in the study taken from the film “The Good, the Bad and the Ugly.” The same video clips were fragmented at four time scales (250, 500, 750, and 1000 ms) and temporally permuted to produce scrambled movies (scrambled). The order of movie and scrambled was counterbalanced across subjects, whereas the fixation condition was always at the beginning. B, For each condition, the source space BLPs in the α and β bands were estimated for each voxel of the brain and over time. C, Pearson's correlation coefficient was used to estimate both the static and transient coupling of 45 nodes. D, Correlation matrices were thresholded between [0:0.05:1], resulting in 21 thresholds and at every time for the dynamic graph analysis. E, For each node, the BC was estimated at each threshold and, after averaging the BC across nodes, the threshold that returns the highest BC was used for the next analyses (red bars). For the static analysis, a range of thresholds was selected for the comparison across conditions (black arrows). The overall procedure was repeated for each frequency band. The β band is represented for illustrative purposes.
Experimental design and statistical analysis
Subjects began the experiment with three 5 min runs in which they were instructed to rest while fixating a small red cross on a black background without any other task to perform (visual fixation). Then, each subject viewed 3 different consecutive movie clips and 3 consecutive scrambled versions, each 5 min long. The order of the movie and scrambled conditions was counterbalanced across subjects, whereas the visual fixation was always presented at the beginning of the experiment to avoid any contamination effect of the audiovisual stimulation over the resting-state functional connectivity measures. Therefore, 36 runs were acquired for each experimental condition (fixation, movie, scrambled). However, three runs (each from a different subject) were excluded because of excessive movement and were not considered for further processing. To ensure a balanced design, one run from the other conditions was also removed. For instance, if a subject moved during visual fixation, for all subsequent analyses, we removed for that subject also one run for natural vision and one for scrambled vision. Therefore, we retained a total of 99 runs (i.e., 33 for each experimental condition).
To ensure the same observation condition as for fixation and to reduce eye movements during the recording session, a small fixation cross was placed at the center of each movie/scrambled clip. The experiment lasted 2 h (including subject preparation and recordings).
Data acquisition
MEG data were recorded using the 153-channel MEG system (Della Penna et al., 2000) installed in a magnetically shielded room at the University of Chieti, Italy. The channels were dc SQUID integrated magnetometers arranged over a whole-head helmet surface. Horizontal and vertical electrooculogram signals, together with the electrocardiogram (ECG) signals, were simultaneously acquired with magnetic signals and were used during preprocessing for offline artifact rejection. The neuromagnetic and electric signals were band-pass filtered at 0.16–250 Hz and sampled at 1025 Hz. The head position relative to the MEG sensors was measured before and after each recording session through a fit of the magnetic field generated by five coils placed on the scalp. The coil positions were digitized by means of a 3D digitizer (3Space Fastrak; Polhemus), together with anatomical landmarks (left and right pre-auricular and nasion) defining a coordinate system. In a separate session, a 3D high-resolution T1-weighted image was recorded from each subject's head (3D MPRAGE pulse sequence, TR = 8.1 ms, TE = 3.7 ms, voxel size = 0.938 mm × 0.938 mm × 1 mm), through a 3 T MR Philips Achieva scanner installed at the University of Chieti, Italy. The MRI anatomical images were coregistered to the head position acquired during the MEG acquisition and then used to build the individual volume conductor.
Data analysis
Overview.
Before providing a detailed description of our methodological approach, for the sake of simplicity, we will summarize the basic analysis steps. After artifact rejection through independent component analysis (ICA), for every subject, we estimated the source-space BLP in each voxel of the brain (Fig. 1B) in the α and β band. Static and dynamic α-/β-BLP connectivity were estimated through the Pearson's correlation coefficient after geometric correction for spatial leakage effects (Fig. 1C) across 45 nodes arranged in seven RSNs (Table 1). This set of nodes is the same as described previously (de Pasquale et al., 2010, 2012, 2016; Betti et al., 2013) and was defined from prior fMRI works as set of functional regions that lead to an accurate RSN classification at the level of individual subjects (Hacker et al., 2013; Baldassarre et al., 2014). Prefrontal and temporoparietal regions are not well covered by this procedure. In fact, these nodes are more variable at the individual subject level and their RSN assignment is also not so consistent (Hacker et al., 2013). We acknowledge that this could be a limitation of the study; in fact, all of the connectivity/graph metrics and related distribution depend on the specific set of nodes adopted in the study. However, one of the main objectives of this study was to test whether centrality measures (topology, dynamics) are preserved (prior hypothesis) or modified (idling hypothesis) in a set of RSNs modulated by the observation of natural scenes, including the visual (VIS), auditory (AUD), language (LAN), and somatomotor (SMN) networks (Hasson et al., 2004; Betti et al., 2013), associated with the control of the visuospatial attention, such as the ventral and dorsal attention networks (VAN and DAN) (Corbetta et al., 2008), metabolically active at rest but suppressed during task performance (i.e., the default mode network, DMN; Raichle et al., 2001); and including a set of cortical hubs previously identified at rest (de Pasquale et al., 2012, 2016).
Node coordinates (MNI) and labels
To investigate whether the audiovisual stimulation preserves (prior hypothesis) or changes (idling hypothesis) the resting connectivity we analyzed the correlation matrices, i) at the level of interaction strength through the analysis of the mean correlation within/across RSNs and ii) at the level of the functional architecture, by associating a graph model, as defined as a set of pairwise relations (edges) between regions (nodes), to the correlation matrices. In turn, the functional architecture was studied through analysis of global topology by using network-based statistics (NBS) (Zalesky et al., 2010) and of brain hubs (i.e., central regions) and their dynamics. To identify brain hubs, we used the binary betweenness centrality (BC), a measure implemented in the Brain Connectivity Toolbox (Rubinov and Sporns, 2010), which requires thresholding of the correlation matrices (Fig. 1D). For each band, thresholds to estimate the static and dynamic BC were selected as those maximizing the average BC across nodes and time (Fig. 1E). Details of each analysis step are provided below.
MEG preprocessing and BLP estimation.
We applied an ICA approach to remove environmental and physiological (e.g., cardiac, ocular) artifacts from sensor space MEG signals and retain the brain independent components (ICs) (Mantini et al., 2011). The sensor maps of the brain IC were scaled to [−100, 100] and were then projected in the source space, represented by a 3D Cartesian grid (4 mm voxel side) by means of a weighted minimum-norm least-squares (WMNLS) estimator implemented in Curry 6.0 (Neuroscan). The noise level used by the WMNLS estimator was set to 6% of the maximum absolute value of the IC map. This value was experimentally selected from independent sets of data (de Pasquale et al., 2010, 2012; Mantini et al., 2011; Marzetti et al., 2013). The individual 3D grid was then projected onto the MNI 152 atlas space through SPM8 so that every voxel centroid was assigned to a set of MNI coordinates. The activity of each voxel in the grid was estimated as the linear combination of the brain IC time courses multiplied by the related weight vectors into the source space. Our analysis was restricted to α- and β-BLP connectivity, the main MEG correlates of fMRI RSNs (Mantini et al., 2007; de Pasquale et al., 2010, 2012; Brookes et al., 2011; Betti et al., 2013; Marzetti et al., 2013). For each voxel in the MNI grid (Fig. 1B), the source space signal was filtered in the α and β bands using separate high-pass and low-pass Chebyshev II IIR digital filters with order 10 (the order was 8 only for the high pass in the α band) and stopband ripple 10 dB. The source-space MEG BLP time series were estimated as pj(t) = 
, where q→j(t′) is the vector activity in voxel j, W = 150 ms is the sliding window size, and t is the time sampled at window shifts corresponding to 20 ms. Because our previous study (Betti et al., 2013) showed that movie watching induces the largest modulations of the interaction strength in these two bands and at frequencies <0.3 Hz, we restricted our analysis to this slow band.
Estimation of BLP connectivity.
The static and dynamic BLP interactions were estimated through the Pearson's correlation coefficient. To minimize spatial leakage effects on the estimation of BLP connectivity, we applied the geometric correction scheme (GCS) whereby leakage from a seed location is modeled on the basis of the forward and inverse models and eliminated before BLP correlation estimation (Wens et al., 2015; Wens, 2015). Specifically, given the source-projected brain IC maps Aic and the associated row uic of the unmixing matrix (1 ≤ ic ≤ Nic), a seed location r0 and its lead field weights L(r0) (the lead field being also projected onto the MNI 152 atlas brain), the GCS provided a set of corrected brain IC maps obtained as
where the superscript + denotes pseudoinversion. These maps were used for the reconstruction of the target sources at r ≠ r0 before computing their BLP correlation with the original, uncorrected seed time series. Equation 1 was derived as a particular case of the GCS (Wens et al., 2015) corresponding to the linear inverse operator W = Σic Aic uic associated with our source reconstruction pipeline (Mantini et al., 2011). Because GCS may be affected by local miscorrection effects mainly due to seed mislocalization (Wens et al., 2015), we masked out in subsequent analyses pairwise correlation values between nodes closer than 45 mm, which is compatible with simulation results reported previously (Wens et al., 2015). Accordingly, these values were removed from all the next analyses. For each run, the static connectivity across nodes was computed as the mean Pearson's correlation coefficient estimated over the filtered, leakage-corrected BLP time series in eight nonoverlapping windows lasting 40 s. By averaging the interaction matrices across runs, we obtained for each experimental condition a group-level static correlation matrix (Fig. 1C). The dynamic interaction matrices were instead obtained using windows lasting 10 s, sliding by 200 ms time. This duration was retrieved in our previous study (Betti et al., 2013) showing that naturalistic viewing influences internodal BLP interaction for fluctuations around 0.1 Hz (the reciprocal of 10 s). Importantly, the correlation time series during naturalistic and time-scrambled movie segments were aligned to the beginning of the stimuli, for each run.
Analysis of interaction strength.
For each subject, the correlation values were averaged, network by network, across its possible node pairs (within-network connectivity; e.g., within-VIS) or across its nodes and those belonging to other networks (between-network connectivity; e.g., VIS-DAN; VIS-AUD; VIS-DMN). The between-network coupling was then averaged across all network pairs. The modulation of the correlation induced by the audiovisual stimuli was quantified by applying the Fisher z transformation to these values and then using repeated-measures ANOVA with condition (fixation, movie, scrambled) and coupling type (within- and between-networks) as main factors separately for the α and β bands. The Duncan test was used for post hoc comparisons.
Analysis of network architecture.
To investigate putative changes of the intrinsic functional architecture induced by the audiovisual stimulation, the correlation matrices obtained for every subject and experimental condition were converted into graphs and then analyzed with the NBS toolbox separately for each band. NBS is a statistical nonparametric technique that operates directly on raw connectivity values and seeks to identify potential connected structures formed by a set of suprathreshold links (graph components) (Zalesky et al., 2010). For each comparison between fixation and movie/scrambled, changes in graphs components were tested by using a primary (t statistic) threshold of 3. Permutation testing (n = 5000) were then used to ascribe a p-value. Each component identified by NBS satisfied p < 0.05. If no significant changes were found using such threshold, the NBS statistics was repeated over a range of decreasing thresholds (from 3 to 2.2). For the graph visualization, we used the MATLAB toolbox BrainNet Viewer (Xia et al., 2013).
The node centrality at rest and during task was analyzed using the measure of BC implemented in the Brain Connectivity Toolbox (Rubinov and Sporns, 2010). BC is related to the number of times a node acts as a bridge between the shortest paths of any two nodes for all the possible pairs considered. The static binary BC was estimated from the group-level voxelwise functional connectivity matrices obtained for each experimental condition (see Fig. 1C, left) according to the following formula:
where σij is the total number of shortest paths from node i to node j and σij(v) is the fraction of those paths passing through v (Rubinov and Sporns, 2010). To compare the pattern of BC across conditions, this measure was normalized to the range [0,1] through the factor 1/[(N − 1)(N − 2)] where N is the number of nodes in the network.
To binarize the correlation matrices we selected a range of thresholds. This strategy avoids the selection of a single common threshold, which may bias the comparison across conditions (Fornito et al., 2013). To identify this range, we binarized the group-level connectivity matrices at each threshold between 0 and 1 with step 0.05 (Fig. 1D) and then estimated the binary BC for each node as a function of the threshold (see Fig. 1E, left). After averaging the BC across nodes separately for each experimental condition and band, we identified a set of thresholds leading to a BC higher than the mean + 3 SDs. The union of these threshold intervals defines a range of thresholds for each band (for the α band: [0.60–0.70]; for the β band: [0.55–0.65]; Fig. 1E, left, red bars). In Figure 1E, the black arrows mark the threshold interval obtained in the β band. Our results report the average BC in that range. To assess the significance of the BC in each node, we estimated a population of 100 random graphs preserving the node degree distribution (Maslov and Sneppen, 2002) for each selected threshold. These BC values were averaged across thresholds and nodes and constituted the significance level for BC (see Fig. 5A, colored solid lines). Those nodes with the highest BC and above the significance level in all the experimental conditions will be referred to as core network nodes (see Fig. 5B).
Dynamics of the core network.
To analyze the dynamics of centrality over time, we first converted into graphs a series of covariance matrices obtained for each subject, run, and time sample (Fig. 1C, right). Then, we dynamically estimated the binary BC every 200 ms run by run separately as follows:
where σij(t) is the total number of shortest paths from node i to node j at time t and σij(v, t) is the number of those paths passing through v at time t.
The dynamic analysis not only requires the BC estimation over time for each node and run, but it also needs to be replied for a population of 100 random graphs. To reduce the computation burden, here we selected a single binary threshold. As described previously (de Pasquale et al., 2016), we assumed the criterion of the highest average graph centrality to identify this threshold. Therefore, we first estimated the binary BC for every run and every node as a function of time over the 0–290 s interval (with step Δt = 200 ms) and threshold ranging between 0 and 1. BC values were then averaged across time and runs separately for every experimental condition (Fig. 1E, right). After averaging the BC across nodes, we selected the final thresholds that returned the highest BC. This threshold was the same across conditions: 0.85 and 0.80, respectively, for the α and β bands (Fig. 1E, right, red bars, for the β band).
To determine whether both viewing intact and time-scrambled movie sequences influenced the spectral content of BC with respect to rest, we estimated the power spectrum density (PSD) by means of the Welch's PSD estimator using 512 points (0.0098 Hz frequency bin; 0% overlap) and Hamming window filters. Then, the PSD was integrated over 5 frequency bands (0–0.1 Hz; 0.1–0.2 Hz; 0.2–0.3 Hz; 0.3–0.4 Hz; 0.4–0.5 Hz) and two-way repeated-measures ANOVA was run with band (five levels) and condition (three levels: fixation, movie, scrambled) as main factors. The ANOVA was carried out for every hub region comprised in the core network (based on Fig. 5B) and separately for the α and β bands.
Joint dynamics of BC.
Our previous study showed that, at rest, the centrality of nodes belonging to a core network fluctuates over time (de Pasquale et al., 2016) and they exhibit some periods of overlap (Fig. 2A). Here, we addressed whether task performance influences these joint dynamics of the BC. To this aim, we first binarized the centrality time courses as follows. For every hub in the core network, we estimated the mean across time of the BC of the corresponding random graph. For each experimental condition, the average over the time points in which the BC was above this value returns a threshold used to binarize the corresponding BC time course for every core hub region. These time courses were then concatenated over runs and subjects.
Dynamic fluctuations of BC during rest. A, Time course of BC in the β band during fixation in a representative subject taken from a subset of hub regions of the core network. B, Binarized time course of BC showing that there are epochs during which BC of hubs belonging to the core network are temporally overlapping; that is, epochs of BC higher (white) and lower (black) than the threshold often occur simultaneously for hub pairs. Red rectangles indicate example of epochs characterized by high (TPi) and low (TNj) joint BC for the hub pair right mPFC/right PUT. The rectangle in yellow indicates an epoch of absence of an hub and presence of another (FPk).
To assess whether hub regions in the core network showed true joint dynamics over time and not random joint fluctuations (e.g., as in the case of a hub with BC always above the threshold and a second hub with BC fluctuating across the threshold), we quantified whether a given hub was highly central in the same epochs of high centrality of another hub, by means of the logical operator true/false (TF) rate as follows:
where TPR represents the true-positive rate and FPR represents the false-positive rate. These quantities were computed as follows:
where TP = presence of Hub1 and Hub2; TN = absence of Hub1 and Hub2; FP = (presence of Hub1 and absence of Hub2 + absence of Hub1 and presence of Hub2)/2; and FN = FP.
Figure 2B depicts an example of an epoch of high joint BC (TPi), low joint BC (TNj), and an epoch of presence of a Hub and absence of another (FPk) for two representative hub regions.
The computation of TF was performed for every pair of hubs in the core network; therefore, for each experimental condition and band, we obtained a TF matrix in which values around 1 represent chance level. To avoid bias due to overcorrection effects, network nodes closer than the mask distance, such as right and left mPFC (DMN) were not considered simultaneously in the same TF matrix; that is, separate computations of TF were performed including these nodes separately.
TF matrices were thresholded to highlight the significant pairs of hubs. Significance was obtained by computing a reference null TF matrix obtained using every node outside of the core network following the same procedure described above for the hubs in the core network. The mean value of each null TF matrix obtained for each condition was used to threshold the corresponding TF matrices for the hub regions separately for the α and β band.
The similarity of the hub TF matrices across conditions was evaluated by means of the correlation coefficient and tested using the Mantel test. Finally, a control analysis tested whether two given hubs were highly central in the same epochs of high joint power; that is, whether the joint power fluctuations could predict joint centrality fluctuations. Therefore, we computed TF matrices according to Equations 4 and 5, where, for each hub pair: TP = presence of joint BC peaks and joint power peaks; TN = absence of joint BC peaks and joint power peaks; FP = (presence of joint BC peaks and absence of joint power peaks + absence of joint BC peaks and presence of joint power peaks)/2; and FN = FP.
Construction of low-level feature models.
Audiovisual stimuli can produce evoked activity that might affect the power correlation, especially in the α band. As a control analysis, for each 5 min movie and scrambled block, we extracted a set of low-level visual/auditory feature models and estimated the putative dependence of the α-BLP connectivity on these movie features. Low-level visual features included models related to luminance and contrast applied to the video image sequence and were computed automatically using MATLAB. Specifically, for each frame, the luminance was calculated from the average of the linear RGB components and the contrast, assigning one to edges in the input image and 0 s elsewhere. In addition, for each time-scrambled block, a model that accounted for the transient scene change was extracted by a human rater. This model consisted of a binary function that was 1 when a scene change was present and 0 elsewhere. Because scene segmentation in naturalistic movie sequences can be driven by top-down processing of conceptual features that are meaningful and correlated with ecologically valid relevant features of the action (Zacks and Swallow, 2007; Zacks et al., 2007), we did not consider this model within the set of low-level visual features. Moreover, our previous study (Betti et al., 2013) already showed a statistically significant relationship between the detection of event boundaries in the movie, as judged by an independent group of observers, and the α-BLP correlation time course, likely reflecting longer time-scale (tens of seconds) adjustments in functional connectivity induced by specific high-level sensory or cognitive events.
Low-level auditory feature included the root-mean-square (RMS) energy extracted from the sound track as in (Lahnakoski et al., 2012, 2017). Specifically, the slow fluctuations of the energy curve are very close to the envelope, making this signal an ideal candidate for the comparison with the slow fluctuations of the BLP connectivity.
To render the low-level models time courses consistent with the computed power correlation time course, we integrated the features models over time using the same parameters as in the correlation analysis (time window = 10 s sliding every 200 ms).
Relationship between BLP connectivity and the estimated features.
In each node pair, we linked the α-power correlation time courses with the time courses of each of the extracted features model using Pearson's correlation separately for each movie/scrambled block. This analysis led to a set of correlation matrices (r values). After Fisher z transformation of the r values, their significance was assessed via a two-tailed t test. p-values were corrected for multiple comparisons using the false discovery rate (FDR).
Results
Reduction of α-BLP but not β-BLP interaction strength during audiovisual stimulation
As a first step, we focused on changes of BLP coupling induced by the audiovisual stimulation to ultimately understand how this relates to modulation of network properties (topology and hubs). On average, audiovisual stimulation reduced the α-BLP connectivity both within and between networks. This effect was quantified by extracting the Fisher z-transformed correlation values in each subject network by network and averaging these values within networks and across all possible pairs of nodes from every possible pair of networks considered (between networks).
Specifically, in the α band, we found a significant main effect of condition in the DAN (F(2,64) = 4.43; p = 0.016, pη2 = 0.12), VIS (F(2,64) = 3.62; p = 0.032, pη2 = 0.11), SMN (F(2,64) = 3.48; p = 0.036, pη2 = 0.098), and a trend toward the significance in the AUD (F(2,64) = 2.97; p = 0.058, pη2 = 0.085). These effects were explained by a stronger α-BLP connectivity during fixation than movie or scrambled conditions (all p-values < 0.05). Moreover, in the α band, there was also a significant coupling type by condition interaction in the DAN (F(2,64) = 6.88 p = 0.002, pη2 = 0.18), VIS (F(2,64) = 4.69 p = 0.013, pη2 = 0.13), and a trend toward significance in the LAN (F(2,64) = 3.06 p = 0.053, pη2 = 0.08). In all cases (except LAN, in which between connectivity was reduced more than within), this effect was due to a reduction of both within- and between-RSNs connectivity in movie and scrambled compared with fixation (all p-values <0.05) (Fig. 3A).
Functional connectivity in the α and β bands. Pairwise BLP correlation in both bands during fixation, movie, and scrambled. The histograms show mean across subjects of the Fisher z-transformed correlation values averaged within and between all possible pairs of networks. Error bars indicate ± SEM. Asterisks denote significance at *p < 0.05, **p < 0.01. Duncan test was used for post hoc comparisons.
Interestingly, in the β band, no significant main effect of condition was found across networks (all p-values >0.05). However, there was a significant coupling type by condition interaction in DAN (F(2,64) = 6.17 p = 0.003, pη2 = 0.16), VIS (F(2,64) = 6.18 p = 0.003, pη2 = 0.16), LAN (F(2,64) = 4.32 p = 0.017, pη2 = 0.12), SMN (F(2,64) = 6.37 p = 0.002, pη2 = 0.16), and DMN (F(2,64) = 3.82 p = 0.027, pη2 = 0.11). In some networks (DAN, VIS, LAN), the effects were similar to those found in α-BLP. However, in other networks (SMN, DMN), the reduction of connectivity in scrambled was not significant (Fig. 3B). In other words, compared with rest, the modulation of BLP connectivity was more robust and consistent in α-BLP than in β-BLP. In β-BLP, the variations in connectivity seem to be more sensitive to the type of movie. Finally, in each network, there was a significant main effect of coupling type both in the α band and the β band (all p-values p <0.001; all pη2 values <0.73), which was accounted for by a stronger within- than between-network connectivity, consistent with the presence of distinct α-/β-BLP correlation networks of brain regions.
In summary, these results indicate that the audiovisual stimulation, both movie and scrambled conditions, strongly reduces the strength of α-BLP functional connectivity both within and between networks. Such a reduction during the observation of temporally scrambled movie scenes extends our previous results on regular movie scenes (Betti et al., 2013). Both results support the idea that audiovisual stimulation suppresses resting α-BLP functional connectivity consistently with the idling hypothesis. However, the connectivity strength in β appears on average to be more variable.
Preservation of large-scale network topology of β-band connectivity during movie observation
Next, using NBS, we statistically compared changes in resting α-/β-BLP connectivity topology induced by the audiovisual stimulation separately for natural and scrambled vision. Figure 4 shows the graph components that were significantly modulated (decreased) and the bar plot shows percentage of total number of within/between connections comprised in the altered component. In the α band, a very high proportion of between network connections decreased in both audiovisual conditions (gray lines on brain). A significant, albeit much lower, number of connections also decreased in the scrambled condition within the visual network (colored lines). In the β band, we found an intermediate number of decreased connections both within (mainly in the visual network) and between networks. Interestingly, there was no significant change of topology for the natural movie condition (Fig. 4).
Reorganization of intrinsic network topology during audiovisual stimulation. NBS shows task-induced between-networks decreases of graph components (gray lines) in the α and β bands during the audiovisual stimulation with respect to fixation. In the β band, the comparison movie versus fixation does not change the overall topography and for this reason is not shown. Only few within-network connections decreased (colored lines). Histograms show the percentage of decreased within- and between-RSNs connections for each band and condition.
The overall maintenance of within-network structure was also quantified through the Mantel test (Glerean et al., 2016) using 10,000 permutations network by network. Specifically, in each subject, we measured the correlation between the connectivity submatrices representing each RSN obtained at rest and during the audiovisual stimulation. Results expressed as percentage of runs displaying significant correlation indicated a strong similarity of within-network topology among fixation, movie, and scrambled conditions in both bands (e.g., β band within-VIS, 100% and 100% out of 33 runs for movie and scrambled, all p < 0.05 Bonferroni corrected; Table 2).
Maintenance of within-network structure
Overall, we can conclude that audiovisual stimulation did not change within-network connection topology (even though it decreased the overall BLP connectivity strength), but it decreased between-network α-BLP connections in both audiovisual conditions, hence causing a change in the overall topology. In contrast, between-network connections in the β band were significantly decreased only in the scrambled condition. Therefore, the number and structure of functional connections in the β band were not significantly altered between visual fixation and natural movie observation.
As a control, we investigated whether differences in functional connections between movie and scrambled conditions reflected differences in signal-to-noise ratio compared with fixation. For each band and node, we estimated the relative (as percentage) mean (averaged across time) power changes during movie and scrambled with respect to fixation. In the β band, we observed a similar power reduction between movie and scrambled (p = 0.31, paired t test). The same result was obtained in the α band (p = 0.36, paired t test). These results indicate that natural and scrambled movie sequences yielded similar α/β power decrements, yet different patterns of connectivity.
Regions of high centrality (hubs) during audiovisual stimulation and rest
The previous results show frequency-specific modulations of between-network connectivity during movie observation. In previous work, we showed that between-network interactions depend on a set of cortical regions (dynamic core network) that are highly central both in static and dynamic analyses (de Pasquale et al., 2016). Here, we consider the effects of audiovisual stimulation on the dynamic core network. We define as “hub” a node with BC above a threshold computed as the mean of the BC values for each condition obtained from a population of random graphs (Fig. 5A, solid line).
Measures of BC during rest and movie viewing. A, BC values in α-BLP and β-BLP for fixation (cyan), movie (magenta), and scrambled (green). These values are compared with the averaged BC obtained from a set of equivalent random graphs (solid lines). The patterns of node centrality vary between the α and β band. B, BC consistency (right y-axis, white bars) in each band defined as the number of times across the three experimental conditions in which the BC in each node exceeds the random graph BC threshold. The product BC consistency × BC values averaged across the three conditions (left y-axis, black bars) provides a core network for each band, defined as a group of nodes with consistently high centrality (red bars). C, For each band, the spatial location of the hubs form the core network.
To identify cortical regions of high centrality, we considered the product between BC consistency and mean BC across experimental conditions. Specifically, “consistency” is defined as the number of times across fixation, movie, and scrambled conditions in which the BC of a node exceeded the threshold for each condition (Fig. 5B, white bars). The product between BC consistency and mean BC across conditions (Fig. 5B, black bars) returned a few regions highly central both at rest and during the audiovisual stimulation (Fig. 5B, red bars). In the α band, these regions included right posterior intraparietal sulcus (PIPs) (DAN), right ventral frontal cortex (VFC) (VAN), right putamen (PUT) (SMN), right V4 (VIS), and left mPFC/mPFC2 (DMN). In the β band, a similar set of regions was found: right PIPs, right VFC, right PUT, and left mPFC/mPFC2 (as in the α band), but also right supplementary motor area (SMN) and right mPFC (DMN). These regions approximately overlap with the dynamic core network identified previously in the resting state (de Pasquale et al., 2016). Because they were highly central and consistent across conditions, we refer to these two groups of nodes, one for each band, as the “core network.” Figure 5C illustrates the spatial location of these nodes, which appear to be strongly overlapping in α and β. In summary, a set of hub nodes identified at rest did not change when examined during audiovisual stimulation both in α and β bands even though connection topology either changed (α and β band, in scrambled) or remained the same (β band in movie) during audiovisual stimulation compared with fixation.
Similarity of the dynamics of core network in the β band at rest and during movie observation
In the resting state, the centrality of nodes belonging to the core network fluctuates over time; that is, they alternate periods of weak and strong centrality with periods of joint fluctuations, as suggested previously (de Pasquale et al., 2016). Here, we tested whether the joint fluctuations of centrality among nodes of the core network are more similar during rest and the observation of movie scenes compared with the observation of time-scrambled scenes. This would be consistent with resting integration dynamics being a temporal prior for natural vision. Alternatively, centrality fluctuations in the resting state may reflect just idling spatiotemporal structures that are reorganized during audiovisual stimulation regardless of whether the stimuli are natural or synthetic.
As a first step, to determine the influence of audiovisual stimulation on BC dynamics, we measured the spectral structure of BC fluctuations in the different conditions by computing, for each node in the core network, the PSD of the BC time course. No significant differences were found between conditions, indicating that the BC fluctuations occur with a similar temporal structure (all p-values >0.05). When we binarized the BC, we confirm prominent joint fluctuations of centrality in core network nodes (see Fig. 2B for a sample of previously identified hub nodes).
The analysis performed through logical operators showed that several node pairs in the α and β bands display significant joint BC dynamics. Finally, we tested whether/how the audiovisual stimulation affected joint dynamics in the core regions by measuring the similarity of the TF patterns between rest and task. In the α band, the TF pattern, and thus the joint dynamics of centrality, differed between fixation and movie, fixation and scrambled (all p-values >0.05). In contrast, in the β band, such fluctuations were more similar between rest and movie (rest vs movie: p-values <0.05) than between rest and scrambled conditions (p-values >0.05) (Fig. 6). Importantly, for each pair of nodes, there was no significant relationship between joint BC and joint power fluctuations (Fig. 6, BC-power).
Matrices of TF positive rate of the joint BC fluctuations in α-BLP (A) and β-BLP (B). These matrices quantify how the dynamics of a hub predicts another hub in each condition. Black dots display significant node pairs. Similarity of BC dynamics across conditions was tested with the Mantel test. TF matrices testing if the occurrence of joint peaks of power could predict joint peaks of BC are also displayed as a control. Neighbor hubs in the core network (e.g., right/left mPFC) were not considered jointly in the same TF matrix (i.e., separate computations of TF were performed including these nodes separately).
As a control analysis, we tested whether BC fluctuations in the core network could predict fluctuations of centrality in sensory regions and vice versa. In fact, in our previous study, we showed that the dynamics of the visual regions depends on the sensory inputs (Betti et al., 2013). Results showed no significant relationship between fixation and both movie and scrambled conditions (p-values >0.05) (Fig. 7).
Matrices of TF positive rate of the joint BC fluctuations for hubs in the core network and for nodes in the visual networks (non-hubs) during fixation, movie, and scrambled in the β band. Black dots mark significant node pairs.
Because audiovisual stimuli can produce evoked activity that might affect the power correlation, especially in the α band, and thus the connectivity-based estimation of centrality, as a control analysis, we extract a set of low-level visual/audiovisual feature models and estimated the putative dependence of the α-BLP connectivity on these movie features for each movie/scrambled block and node pair. Results showed no significant relationship between α-BLP connectivity and feature models for each node pair both for the movie and scrambled condition (all p-values > 0.05, FDR corrected. Correlation values ranged from −0.33 to 0.23. Mean values were as follows: movie-luminance = ∼0; movie-contrast = ∼−0.07; movie-RMS = ∼0.03; scrambled-luminance = ∼−0.02; scrambled-contrast = ∼−0.1; and scrambled-RMS = ∼−0.05).
In summary, these results show that, in the β band, not only is a group of cortical hubs involved at rest and during movie observation, but also that their dynamics, measured as joint fluctuations of centrality (BC), are highly similar. This result does not extend to artificial movie stimuli.
Discussion
We studied changes in static and dynamic functional connectivity measured with BLP correlation at rest and during the observation of normal or time-scrambled movie clips. In the α band, audiovisual stimulation, compared with rest, reduced BLP connectivity, reorganized the overall topology (especially the between-network connections), and changed the dynamics of nodal centrality. In contrast, the β-band intrinsic topology and dynamics of integration at rest were similar to those measured during the observation of normal movie sequences, but altered for time-scrambled movies.
Different roles for α-BLP and β-BLP functional connectivity
α-/β-BLP connectivity is the main physiological correlate of fMRI connectivity at rest (Mantini et al., 2007; de Pasquale et al., 2010, 2012; Brookes et al., 2011; Hipp et al., 2012; Betti et al., 2013; Marzetti et al., 2013). The behavioral significance of intrinsic connectivity remains unclear. This study suggests key differences in the functional significance of α-BLP and β-BLP connectivity. Although α-BLP connectivity strength and topology were significantly decreased during audiovisual stimulation, there was an overall preservation of β-BLP connection topology and variable changes in connectivity strength during movie observation.
EEG/MEG studies during visual attention or visual processing associate increments of α-oscillations/power to cortical inhibition and decrements of α-oscillations/power to enhanced cortical excitability and information processing (Thut et al., 2006; Siegel et al., 2008; Jensen et al., 2012; Sauseng et al., 2013). Moreover, notwithstanding differences in spatial scale, desynchronization of α-BLP activity during audiovisual stimulation recalls decreases of spatial correlation of local field potentials observed in monkey and cat primary visual cortex during visual stimulation (Smith and Kohn, 2008; Nauhaus et al., 2009, 2012; Sato et al., 2012). Overall, these studies are consistent with the idea that active processing and sensory stimulation decrease correlated noise in the cortex (Smith and Kohn, 2008; Nauhaus et al., 2009, 2012; Sato et al., 2012). Our study shows that audiovisual stimulation induces, not only reductions of within/between network α power connectivity, likely akin to cortical noise, but also modulations of the overall topology, especially between networks, whereas within-network topology is largely preserved. Our interpretation is that the reduction of α-rhythm resting connectivity reflects a relative breakdown of idling states between networks that are reconfigured in a task-specific manner, as shown by fMRI and MEG studies (Nir et al., 2006; Cole et al., 2014; Krienen et al., 2014; Spadone et al., 2015; Gratton et al., 2016; Kim et al., 2017).
An incomplete but intriguing picture is emerging from this body of work. Strong interhemispheric connections at rest weaken during visual processing, as shown in fMRI connectivity (Spadone et al., 2015) and α-BLP correlation (this study and Betti et al., 2013). Correspondingly, there is a modulation of intrahemispheric connections between networks, as manifested in fMRI increases (Cole et al., 2014; Spadone et al., 2015; Kim et al., 2017), β-/γ-BLP correlation increases and α-BLP correlation decreases (this study and Betti et al., 2013). Disruption of interhemispheric connections is the main effect of stroke on fMRI resting connectivity which correlates with multiple cognitive deficits (Siegel et al., 2016). Therefore, an important link may exist between fMRI and α-BLP interhemispheric connectivity, their integrity at rest, and their decrement during normal task processing.
Changes in α-BLP connectivity cooccurred with a relative invariance of β-BLP topology and more variable changes in BLP connectivity strength during movie viewing. This is in contrast to significant changes in connectivity patterns (within/between networks) observed during viewing of scrambled movies. Changes during scrambled viewing may represent different factors, but more likely represent a desynchronization of visual cortex induced by the high temporal frequency of the visual stimuli (range 250–1000 ms duration). These changes were mainly limited to the visual network, which is consistent with the sensitivity of early visual areas to higher temporal frequencies (Hasson et al., 2008).
In summary, the intrinsic β-band connectivity pattern resembles the one during natural visual stimulation. This similarity in β-band connectivity is reinforced by the analysis of hub dynamics.
Dynamics of core network in the β band are similar between rest and natural vision
A second notable result was that a set of cortical regions appeared consistently central across conditions. Moreover, joint fluctuations of centrality across these hubs in the β band were significantly similar during visual fixation and movie, not scrambled, viewing. This result was specific to the β band because no significant similarity was observed between rest and both movie conditions in the α band. Control analyses indicated that joint fluctuations of centrality were not due to power changes nor did they cooccur with centrality fluctuations in early visual regions.
Our interpretation is that, in contrast to α-band connectivity, which reflects changes in the inhibitory tone of cortical circuitry, β-BLP topology and its dynamics underlie intrinsic signals that are more task related, perhaps reflecting the timing of visual stimulation or semantically related information. Moreover, β-BLP connectivity is not dependent on the visual input; for example, evoked responses. Therefore, β-band connectivity during natural vision seems intrinsically driven. Prior fMRI work on temporal sensitivity to natural scenes showed that early visual cortex is sensitive to short temporal patterns, whereas higher-order visual regions tend to prefer longer temporal stimuli that are more consistent and predictable in time (Hasson et al., 2008). Another study showed modifications of functional connectivity within visual occipital cortex and of (directional) interactions between visual and dorsal attention regions during a visual attention task compared with rest. Moreover, in higher-order dorsal attention regions, the strength and directionality of connections at rest predicted those measured during visuospatial attention (Spadone et al., 2015). As noted above, in contrast to α-band connectivity and dynamics, which are reorganized going from rest to visual stimulation, β-band topology and its dynamics appear to be preset at rest for natural vision. This rather surprising, yet insightful, result supports the functional significance of temporal dynamics of intrinsic activity in associative cortical regions.
A stronger similarity in connectivity for visual fixation and natural visual stimuli compared with scrambled visual stimuli is also consistent with recent work showing that resting connectivity in human visual cortex is more similar to that induced by natural scenes than synthetic stimuli (Wilf et al., 2017). Similarly, at level of microcircuits in animals, the mean and variability of neural responses at rest approximate the statistics of natural stimuli, not synthetic ones (Fiser et al., 2004; Berkes et al., 2011). However, a greater understanding of rest- and movie-related connectivity patterns will require a more in-depth analysis of the semantics contained in movie stimuli (Huth et al., 2012).
The similarity between topology and its dynamics during natural visual stimuli and intrinsic patterns is consistent with the idea that brain activity adapts to the statistical properties of the environment to which it is exposed in the course of development and experience (Simoncelli and Olshausen, 2001). Repeated visual stimulations produce reverberation of sensory evoked cortical activity in spontaneous activity (Han et al., 2008); therefore, repeated exposure to natural stimuli may shift, through learning mechanisms (Hebb, 1949), not only in space but also in time-intrinsic functional connectivity patterns.
Why are rest and natural vision topology and dynamic integration similar?
One influential hypothesis is that the function of intrinsic activity is to provide a spatial prior in connectional space from which task-evoked patterns are both generated and constrained (Raichle, 2011). However, behavior and cognitive states also unfold over time (Dehaene, 1993; Busch and VanRullen, 2010). For instance, extracting the meaning of a scene in a movie takes time, from seconds to tens of seconds (Kurby and Zacks, 2008). The slow-varying timescale of natural stimuli resembles that of neural fluctuations in associative cortex (Hasson et al., 2008; Honey et al., 2012). It is also possible that spontaneous fluctuations of integration might be tuned to the slow timescales distinctive of naturalistic stimuli and cognition.
Our current findings and previous ones (de Pasquale et al., 2012, 2016) show that fluctuations of centrality in hubs occur on a timescale of seconds to tens of seconds. The slow dynamics of integration might then represent a suitable candidate to encode the slow-varying temporal structure of naturalistic scenes. Through fluctuations of centrality, hubs may accumulate and integrate information across multiple systems. In the course of development and through experience, task-relevant temporal patterns of interaction become preserved in spontaneous activity. This is metabolically expensive (Attwell and Laughlin, 2001; Lennie, 2003), but computationally efficient. During movie observation, joint centrality fluctuations, which reflect primarily intrinsic dynamics, may anticipate and thus precompute the temporal structure of natural stimuli. This hypothesis is consistent with active predictive models of brain function emphasizing a role of the spontaneous activity in maintaining active representations (Engel et al., 2001). It also predicts a relationship between temporal anticipation of natural stimuli and slow-varying fluctuations of higher cortical areas. These higher-order regions may integrate and compare over time discrete snapshots of sensory information taken in sensory cortex (VanRullen and Koch, 2003; Baumgarten et al., 2015; VanRullen, 2016). However, further studies are needed to test this hypothesis more directly.
Footnotes
This work was supported by European Community's Seventh Framework Programme 5 (FP7/2007-2013), Grant Agreement HEALTH-F2-2008-200728 BrainSynch, and by the Human Connectome Project (Grant 1U54MH091657-01). V.B. was supported by a fellowship from the University of Chieti. M.C. was supported by the National Institutes of Health (National Institute of Mental Health Grant R01 MH096482-01 and National Institute of Neurological Disorders and Stroke Grant R01 NS095741).
The authors declare no competing financial interests.
- Correspondence should be addressed to Dr. Viviana Betti, Department of Psychology, Sapienza University of Rome, Via dei Marsi, 00185 Rome, Italy. viviana.betti{at}gmail.com
