# Lag threads organize the brain’s intrinsic activity

See allHide authors and affiliations

Contributed by Marcus E. Raichle, February 27, 2015 (sent for review December 16, 2014; reviewed by Gyorgy Buzsáki, Olaf Sporns, Jonathan D. Victor, and Rafael Yuste)

### This article has a Correction. Please see:

- Correction for Mitra et al., Lag threads organize the brain’s intrinsic activity - December 22, 2015

## Significance

It is well known that slow intrinsic activity, as measured by resting-state fMRI in a variety of animals including humans, is organized into temporally synchronous networks. The question of whether intrinsic activity contains reproducible temporal sequences has received far less attention. We have previously shown that human resting-state fMRI contains a highly reproducible lag structure. Here, we demonstrate that this lag structure is of high dimensionality and consists of multiple highly reproducible temporal sequences, which we term “lag threads.” Moreover, we demonstrate that the well-known zero-lag temporal correlation structure of intrinsic activity emerges as a consequence of lag structure. Thus, lag threads may represent a fundamental and previously unsuspected level of organization in resting-state activity.

## Abstract

It has been widely reported that intrinsic brain activity, in a variety of animals including humans, is spatiotemporally structured. Specifically, propagated slow activity has been repeatedly demonstrated in animals. In human resting-state fMRI, spontaneous activity has been understood predominantly in terms of zero-lag temporal synchrony within widely distributed functional systems (resting-state networks). Here, we use resting-state fMRI from 1,376 normal, young adults to demonstrate that multiple, highly reproducible, temporal sequences of propagated activity, which we term “lag threads,” are present in the brain. Moreover, this propagated activity is largely unidirectional within conventionally understood resting-state networks. Modeling experiments show that resting-state networks naturally emerge as a consequence of shared patterns of propagation. An implication of these results is that common physiologic mechanisms may underlie spontaneous activity as imaged with fMRI in humans and slowly propagated activity as studied in animals.

Spontaneous (intrinsic) neural activity is ubiquitously present in the mammalian brain, as first noted by Hans Berger (1). Spontaneous activity persists in all physiological states, although the statistical properties of this activity are modified by level of arousal and ongoing behavior (2⇓⇓⇓⇓⇓–8). Invasive studies in animals using diverse techniques—for example, local field potentials (9⇓–11), voltage-sensitive dyes (12⇓⇓–15), and calcium imaging (4, 16, 17)—have demonstrated richly organized intrinsic activity at multiple temporal and spatial scales. The most used technique for studying whole-brain intrinsic activity in humans is resting-state functional magnetic resonance imaging (rs-fMRI). Biswal et al. first reported that slow (<0.1 Hz) spontaneous fluctuations of the blood oxygen level-dependent (BOLD) signal are temporally synchronous within the somatomotor system (18). This basic result has since been extended to multiple functional systems spanning the entire brain (19⇓⇓–22). Synchrony of intrinsic activity is widely referred to as functional connectivity; the associated topographies are known as resting-state networks (RSNs) (23) and, equivalently, intrinsic connectivity networks (24).

Almost all prior rs-fMRI studies have used either seed-based correlation mapping (25) or spatial independent components analysis (sICA) (26). Critically, neither or these techniques provide for the possibility that activity within RSNs may exhibit temporal lags on a time scale finer than the temporal sampling density. However, we recently demonstrated highly reproducible lags on the order of ∼1 s by application of parabolic interpolation to rs-fMRI data acquired at a rate of one volume every 3 s (*SI Appendix*, Fig. S1) (27). Moreover, this lag structure can be modified, with appropriate focality, by a variety of task paradigms (27).

Investigations of rs-fMRI lag structure previously have been limited by the concern that observed lags may reflect regional differences in the kinetics of neurovascular coupling rather than primary neural processes (28, 29). However, our previous dimensionality analysis demonstrated that there are at least two independent lag processes within the brain (27). The neurovascular model can account for only one of these. Hence, there must be at least one lag process that is genuinely of neural origin. We have since made significant methodological improvements (*Theory* and Fig. 1) that enable a more detailed characterization of lag structure in BOLD rs-fMRI data. We report our results in two parts.

In part I, we present an expanded view of the lag structure within the normal adult human brain derived from BOLD rs-fMRI data in 1,376 individuals. Specifically, we show that at least eight orthogonal lag processes can be reproducibly demonstrated. We refer to these processes as “threads” by way of analogy with modern computer programming practice in which single applications contain multiple, independent thread sequences.

In part II, we investigate the relation between lag threads and zero-lag temporal correlations—that is, conventional, resting-state functional connectivity. We find that, although there is no simple relation between lag and zero-lag temporal correlation over all pairs of voxels, apparent propagation is largely unidirectional within RSNs. We also show that the zero-lag temporal correlation structure of rs-fMRI arises as a consequence of lags, whereas the reverse is not true. These results suggest that lag threads account for observed patterns of zero-lag temporal synchrony and that RSNs are an emergent property of lag structure.

## Theory

We define the lag between two fMRI time series by computing the cross-covariance function at intervals of one frame and identifying the local extremum using parabolic interpolation (*SI Appendix*, Fig. S1). This analysis assumes the existence of a single temporal lag between regions. The validity of this assumption depends on the fact that BOLD fMRI time series are aperiodic (30, 31) (see *SI Appendix*, *Lagged Cross-Covariance Curves Exhibit a Single Peak* for additional discussion of this point). Measured lags at the group level (i.e., averaged over individuals) typically assume values in the range ±1 s. Apparent propagation is inferred on the basis of observed lag between two time series. This formulation makes no assumptions regarding the path over which the activity “propagates” between regions. Thus, “propagation,” as defined here, entails lags on the order of ∼1 s in activity over spatial scales on the order of centimeters.

As an aid to understanding the methodology, we describe our approach to characterizing lag structure using a simple illustrative model containing three orthogonal lag processes (threads) propagating through six nodes (Fig. 1). Apparent propagation, as defined here, is shown using synthetic time series with “1/f” spectral content duplicated from real BOLD rs-fMRI data (31) (see *SI Appendix*, *Simulating Synthetic BOLD fMRI Time Series* for further detail) propagating through six nodes (Fig. 1*A*). The superposition of the three thread processes is shown in Fig. 1*B*. Analysis of the superposed time series observed at the six nodes (using the procedure illustrated in *SI Appendix*, Fig. S1) yields the time delay matrix shown in Fig. 1*C*; we call this matrix *TD*. Having computed *TD*, we can compute a mean for each column, using the previously described projection strategy for computing BOLD rs-fMRI lag topographies (27). In the case of a single lag process, the projection strategy is sufficient to recover the ordering (*SI Appendix*, Figs. S2 and S3). However, as shown in Fig. 1*D*, in the case of multiple superposed lag processes, the column-wise projection generates an oversimplified approximation of the dynamic system.

In the present work, to recover lag processes in multidimensional time series, we use principal components analysis (PCA). PCA assumes linear superposition of components. The validity of this assumption is discussed in *SI Appendix*, *Validity of Applying PCA to Recover Lag Thread Topographies*. In each column, *TD*, the vector corresponding to column *i*, is a lag map of the system with reference to time series *i*. That is, the first column of *TD*_{z}, constructed by subtracting the mean of each column from *TD*. The columns of *TD*_{z} are zero-centered lag maps. Application of PCA to *TD*_{z} recovers the eigenspectrum representing the number of lag threads present in the system. Fig. 1*E* shows that precisely three nonzero eigenvalues are found in this illustrative case. The eigenvectors corresponding to these nonzero eigenvalues can be used to recover the topography of the lag threads; the node diagrams above the nonzero eigenvalues in the *Lower Right* panel of Fig. 1 illustrate the recovered lag processes. In the case of no delays (*SI Appendix*, Fig. S2) or only a single set of delays (*SI Appendix*, Fig. S3), PCA finds zero or one nonzero eigenvalue. Thus, *TD* analysis is sufficient to assess the number of lag threads in the system. Although Fig. 1 illustrates *TD* and *TD*_{z} as square matrices (i.e., the number of voxels in each lag map is equal to the number of lag maps), lag thread computation is algebraically well defined also when the number of voxels greatly exceeds the number of lag maps.

To increase the signal-to-noise ratio (SNR) in real BOLD rs-fMRI data, we produced (6 mm)^{3} voxel resolution lag maps from time series extracted from 330 (15 mm)^{3} cubic regions of interest (ROIs), uniformly distributed throughout gray matter (see *SI Appendix* for further detail).

## Methods

A large data set (*n* = 1,376) was obtained from the Harvard-MGH Brain Genomics Superstruct Project (32) (Table 1). The 1,376 subjects were randomly divided into two groups of 688 subjects to test the reproducibility of our analyses. Please see *SI Appendix* for further details regarding preprocessing and computational methods.

## Results

### Part I.

#### Existence and reproducibility of lag threads.

Fig. 2 shows the topography of four lag threads derived from real BOLD rs-fMRI data acquired in the first group of 688 subjects. Blue hues indicate regions that are early—that is, “sources” of propagated BOLD activity—and red hues indicate regions that are late—that is, “destinations” of propagated BOLD activity. The range of lag values in each thread is ∼2 s. The threads exhibit a high degree of bilateral symmetry. Interestingly, although specific anatomical structures are often prominent sources or destinations within threads, these topographies do not respect RSN boundaries (Fig. 2 and Movies S1–S4). BOLD rs-fMRI signals propagate in lag threads both within and across RSNs. We have previously reported this principle in relation to the distribution of lag values over pairs of nodes in the brain (27).

Fig. 3*A* shows the first 20 eigenvalues derived by spatial PCA of the lag threads derived from the first 688-subject group (the first four threads of which are illustrated in Fig. 2). In Fig. 1*E*, only nonzero eigenvalues correspond to lag threads, making it easy to infer the dimensionality of the system. In real data, the presence of noise means that all eigenvalues are nonzero; hence, dimensionality must be estimated (33). Using an information criterion, we estimated the dimensionality in Fig. 3*A* to be 8 (see *SI Appendix* for further detail). To explore the reproducibility of our results, we applied the same calculations to a second, separate group of 688 subjects and obtained eigenspectrum and dimensionality estimates indistinguishable from the results shown in Fig. 3*A*. Reproducibility of lag thread topographies across the two groups of 688 subjects is illustrated in Fig. 3*B*. It is evident that lag threads are highly reproducible across groups. The full topographies of the first eight lag threads are reported in *SI Appendix*, Figs. S7–S14.

#### Multiple lag threads in a selected number of regions.

The high dimensionality of the lag system shown in Fig. 3 theoretically could reflect our specific choices of ROIs and lag map resolution (*SI Appendix*, Fig. S5). To test this possibility, we constructed 17 ROIs by thresholding the eight thread maps to define prominent sources and destinations. These ROIs are shown in Fig. 4*A*. Many of these regions have been previously identified as critical nodes that organize the brain’s ongoing activity (34). The corresponding 17 × 17 time delay matrices in the two groups of 688 subjects are each shown in Fig. 4*B*. Excellent reproducibility is evident (Fig. 4*C*). Moreover, the maximum likelihood dimensionality estimate in both groups is 8 (Fig. 4*D*), precisely the same result shown in Fig. 3*A*. Fig. 4 demonstrates that high dimensional lag structure can emerge from a relatively small set of ROIs as long as key regions of the brain are represented. Seed-based lag maps for these ROIs are shown in *SI Appendix*, Figs. S15–S24.

### Part II.

#### Lag threads in relation to zero-lag temporal correlation.

Having demonstrated the existence and reproducibility of multiple lag threads in human rs-fMRI data, we next investigated the relation of lag threads to zero-lag temporal correlation structure. To this end, we examined shared patterns of propagation across lag threads. We did this by defining a matrix, *k* is the number of lag threads (*SI Appendix*, Eq. **S9**). The corresponding 6526 × 6526 (voxels *A*) reveals commonalities in signal propagation across lag threads (see *SI Appendix* for additional discussion). Critically, the voxel-wise correlation structure across lag threads (Fig. 5*A*) resembles the conventional zero-lag temporal correlation structure of BOLD rs-fMRI (Fig. 5*B*) (RSN membership as in ref. 35). We quantitatively confirmed the similarity between the correlation structures in Fig. 5 *A* and *B* by computing the Pearson correlation between the unique values in each matrix (*r* = 0.41). Thus, there is a correspondence between common patterns of propagation across lag threads (Fig. 5*A*) and the zero-lag temporal structure of rs-fMRI.

The similarity between the matrices shown in Fig. 5 *A* and *B* raises the question of whether, across pairs of voxels, lag is directly related to zero-lag temporal correlation. Hypothetically, voxel pairs exhibiting shorter temporal lags could exhibit higher temporal correlation at zero lag. However, a 2-dimensional histogram of zero-lag Pearson *r* versus lag for all voxel pairs (Fig. 5*C*) demonstrates no systematic relation. A second possibility is that regions belonging to the same RSN are iso-latent (i.e., exhibit the same lag value) in each lag thread. For example, the entire default mode network could be early in one thread but late in another. However, the topographies in Fig. 2 show that this is not the case. Additionally, we have previously shown that RSNs are not iso-latent. Rather, the range of intra- and inter-RSN lag values is the same and no RSN is either early or late as a whole (figure 9 in ref. 27). Therefore, the similarity between the cross-thread correlations (Fig. 5*A*) and zero-lag BOLD rs-fMRI temporal correlations (Fig. 5*B*) requires a more subtle explanation. To provide the basis for this explanation, we introduce the concept of “lag thread motifs.”

#### Lag thread motifs.

We illustrate the concept of a lag thread motif using a four-voxel model system with two lag threads (Fig. 6). The overall pattern of propagation between the four voxels is different in the two threads (Fig. 6*A*). However, the sequence of propagation through voxels 1 and 2 is identical. Voxels 1 and 2, therefore, constitute a lag thread motif: a set of regions in which the sequence of propagation is the same across lag threads. There are no other motifs in Fig. 6. The patterns of propagation shown in Fig. 6*A* are realized in Fig. 6*B* using synthetic time series with 1/f spectral content duplicating that of real BOLD rs-fMRI data (31) (*SI Appendix*, *Simulating Synthetic BOLD fMRI Time Series*). Again, although the overall pattern of propagation differs between threads, the sequence of propagation between the first two voxels is preserved (dotted arrows in the dark red boxes).

We asked whether the thread motif model can explain the findings in Fig. 5, specifically, the similarity between the cross-thread correlations (Fig. 5*A*) and zero-lag temporal correlations (Fig. 5*B*), and the absence of a systematic relation between zero-lag Pearson *r* versus lag for all voxel pairs (Fig. 5*C*). We explore this question in Fig. 7, which presents a simulation experiment based on the model presented in Fig. 6, but scaled up to include 30 voxels, eight orthogonal lag threads, and two thread motifs (see *SI Appendix* for details concerning generation of orthogonal model lag threads; *SI Appendix*, Fig. S6 shows an explicit description of the model). Motif 1 propagates through voxels 1–5; motif 2 propagates through voxels 6–10. Fig. 7*A* shows the voxel-wise correlation across simulated lag threads (paralleling Fig. 5*A*). Voxels sharing a thread motif necessarily exhibit perfectly correlated lag sequences (diagonal blocks labeled “1” and “2”). Fig. 7*B* shows the 30 × 30 zero-lag temporal correlation matrix (paralleling Fig. 5*B*) computed on the basis of the synthetic 1/f time series representing the lag threads (see *SI Appendix* for additional details). Thus, thread motifs are sufficient to induce lag thread correlation structure, as in real rs-fMRI data (Fig. 5*A*). Moreover, the matrices shown in Fig. 7 *A* and *B* (synthetic data) exhibit the same similarity as the matrices shown in Fig. 5 *A* and *B* (real data). This similarity suggests that the existence of shared lag thread motifs is sufficient to explain zero-lag temporal correlations. We note that this model depends on the 1/f spectral content of BOLD rs-fMRI time series. The existence of an association between lag structure and temporal correlation requires that the underlying time series exhibit some degree of autocorrelation; lag structure would be dissociated from correlation structure in a system in which the signals were comprised of infinitely narrow impulses or white noise. Finally, Fig. 5*C* indicates that there is no systematic relationship in real BOLD rs-fMRI data across voxel pairs between conventional zero-lag Pearson *r* and lag. This feature is present also in our simulation (Fig. 7*C*). The low correlation between conventional zero-lag Pearson *r* and lag in the synthetic data (*r* = –0.04) confirms that thread motifs need not introduce a systematic relation between these quantities. Therefore, the results of Fig. 7 (simulation) suggest that the results shown in Fig. 5 (real data) can be explained if thread motifs correspond to conventional RSNs—in other words, if intra-RSN sequences of propagation are preserved across threads.

#### RSNs correspond to lag thread motifs.

The simulation in Fig. 7 suggests a model in which conventionally defined RSNs correspond to thread motifs and implies two testable predictions. First, if the sequence of propagation is preserved within thread motifs, it follows that the dimensionality of intra-RSN lag structure should be 1 (see *SI Appendix*, Fig. S3 for further explanation). Second, although the simulation in Fig. 7 contains no systematic relation, over all voxel pairs, between lag and zero-lag temporal correlation (Fig. 7*C*), if we examine the relationship between zero-lag Pearson *r* and lag considering only voxels within a motif, a substantial negative correlation emerges (*r* = –0.75; Fig. 7*D*). The basis for this relation is that, within a single thread motif, more nearly synchronous time series must be more correlated. However, in general, threads propagate in multiple directions outside of motifs (e.g., as in Fig. 6). Consequently, relations of the type shown in Fig. 7*D* are obscured in Figs. 5*C* and 7*C*, because the fraction of intramotif voxel pairs is a small fraction of all voxel pairs. Therefore, if RSNs correspond to lag thread motifs, voxel pairs within an RSN should exhibit a substantial negative correlation between zero-lag Pearson *r* and lag (Fig. 7*D*).

We test these predictions in Fig. 8 *A* and *B*. Fig. 8*A* shows the maximum likelihood dimensionality of the temporal lag structure calculated within RSNs as defined in ref. 35. As predicted by the thread motif model, the maximum likelihood dimensionality is 1 for all RSNs except the sensorimotor network (SMN). Fig. 8*B* shows the Pearson correlation derived from a scatter plot of zero-lag temporal correlation versus lag over all intranetwork voxels (as in Figs. 5*C* and 7*C*). Again, as predicted by the thread motif model, there is a substantial negative correlation in every RSN except the SMN. Therefore, although lag threads represent various patterns of propagation with generally reciprocal signaling across regions, within each RSN, BOLD rs-fMRI signal propagation is largely unidirectional.

The apparently anomalous dimensionality result obtained with respect to the SMN (Fig. 8 *A* and *B*) highlights an interesting point concerning the correspondence between RSNs and thread motifs. Conventional BOLD rs-fMRI analyses generally agglomerate primary somatomotor and somatosensory areas into a single RSN (22, 36) (see also Fig. 5*B*). Separation of motor and sensory areas into distinct parcels has only recently been achieved using a boundary mapping technique (37). In contrast, lag threads sharply distinguish primary sensory versus primary motor cortices (Fig. 8*A*). Primary motor cortex is earlier than primary sensory cortex in most threads (Fig. 2, threads 3 and 4), but the ordering is reversed in other threads (*SI Appendix*, Figs. S11 and S12). This feature is also illustrated in Fig. 8*A*. Consequently, the observed dimensionality of lag structure in the SMN (as conventionally defined) is 2 (Fig. 8*A*). We verified that the dimensionality of lag structure in the separated motor and sensory components of the SMN is 1 in both cases.

#### Zero lag temporal synchrony emerges from lag structure.

The previous results suggest that RSNs correspond to lag thread motifs—that is, that the sequence of propagation within RSNs is largely preserved across lag threads. This finding raises the possibility that zero-lag temporal synchrony (i.e., conventional functional connectivity) within RSNs emerges from lag structure. To test this hypothesis, we converted each of the first eight lag threads extracted from real data (Fig. 2) into time series with the same spectral content as BOLD rs-fMRI (as in Fig. 6*B*; see *SI Appendix*, *Simulating Lag Threads* for further details). We then superposed these time series, weighted in proportion to their respective eigenvalues (Fig. 4*A* and *SI Appendix*, Eq. **S9**), to reconstruct synthetic BOLD rs-fMRI data with appropriate spectral content and imposed structure derived only from lag threads. The zero-lag temporal correlation matrix computed from the reconstructed time series is shown in Fig. 8*C*. This matrix is strikingly similar to the zero-lag temporal correlation matrix computed from real BOLD rs-fMRI data (Fig. 5*B*). Fig. 8*D* shows a scatterplot of the real (Fig. 5*B*) versus reconstructed (Fig. 8*C*) zero-lag temporal correlation values. Light blue and dark blue dots in Fig. 8*D* represent inter-RSN and intra-RSN correspondence, respectively. The scatterplot quantitatively demonstrates substantial agreement between the zero-lag temporal correlation structure of real and reconstructed data (*r* = 0.58). The model better predicts intra-RSN correlation structure (*r* = 0.62) versus inter-RSN correlation structure (*r* = 0.24); the implication of this difference is at present not understood. Nevertheless, Fig. 8 *C* and *D* suggests that intra-RSN synchrony RSNs is an emergent property of lag structure.

Fig. 8 demonstrates that zero-lag temporal correlation can arise from patterns in lag. Fig. 9 demonstrates that lag structure is not uniquely determined by zero-lag temporal correlation structure. To illustrate this point, we generated synthetic, multidimensional time series with spectral content and second order statistics (i.e., zero-lag correlation structure) matched to that of the real BOLD rs-fMRI data (Fig. 9 *A* and *B*; see *SI Appendix*, *Zero-Lag Temporal Correlation Need Not Specify Lag Structure* for details). Fig. 9*C* shows that the average zero-lag temporal correlation matrix, in both synthetic (blue) and real (red) data, converges to the average structure derived from real data. Fig. 9*D* concerns lag and demonstrates a very different result. In particular, the average time delay matrix computed over synthetic data rapidly converges to the all zeros matrix (blue line in Fig. 9*D*), whereas the *TD* matrix computed over real data converges to a consistent, nonzero delay structure (red line in Fig. 9*D*). Fig. 8*C* demonstrates that lag structure can be used to reconstruct zero-lag temporal correlations, whereas Fig. 9 shows that the reverse does not hold. The implication of these results is that lag structure represents the more fundamental level of organization in rs-fMRI.

## Discussion

### Summary of Findings.

The structure of human intrinsic brain activity, as imaged with resting-state BOLD fMRI, has been understood predominantly in terms of zero-lag, temporal synchrony within widely distributed functional systems (RSNs). We previously demonstrated that interregional lags are reproducibly present in BOLD rs-fMRI data and that these lags are not attributable to hemodynamic factors (27). We have substantially expanded on our previous findings here. In part I, we demonstrated that lag threads in human rs-fMRI exhibit multiple, highly reproducible patterns of propagated activity (lag threads) in BOLD rs-fMRI data (Figs. 2–4). We also showed that there are most likely eight lag threads in our current analysis and that eight lag threads can emerge from a small set of key ROIs (Figs. 3 and 4).

In part II, we investigated the relation between lag threads and zero-lag temporal correlations in BOLD rs-fMRI. To this end, we examined common patterns of propagation across lag threads by computing voxel-wise correlations across eight lag threads (Fig. 5*A*). We found that voxel-wise lag–thread correlations and BOLD rs-fMRI zero-lag temporal correlations exhibit a similar structure (Fig. 5). To explain this similarity, we hypothesize the existence of lag thread motifs—that is, sequences of propagation through subsets of regions that are shared across multiple threads (Fig. 6). Simulation experiments showed that zero-lag temporal synchrony within RSNs naturally emerges as a consequence of lag thread motifs (Fig. 7). We also demonstrated that conventionally defined zero-lag RSNs very likely correspond to lag thread motifs (Fig. 8 *A* and *B*). Finally, we reproduced, to a fair approximation, the zero-lag temporal correlation structure of BOLD rs-fMRI using synthetic time series with imposed structure derived only from lag threads (Fig. 8 *C* and *D*). The reverse relation does not hold—that is, zero-lag temporal correlation structure does not determine a unique lag structure (Fig. 9). Hence, temporal synchrony can be understood as a consequence of BOLD rs-fMRI lag structure.

It is important to note that, although the thread motif model provides a basis for understanding how synchrony arises from patterns in lag, this model is in no way imposed in the reconstruction in Fig. 8*C*; BOLD time series were reconstructed purely on the basis of lag threads calculated from real data. We conclude that lag threads represent a fundamental organizing property of the brain’s intrinsic activity.

### Physiology of Propagated Activity.

A prominent concern regarding lags in human rs-fMRI has been that all observed phenomenology could be attributable to regional variations in the kinetics of neurovascular coupling (28, 29, 38⇓–40). Such regional differences can be represented as an ordered sequence, as illustrated in Fig. 1. Importantly, this set of temporal shifts can account for only a single lag thread (27). Delayed signals in large venous structures most likely contribute to the topography of lag thread 5 (*SI Appendix*, Fig. S11). However, the existence of at least eight lag threads demonstrates that regional differences in neurovascular coupling account for a minor component of BOLD rs-fMRI lag structure.

Our results raise the question of what physiological mechanisms might underlie propagation of slow spontaneous activity over the whole brain. Propagation of low-frequency (<1 Hz) spontaneous activity has been extensively described in the rodent brain using various modalities, including whole-cell recordings (9), local field potentials (5, 41, 42), voltage-sensitive dyes (11, 12, 43), and calcium imaging (16). We have previously shown that spontaneous BOLD signal fluctuations correspond to low-frequency (<1 Hz) local field potentials, also known as slow cortical potentials (SCPs), which represent slow endogenous changes in excitability (44, 45). Thus, we speculate that propagation of activity in the BOLD signal is likely to represent propagation of slow changes in neuronal excitability.

Propagated changes in neuronal excitability have been previously described in terms of UP/DOWN states (UDSs). UDSs are slow (<1 Hz), spontaneous, subthreshold changes in neuronal membrane potential. These membrane potential fluctuations are effectively synchronous at a submillimeter spatial scale (15, 46, 47) but exhibit multiple, complex patterns of propagation over larger spatial scales (3, 4, 9, 15, 41, 48, 49), spanning thalamus (41), striatum (50), and cortex (41). Although UDSs were initially associated with anesthesia and slow wave sleep, it is now known that UDSs persist and propagate during quiet wakefulness (9, 11, 15). The lags we found (on the order of ∼1 s over the whole brain) in human spontaneous activity are comparable to UDS propagation delays in rodents (9, 11, 12, 15, 41, 43, 51). There also are intriguing correspondences between the directionality of lags in the BOLD signal and lags in UDSs in rodents. Reports by Hahn et al. (9) and Sirota et al. (42) both document that slow fluctuations in somatosensory neocortical areas lead activity in hippocampus by less than 1 s. These findings agree with lags between hippocampus and somatosensory cortex in a seed-based lag map derived from entorhinal cortex (*SI Appendix*, Fig. S17).

Thus, UDSs might underlie BOLD rs-fMRI signal fluctuations (6, 44, 52). However, some features of these two phenomena are discrepant. In particular, UDSs are generally periodic with frequency content in the range of 0.5–0.8 Hz (9, 15, 46), whereas the resting-state BOLD fMRI signal is aperiodic and dominated by frequencies ≤0.1 Hz (31). Whether or not UDSs are responsible for the present results, we note that there is no consensus regarding the mechanisms underlying slowly propagated activity. Proposed explanations include shifts in excitatory–inhibitory balance (51), thalamo–cortical interactions (41), astrocytic gliotransmission (52, 53), and metabolic neuromodulators such as adenosine (51, 54). Future work is required to definitively elucidate the physiologic mechanisms underlying propagation in the BOLD rs-fMRI signal.

### Topography of Lag Threads.

Each lag thread (Fig. 2, Movies S1–S4, and *SI Appendix*, Figs. S7–S14) exhibits a unique, highly reproducible (Fig. 3*B*) topography. It is evident that these topographies generally are bilaterally symmetric, as are most functional systems and RSNs. Although fragments of functional systems can be observed in the lag threads (e.g., frontopolar cortex components of the frontoparietal control network in thread 1), the contours of the thread maps do not correspond to the topographies of RSNs. This point relates to our previous demonstration that RSNs as a whole are not iso-latent; rather, each RSN contains a wide range of latency values and includes both early and late nodes (27).

Several topographic features of the lag threads suggest functional properties. Because these points are speculative, we consider only the first two lag threads. In thread 1, the earliest regions are brainstem, thalamus, hippocampus, and putamen; late areas include frontopolar cortex and central insula. This sequence suggests a “bottom–up” process in which activity begins subcortically and propagates to progressively higher order areas of the cerebral cortex (Fig. 2 and Movie S1). Another striking finding is that the putamen is early, whereas the caudate is late. See Movie S1 to visualize activity propagation from posterior putamen to the head of the caudate. As far as we are aware, this finding has not been previously described. In thread 2, thalamus, hippocampus, and brainstem are late with respect to cerebral cortex—that is, opposite to their temporal position in thread 1. This sequence suggests a “top–down” process in which activity propagates from higher to lower centers (Fig. 2 and Movie S2). However, other features—for example, late cerebellum, early putamen, and late caudate—are common to threads 1 and 2. The functional significance of specific lag thread topographies is unknown and requires further study.

### The Role of Lag Threads in Integration Versus Segregation.

Activity in the brain, at various spatial scales, has been discussed in terms of two fundamental concepts: synchrony (8, 13, 55) and lagged propagation (3, 4, 41, 56, 57). Taken to their logical extremes, synchrony and lag are opposed in a simple system: A perfectly synchronous system contains no lags, and a system with a single set of lags is not synchronous (58) (*SI Appendix*, Fig. S3). The fact that the brain’s spontaneous activity exhibits both of these properties may be a manifestation of the dual functions of neuronal segregation and integration (27, 59). Conventional zero-lag resting-state functional connectivity analysis has provided a powerful tool for using synchronicity to map spatially distinct functional areas (19⇓⇓⇓–23, 37). However, functional parcellations do not explain how spatially segregated modules in the brain become integrated (59). Lag threads demonstrate that spontaneous activity exhibits apparent propagation both within and between spatially segregated RSNs. Therefore, lag threads may explain how spatially segregated networks can be integrated over a time scale of seconds.

Conversely, lag threads pose their own problem: If spontaneous activity is characterized by a lag structure, how does synchrony arise? Our results suggest that lag thread motifs provide an answer. Preservation of lag sequencing within certain regions of the brain (i.e., RSNs) across multiple threads gives rise to zero-lag synchrony within these systems (Fig. 8 *C* and *D*). Thus, the lag thread motif model unifies the coexistence of synchrony (spatial segregation) and lags (temporal integration) in the brain’s spontaneous activity.

The physiological functions served by lag threads remain unknown, but previous work sheds some light on this matter. We have shown that the lag structure of rs-fMRI is focally modulated, in humans, following the performance of a motor task (27), suggesting that the lag structure of intrinsic activity may be involved in learning and memory. Indeed, the spatiotemporal structure of UDSs (discussed in *Physiology of Propagated Activity*), a potential correlate of lags in the BOLD rs-fMRI signal, has been linked to consolidation and plasticity mechanisms (9, 15, 41, 60). Additional support for this perspective comes from studies of neurodevelopment showing that precise patterns of propagated intrinsic activity are essential for fine-tuning synaptic connections (61⇓–63). It is believed that persistence of this principle into adulthood supports the brain’s capacity for lifelong plasticity (63⇓⇓–66). A second hypothesis is that the cortex, like the spinal cord, acts as a central pattern generator (67, 68) and that patterns of propagated intrinsic activity represent neuronal programs that are recruited to perform tasks (4, 5, 69). Thus, lag threads may form the basis for activity sequences that naturally play out in responses to events.

Regardless of their specific functions, the reproducibility of lag thread phenomenology suggests that this organizational feature is essential to normal brain physiology and function. We hypothesize that perturbed lag thread structure may underlie some neuropathological conditions. If so, these conditions may not manifest as altered conventional functional connectivity, as changes in lag thread structure (for instance, altered thread hierarchy in Fig. 3*A*) may not change zero-lag temporal correlations. Hence, an understanding of the physiologic functions of lag threads may lead to better understanding of the brain in health and disease.

### Limitations and Future Directions.

The SNR of BOLD rs-fMRI is limited. Accordingly, extensive averaging over very large subject groups was required to obtain stable lag estimates at (6 mm)^{3} voxel resolution using (15 mm)^{3} reference ROIs. We are optimistic that future improvements in BOLD fMRI (70), for instance, increased temporal resolution, will allow detection of lag threads in smaller populations, provided that voxel-wise SNR remains adequate and preprocessing strategies effectively remove artifact. Alternative approaches to studying rs-fMRI lags, at coarser resolution, but with less sensitivity to SNR limitations, include deriving lag structure from selected ROIs, as in Fig. 4, and computing lag projections, as in ref. 27.

A second caveat is that the presently reported correspondence between lag thread motifs and RSNs (Fig. 8 *A* and *B*) reflects a specific RSN parcellation (*SI Appendix*, Fig. S5) (35), although the inferences derived from this result most likely depend only minimally on the details on any particular parcellation scheme (21, 22). We note that the correlation-based results (Figs. 5 *A* and *B* and 8 *C* and *D*) are parcellation independent.

Third, as lag threads are simply principal components of lag structure, they formally constitute only a basis set for lagged activity. Consequently, the sign of the lag threads are, by definition, undetermined by a factor of ±1. Moreover, the assumption of linear superposition in PCA implies that topologically complex or nonlinear temporal sequence topographies cannot be recovered. However, we did find that kernel PCA, a nonlinear technique, recovers lag thread topographies (*SI Appendix*, Fig. S28) quite similar to those shown in Fig. 2. Additionally, the sign and topographies of seed-based lag maps (*SI Appendix*, Figs. S15–S24) are uniquely determined. We used these maps to demonstrate that lag thread topographies reasonably separate seed-based lag maps into common clusters and that the sign of each lag thread has most likely been correctly assigned (*SI Appendix*, *Validity of Applying PCA to Recover Lag Thread Topographies*).

Fourth, there is an ambiguity concerning voxels with lag values near zero in each lag thread. One possibility is that these voxels are in the middle of the temporal sequence represented by the lag thread. Alternatively, the voxel may not participate in the temporal sequence. At present, we cannot distinguish between these possibilities.

Finally, Fig. 8*C* shows a temporal correlation matrix computed on the basis of reconstructed BOLD rs-fMRI time series. This result reproduces many features of real data (Fig. 5*B*), but the correspondence obviously is imperfect (Pearson *r* = 0.58; Fig. 8*D*). Importantly, we assumed that the spectral content of BOLD rs-fMRI is uniform over gray matter and that lag threads superpose linearly. These assumptions represent an approximation (71), although the extent to which spectral shapes are regionally dependent at frequencies below 0.1 Hz is uncertain (31). We also restricted our reconstruction to only the first eight lag threads deemed significant by maximum likelihood dimensionality analysis. Although the remaining lag threads contribute less individual variance, they may collectively play an important role in shaping correlation structure. In view of these approximations, a more complete model may be expected to provide a closer match between a reconstructed and true correlation structure of BOLD rs-fMRI time series. As our reconstruction relies only on lag threads, we have also excluded other phenomena that may contribute to the coordination of zero-lag correlation structure.

## Acknowledgments

We thank John McCarthy, Manu Goyal, and Carl Hacker for helpful discussion and Dr. Randy Buckner for assistance in obtaining the MRI data, which came from the Brain Genomics Superstruct Project, funded by the Simons Foundation Autism Research Initiative. This work was supported by National Institutes of Health Grants NS080675 (to M.E.R. and A.Z.S.), P30NS048056 (to A.Z.S.), and F30MH106253 (to A.M.).

## Footnotes

- ↵
^{1}To whom correspondence may be addressed. Email: marc{at}npg.wustl.edu or anishmitra{at}wustl.edu.

Author contributions: A.M., A.Z.S., and M.E.R. designed research; A.M., A.Z.S., T.B., and M.E.R. performed research; A.M., A.Z.S., and M.E.R. contributed new reagents/analytic tools; A.M., A.Z.S., T.B., and M.E.R. analyzed data; and A.M., A.Z.S., T.B., and M.E.R. wrote the paper.

Reviewers: G.B., NYU Neuroscience Institute; O.S., Indiana University; J.D.V., Weill Cornell Medical College; and R.Y., Columbia University.

The authors declare no conflict of interest.

See Commentary on page 5263.

This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1503960112/-/DCSupplemental.

## References

- ↵.
- Berger H

- ↵
- ↵
- ↵.
- Ikegaya Y, et al.

- ↵.
- Luczak A,
- Barthó P,
- Marguet SL,
- Buzsáki G,
- Harris KD

- ↵
- ↵
- ↵.
- Buzsáki G,
- Draguhn A

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Petersen CC,
- Hahn TT,
- Mehta M,
- Grinvald A,
- Sakmann B

- ↵
- ↵.
- Miller JE,
- Ayzenshtat I,
- Carrillo-Reid L,
- Yuste R

- ↵
- ↵.
- Buckner RL,
- Krienen FM,
- Castellanos A,
- Diaz JC,
- Yeo BT

- ↵.
- Choi EY,
- Yeo BT,
- Buckner RL

- ↵
- ↵.
- Yeo BT, et al.

- ↵
- ↵.
- Seeley WW, et al.

- ↵.
- Biswal BB, et al.

- ↵.
- Beckmann CF,
- DeLuca M,
- Devlin JT,
- Smith SM

- ↵.
- Mitra A,
- Snyder AZ,
- Hacker CD,
- Raichle ME

- ↵
- ↵
- ↵
- ↵
- ↵.
- Buckner R, et al.

- ↵.
- Minka TP

*Advances in Neural Information Processing Systems*(MIT Press, Cambridge, MA), Vol 13, pp 598–604 - ↵.
- van den Heuvel MP,
- Kahn RS,
- Goñi J,
- Sporns O

- ↵
- ↵.
- Esposito F, et al.

- ↵.
- Gordon EM, et al.

- ↵
- ↵
- ↵
- ↵.
- Sheroziya M,
- Timofeev I

- ↵.
- Sirota A,
- Csicsvari J,
- Buhl D,
- Buzsáki G

- ↵.
- Mohajerani MH,
- McVea DA,
- Fingas M,
- Murphy TH

- ↵.
- He BJ,
- Snyder AZ,
- Zempel JM,
- Smyth MD,
- Raichle ME

- ↵
- ↵
- ↵
- ↵.
- Fucke T, et al.

- ↵
- ↵.
- Plenz D,
- Kitai ST

- ↵.
- Amzica F,
- Steriade M

- ↵.
- Poskanzer KE,
- Yuste R

- ↵
- ↵.
- Fellin T, et al.

- ↵
- ↵
- ↵.
- Massimini M,
- Huber R,
- Ferrarelli F,
- Hill S,
- Tononi G

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Katz LC,
- Shatz CJ

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Brier MR, et al.

## Citation Manager Formats

## Article Classifications

- Biological Sciences
- Neuroscience

## See related content: