Developmental variations in recurrent spatiotemporal brain propagations from childhood to adulthood

Domain: biology
Published: 2026-01-10
Editorial outcome: minor_revision

Abstract
The brain undergoes profound structural and functional development from childhood to adolescence. Convergent evidence suggests that neurodevelopment proceeds in a hierarchical manner, characterized by heterogeneous maturation patterns across brain regions and networks. However, the maturation of the intrinsic spatiotemporal propagations of brain activity remains largely unexplored. This study aims to bridge this gap by delineating spatiotemporal propagations from childhood to early adulthood. By leveraging a recently developed approach that captures time-lag dynamic propagations, we characterized intrinsic dynamic propagations along three axes: sensory-association (S-A), ‘task-positive’ to default networks (TP-D), and somatomotor-visual (SM-V) networks, which progress towards adult-like brain dynamics from childhood to early adulthood. Importantly, we demonstrated that as participants mature, there is a prolonged occurrence of the S-A and TP-D propagation states, indicating that they spend more time in these states. Conversely, the prevalence of SM-V propagation states declines during development. Notably, top-down propagations along the S-A axis exhibited an age-dependent increase in occurrence, serving as a superior predictor of cognitive scores compared to bottom-up S-A propagation. These findings were replicated across two independent cohorts (Human Connectome Project: Development and Nathan Kline Institute-Rockland Sample), emphasizing the robustness and generalizability of these findings. Our results provide new insights into the developmental progression of functional dynamics during youth and their role in supporting cognition.

Introduction
The brain is a dynamic system with continuously evolving neural activity that interacts across networks. Within the dynamic landscape, temporal synchronization between brain networks exhibits rich time-varying features, involving rapid reconfigurations into transient whole-brain states and forming spatiotemporal patterns over time1–5. Accumulating evidence describes recurrent spatiotemporal patterns, detectable using a range of recording and imaging techniques in humans, nonhuman primates, and rodents1,6–8. These dynamic changes may index ongoing cognition and arousal, and reflect the brain’s continuous adaptation to internal and external environment – which is critical for maintaining cognitive flexibility9–11. Importantly, recent work has identified these dynamic fluctuations as a combination of instantaneous standing waves and time-lagged traveling waves propagating across functionally defined networks10. In adults, the topography of these waves prominently features spatiotemporal propagations along sensorimotor-association (S-A) and ‘task-positive’-’task-negative’ axes10. However, it remains largely unknown how these spatiotemporal patterns evolve and mature from childhood to early adulthood and contribute to cognitive development.

During childhood to adolescence, the brain undergoes profound structural and functional changes that shape cognition, emotion, and behavior12. Substantial evidence suggests that the brain development during this period is heterochronic across diverse brain regions. Early maturation primarily occurs in unimodal areas, including sensorimotor, visual and auditory cortices, succeeded by the higher-order association cortices that support advanced cognitive functions such as semantic processing, decision-making, problem-solving, and socio-emotional cognition13–18. Concurrently, functional brain organization also exhibits heterochronic changes during this period. In childhood, functional connectivity is predominantly organized among different unimodal regions, likely to facilitate sensory and motor experiences during the early stages of life. By adolescence, however, there is a notable shift towards an increased prominence of the functional organization along the S-A axis19,20. Despite these advances, studies investigating functional dynamics in youth have largely relied on instantaneous synchronization with zero-lag approaches, which may not capture the whole-brain traveling waves and their propagations across networks10. As a result, how the dynamic propagations develop and interplay with the functional organization during childhood to early adulthood remains poorly understood.

Here, we investigate age-related changes in brain dynamic propagations during childhood and early adulthood. To capture the time-varying nature of propagation across the developmental samples previously unattainable, we leveraged Complex Principal Component Analysis (CPCA), a recently developed method designed to decompose dynamics systems as a mixture of traveling and standing waves10. CPCA is capable of identifying zero-lag and time-lag inter-regional synchrony by detecting recurring spatiotemporal propagations that are representative of prior studies on temporal dynamics (e.g., quasi-periodic patterns)8,10,21,22. Prior work using this method identified oscillating waves in young adults spanning sensorimotor to association networks (S-A)10. Beyond this, they also captured an additional axis of propagation between task-positive and default networks (TP-D)10. Building on these findings, we characterized these propagations and investigated their maturation from childhood to early adulthood. We hypothesized that these spatiotemporal propagations develop progressively from childhood to early adulthood, maturing to resemble adult-like brain dynamics. We recognized that these time-resolved propagations share similar spatial patterns with the cortical gradients of static functional connectivity. To examine this relationship, we further compared the time participants spent in each propagation state with the explanatory power of corresponding functional gradients to reveal how dynamic propagations contribute to the development of static functional gradients during this period. Next, considering that the direction of propagations (e.g., top-down versus bottom-up) might vary between youth and adults23, we examined age-related variations in the prevalence of specific propagation direction during this period. We demonstrated that these propagations can predict individual differences in cognitive performance. Finally, we validated the robustness and generalizability of our findings, including the characterization of dynamic propagations and developmental trajectories, through replication in an independent dataset.

Results
Reliable recurring dynamic propagations were characterized from childhood to adulthood
We first characterized propagation components, which represent distinct spatiotemporal patterns that dynamically evolve brain-wide over time, using CPCA10. Unlike traditional dynamic techniques, which confines dynamic states within a fixed window size, CPCA flexibly tracks recurring dynamic cycles in the phase domain from 0 to 2π (see CPCA framework in Methods and Figure. S1). Whole-brain fMRI signals are decomposed into recurring spatiotemporal patterns, propagating across different brain regions from the starting point (phase 0 radians) to the midpoint (phase π radians), and then completes a full circle back to the initial starting phase (phase 2π or 0 radians). To examine developmental changes, we separately extracted the propagation components in developmental datasets (HCP-D, N=408, 215 females, age range: 8–21 years) and middle-aged adult reference dataset (HCP-A, N=399, 225 females, age range: 36–64 years). The HCP-A dataset was selected as a reference due to its identical acquisition parameters to HCP-D. In both datasets, the first three propagation components explained over 50% of the total variance of resting-state fMRI dynamics. To illustrate the propagations across the cortex at different cycle stages, we visualized the spatiotemporal patterns of each component over the phase 0 to 2π sampled at specific phase points (Fig. 1A). The first three group propagation components in the reference dataset revealed distinct spatiotemporal dynamics. The first propagation spans the cortex between the somatosensory and transmodal association networks (Pattern S-A). The second propagation captures a cycle between the task-positive networks (i.e., frontoparietal and attention networks) and the default mode network (Pattern TP-D). The third component propagates between somatomotor and visual networks (Pattern SM-V). Each of these three spatiotemporal patterns explained over 10% of the variance for both HCP-A and HCP-D (Fig. 1B). Of note, the raw order of these components was determined by the explained variance ratio (EVR), which varies across individuals and between youth and adult cohorts. To enable inter-individual comparability of spatiotemporal patterns for both within the HCP-D cohort and relative to the adult reference sample (HCP-A), we used orthogonal alignment to project the individual-specific components onto a shared dimensional space from adult reference. At the group level, each of the spatiotemporal patterns exhibited high similarity (r > 0.85) across youth and adult cohorts (Fig. 1C, diagonal). At the individual level, the spatiotemporal patterns were replicable across sessions within individuals, achieving high test-retest reliability with discriminability over 0.9 for Pattern S-A and TP-D, and over 0.85 for Pattern SM-V (Fig. 1D). The duration of the recurrent periodic cycle remained consistent within cohorts (HCP-D: 30.04±3.02 seconds, HCP-A: 30.68±3.23 seconds, Figure. S2).

Maturation of spatiotemporal patterns in youth towards adult-like dynamics
Having demonstrated reproducible and reliable spatiotemporal patterns in individuals in youths and adults, we examined their developmental effect by comparing individual-specific propagation patterns obtained from childhood to early adulthood cohort HCP-D with the middle-aged adult reference patterns (Fig. 1E). Overall, the propagation patterns are highly similar between youths and adults (Pattern S-A, mean r=0.93±0.02; Pattern TP-D: r=0.89±0.03; Pattern SM-V: r=0.90±0.04). Notably, the similarity increased significantly with age for all three patterns, controlling for sex and head motion (Pattern S-A: ρ=0.42; Pattern TP-D: ρ=0.32; Pattern SM-V: ρ=0.15, all p<0.05) from childhood to early adulthood, indicating a progressive refinement of functional dynamics toward an adult-like spatiotemporal pattern during this period. To further assess the extent to which these three dynamic propagations collectively explain age-related neurodevelopmental changes during this period, we employed a partial least squares regression (PLSR) model with five-fold cross-validation to estimate individual ages using three propagation patterns and achieved high accuracy (Fig. 1F, ρ=0.80 and Mean Absolute Error = 1.87 years).

Individuals spend more time in Patterns S-A and TP-D but less time in Pattern SM-V from childhood to early adulthood
Next, we examined the developmental effect on the time participants spent in each propagation state. We first defined the intrinsic moment-to-moment dominant state based on the strength of the propagation patterns. The dominant state was determined using a winner-take-all approach, selecting the pattern with the highest propagation strength (temporal magnitude) at each time point. We counted the time points of the dominant state and normalized it by the total number of time points in the scan, referring to this measure as the occurrence ratio (Fig. 2A). Time points where none of the first three patterns was dominant were categorized as the “other patterns” state. We then grouped participants into three developmental stages (children: 8–12 years, adolescents: 13–17 years, early adulthood: 18–21 years) and compared the occurrence ratio of each dynamic state across these developmental periods (Fig. 2C). This approach allows for identifying the stage-specific differences linked to developmental milestones and might capture stage-specific changes tied to key neurodevelopmental transitions. For all three propagations, the occurrence ratio showed significant group differences across developmental stages (F (3,803)=43.84, p<0.01 for Pattern S-A; F (3,803)=35.09, p<0.01 for Pattern TP-D; and F (3803)=35.88, p<0.01 for Pattern SM-V). The occurrence ratio of Pattern S-A significantly increased from childhood to adolescence (t 309=3.33, p<0.01, effect size = 0.38), from adolescence to early adulthood (t 267=2.81, p<0.01, effect size = 0.36), and from early adulthood to middle adulthood (t 494=3.04, p<0.01, effect size = 0.34). Pattern TP-D exhibited no significant change of state occurrence ratio from childhood to early adulthood but a significant increase between early adulthood to middle-aged adulthood (t 494=4.71, p<0.01, effect size = 0.53). In contrast, Pattern SM-V decreased its occurrence ratio starting from adolescence (adolescence to early adulthood: t 267=-3.77, p<0.01, effect size=-0.48, early adulthood to middle adulthood: t 494=-3.25, p<0.01, effect size=-0.37). We also replicated the above analysis using a restricted subsample of the HCP-A cohort (age < 50 years) as well as the HCP-Young Adult (HCP-YA, age 21-35 years, n=892) as an alternative reference. Similar results were observed (Figure. S3). Considering the brain state at a given time point could be a combination of multiple propagation patterns instead of a single ‘winner-takes-all’ pattern. In addition to the ‘winner-takes-all’ strategy, we also used a weighted approach that quantifies the relative contribution (i.e., proportion) of each dynamic pattern by comparing their respective magnitudes at a given point. We then averaged these proportions across all time points to estimate the mean occurrence contribution of each pattern for each individual (Figure. S17).

Spatiotemporal propagations reflect both dynamic and static functional organization
Notably, the spatial patterns of S-A, TP-D, and SM-V closely align with established functional gradients, which represent the low-dimensional axes of static functional connectivity across brain regions24. This spatial consistency between dynamic and static functional organization has been highlighted in previous studies in adults10,25. Here, we first examined their relationship in developmental cohorts (HCP-D) by comparing the occurrence time of dynamic states and the EVR of static functional gradients. Across all three dynamic states, their occurrence ratio showed strong correlations with the EVR of the corresponding static functional gradients (Fig. 3A, Pattern S-A: r=0.46, p<0.01, Pattern TP-D, r=0.44, p<0.01, Pattern SM-V, r=0.52, p<0.01).
Of note, recent findings demonstrated that static functional gradients undergo significant developmental changes from childhood to adolescence19,20. We hypothesize that these developmental changes are supported by the corresponding dynamic states. Specifically, as participants spend more time in a particular dynamic state, that state becomes the dominant gradient in static functional organization. To test this, we examined the age effects on EVR of functional gradients before and after regressing out corresponding dynamic patterns from the fMRI signals (Fig. 3B, C). Notably, regressing out dynamic patterns significantly reduced the age effect on corresponding functional gradients, in particular for the functional gradient and pattern along the S-A axis (ρ=0.30 to ρ=0.14) (Fig. 3C). In addition, this attenuation was also observed, though to a lesser extent, in non-corresponding functional gradients as well (Figure. S15). These findings provide compelling evidence that the time participants spend in dynamic states plays a critical role in the developmental progression in static functional organization from childhood to early adulthood.

Increased time spent in top-down propagation and decreased time spent in bottom-up propagation from childhood to early adulthood
The above analyses revealed that the Pattern S-A became more prevalent with age. However, the analysis did not distinguish directionality of propagations that ascend the S-A axis (i.e., bottom-up) and those that descend the S-A axis (i.e., top-down). To address this, we defined ‘bottom-up’ propagation as the temporal changes from sensory to association regions during the 0 to π phase and ‘top-down’ propagation as the reverse from association to sensory areas during the π to 2π phase (Fig. 4A). We then calculated the time participant spent within each propagation phase bin (π/16) and assessed its relationship with age. Notably, this analysis revealed an age-related increased prevalence of top-down propagation (i.e., π to 2π), along with a decrease in time spent in bottom-up propagation, particularly during the later propagation phase (i.e., π/2 to π) of the bottom-up process as it ascends to the default mode network (Fig. 4B). Additionally, we investigated whether dynamic propagations could predict the total cognitive composite score, a comprehensive measure of crystallized and fluid intelligence26. Using a PLSR model with 5-fold validation with 1000 iterations, we examined the predictive performance (i.e., r-value between actual and predicted scores) for each of four propagation phase stages. Our results showed that both bottom-up and top-down propagation achieved significant prediction performance (r 0-π/2=0.17±0.04, p _permutation=0.053; r π/2-π=0.19±0.03, p _permutation=0.013; r π/2-3π/2=0.20±0.03, p _permutation<0.01; r 3π/2-2π=0.24±0.03, p _permutation<0.01). Given that top-down processing involves goal-directed behavior and selective attention critical for complex cognitive tasks, it appears the top-down propagation (π to 2π stage) predicts cognitive score better than the bottom-up propagations.

Developmental shifts in propagation directionality and cognitive relevance of Patterns TP-D and SM-V
Similar to the analysis for Pattern S-A, we also examined age-related changes in propagation directionality for Patterns TP-D and SM-V. For Pattern TP-D, we observed an increasing occurrence of external-to-internal propagation from attention network to default network (π/2 to 3π/2) and a decreasing occurrence of the reverse internal-to-external propagation from default network to attention network (0 to π and 3π/2 to 2) from childhood to early adulthood (Figure. S4). For Pattern SM-V, the occurrence increased at the propagation phase from the insular and motor regions to visual networks, while it decreased during the phase from the visual to motor networks. When predicting cognitive scores, both Pattern TP-D and SM-V showed weaker predictive performance compared to Pattern S-A (Figure. S4, Pattern TP-D: r 0-π/2 =0.21±0.03, rπ/2-π=0.14±0.03, rπ/2-3π/2 =0.18±0.03, r 3π/2-2π=0.15±0.03; Pattern SM-V: (r 0-π/2 =0.15±0.03, rπ/2-π=0.12±0.04, rπ/2-3π/2 =0.17±0.03, r 3π/2-2π=0.12±0.04).

Validation of dynamic propagations and their developmental effects across independent cohorts
The above analyses revealed developmental effects in dynamic propagations from childhood to early adulthood in HCP-D cohort. To evaluate the generalizability of our findings, we replicated our analysis in an independent cohort from the NKI-RS datasets and selected the participants including children, adolescents, and early adults (N=269, age=6-21). This dataset consists of three runs with distinct acquisition parameters (TR=645ms, 1400ms, 2500ms), with the main figures (Fig. 5) presenting results averaged across acquisition settings for differences in TR. The propagation patterns and developmental effects observed in NKI-RS were consistent with those identified in HCP-D. Specifically, the first three propagation patterns explained over 50% variance of resting state fMRI dynamics, each contributing more than 10% of EVR (Figure. S5B). Across cohorts, the propagation patterns demonstrated high similarity compared to the adult reference (HCP-A) (r=0.88±0.04 for Pattern S-A, r=0.81±0.06 for Pattern TP-D, and r=0.87±0.03 for Pattern SM-V). More importantly, these patterns showed similar significant age-related resemblance to adult reference patterns, with Pattern S-A showing the strongest correlation (ρ=0.46, p<0.05) (Fig. 5A, Figure. S5D provides results for each acquisition). Across developmental stages from childhood, adolescent, early adulthood, to middle-age adulthood in NKI-RS, we observed consistent age-related shifts in occurrence ratios in NKI-RS. Specifically, older individuals exhibited a higher proportion of time spent in Patterns S-A and TP-D with a reduced presence in Pattern SM-V (Fig. 5C and Figure. S6 for each acquisition). The correlation between occurrence ratios of dynamic patterns and EVR of static functional gradients also showed high similarities (Figure. S7 for each acquisition).
Consistent with HCP-D, a decreased occurrence time was observed in the later phase of bottom-up propagation (π/2 to π). The increasing trend in the early phase (0 to π/2) of bottom-up propagation from HCP-D was more pronounced in NKI-RS, while no significant increase in top-down propagation was detected (Fig. 5D and Figure. S8 for each acquisition). This discrepancy might be attributed to the limited data available per individual in the NKI-RS compared to HCP-D. Although we merged occurrence time across three runs to increase data quantity per individual, the NKI-RS dataset have less data per individuals and were collected using varying scan acquisitions, resulting in greater heterogeneity compared to the more standardized HCP-D cohorts. Nevertheless, the overall robustness of the dynamic propagations and their developmental effect was evident across cohorts.

Sensitivity analyses using an alternative parcellation and low-motion subset
To further validate the robustness of our findings to methodological variations, we conducted two sensitivity analyses within the HCP-D cohort: (1) using a different parcellation (Schaefer 400 parcels)27, (2) using a low-motion subset consisting of individuals without substantial head movement (N=287, Mean Framewise Displacement <0.2mm). Consistent findings were observed across both analyses. Specifically, from childhood to early adulthood, individuals spent more time (i.e., occurrence ratio) for Pattern S-A and Pattern TP-D while less time for pattern SM-V (Figure. S10). The occurrence ratio of the propagation states also showed a strong correlation with the EVR of static functional gradients (Figure. S11). The developmental effect on occurrence ratio at different propagation phases were also confirmed across all three propagation patterns (Figure. S12). From childhood to adulthood, participants exhibited a decrease in time spent in bottom-up propagation and an increase in time spent in top-down propagation along the S-A axis. Notably, Pattern S-A showed better predictive performance for the total composite cognitive scores compared to Pattern TP-D and Pattern SM-V (Figure. S13), with the highest prediction performance using the top-down propagation in Pattern S-A (Figure. S12). Collectively, these replications highlight the robustness of our results associated with head motion and scale of cortical parcellation.

Methods
Datasets
We utilized fMRI data from three large-scale datasets in this study: the Human Connectome Project Development (HCP-D) and Aging (HCP-A), and the Nathan Kline Institute-Rockland Sample (NKI-RS). Our primary analyses were conducted using the HCP-D dataset for the youth sample, with the HCP-A dataset serving as the adult reference. These analyses were then replicated in the NKI-RS dataset. Details of demographic information are summarized in Table S1.
HCP datasets are publicly available neuroimaging data from the Lifespan 2.0 Release of HCP project64,65. Resting-state fMRI data were collected across four sites (Harvard, UCLA, University of Minnesota, and Washington University in St. Louis), using 3 Tesla Siemens Prisma scanners. The acquisition protocol is identical for both HCP-D and HCP-A. Specifically, a 2D multi-band gradient-recalled echo (GRE) echo-planar imaging (EPI) sequence was employed (TR=800ms, TE=37ms, flip angle=52 degree, with 2mm isotropic voxel resolution). With this acquisition, during two sessions, two T2*-weighted image runs were acquired, each lasting 6.5 minutes and using different phase encoding directions (AP and PA), resulting in a total of 26 minutes of resting-state fMRI data for each subject. We included the participants whose data had been minimally preprocessed and passed the HCP quality control assessment66. Our analyses included 408 participants (8–21 years) from HCP-D and a group-average of 399 participants (22-64 years) from HCP-A as adult reference. Additionally, we included HCP Young Adult cohort (N=892, age 21-35, TR=720ms, TE=33.1ms, flip angle=52 degree, with 2mm isotropic voxel resolution) to further validate the replication of the spatiotemporal patterns observed in the HCP-A adult reference group.
The NKI-RS dataset67 is part of a community-ascertained lifespan dataset. Resting-state fMRI data acquired in three acquisitions, each with distinct temporal and spatial resolutions: a 9.7 minute scan with TR=645ms at 3-mm isotropic resolution (900 time points, TE=30ms, flip angle=60 degree), a 9-minute scan with TR=1400ms at 2-mm isotropic resolution (404 time points, TE=30ms, flip angle=65 degree), and an 5min scan with TR=2500ms at 3-mm isotropic resolution (120 time points, TE=80ms, flip angle=80 degree). The details of the acquisition can be found at (https://rocklandsample.org/for-researchers-rockland-sample-i-2011-2022-neuroimaging-protocol). After the preprocessing and quality control, individuals aged 6–21 years (N=298 for TR=645ms; N=288 for TR=1400ms; N=282 for TR=2500ms) were included as the youth sample. We also included participants with age 22–64 (N=526 for TR=645ms, N=533 for TR=1400ms, N=533 for TR=2500ms) as an early adult group in replication analysis. A total of 269 young participants (age 6–21 years) and 487 middle-age adults (age 22–64 years) with data from all three runs were included for replication in Fig. 5.

MRI data preprocessing
HCP-D and HCP-A datasets. Resting-state fMRI data were preprocessed through the HCP pipeline66. Briefly, the minimal preprocessing of fMRI data included distortion correction, head movement correction, bias field correction, intensity normalization, registration to a mixture of volume and surface “grayordinate” CIFTI template using multimodal surface mapping (MSMAll), high-pass (>0.009 Hz) filtering, ICA-FIX denoising, and spatially smoothed with 2-mm FWHM Gaussian kernel. In our study, we further applied global signal regression and temporal filtering to the HCP minimally preprocessed data using the conventional resting-state fMRI low-frequency band-pass range of 0.01-0.08 Hz. Finally, the time series were averaged within each parcel using Schaefer 17-network 200 parcellations27. To assess the robustness of the findings, Schaefer 17-network 400 parcellation was also employed in our replication analysis27. We further replicated the main results using a subset of individuals (N=287 with low head motion (mean framewise displacement <0.2mm)68.
NKI-RS. Both replication datasets were preprocessed using fMRIPrep69 and XCP-D pipelines70. Similar to the HCP pipeline, the functional MRI data pipeline of fMRIPrep and XCP-D includes head motion correction, registration to the “grayordinate” CIFTI template. Nuisance regressors include 36 parameters (six head motion parameters, white matter, cerebrospinal fluid, and global signal, along with their squired terms, derivatives, and squared derivatives)71. The data were temporally filtered with the same low-frequency range of 0.01-0.08 Hz and parcellated using Schaefer 17-network 200 parcellation27.

Dynamic spatiotemporal pattern analysis
To extract the propagation patterns of resting state fMRI signals, we adopted the CPCA approach10,72,73. This method converted the fMRI data into analytic signals followed by a dimensionality reduction technique to decompose the high-dimensional analytic data into low-dimensional principal components. Specifically, for each parcel, we first computed the analytic signal xa(t) using Hilbert transform as follows:
xa(t) = x(t) + i * H[x(t)] (1)
The real part x(t) of the analytic signal represent the original fMRI time series, while the imaginary part corresponds to the Hilbert-transformed signal, which is phase-shifted signals by π/2 radians (Figure. S1A). Next, we applied fast randomized singular value decomposition (SVD: https://github.com/facebookarchive/fbpca) on the complex-valued data matrix Xa (size: the number of time points by the number of parcels), decomposing into ten components.
Xa = UΣVH (2)
The left singular vectors U contain the orthogonal temporal vectors (size: the number of time points by the number of components), Σ is a diagonal matrix which holds the singular values, and the right singular vectors V consist of the orthogonal spatial vectors (size: the number of components by the number of parcels) with H as Hermitian transpose. These vectors are complex-valued signals which can be represented in terms of amplitude and phase components using Euler’s formula. The amplitude quantifies the strength of the components, and the phase specifies the relative delay, expressed in radians from 0 to 2π. To capture a full cycle of cortical propagations, we extracted the spatiotemporal patterns across 32 equally spaced phase segments (i.e., 2π/32), yielding 32 spatial maps at each phase segment from 0 to 2π for each component. The selection of 32 phase bins achieves a balance between computational efficiency and temporal resolution in the phase domain. Specifically, we first averaged the temporal vector U for each component within each of 32 phase segments to generate the resampled temporal vector (32 by the number of components). The spatiotemporal patterns were then computed by performing a component-wise outer product between the resampled temporal vector and the corresponding spatial vector. In Fig. 1A, we demonstrated and visualized spatial maps at phase 0, 5π/32, 10π/32, 15π/32, 20π/32, 25π/32, and 30π/32 for the first three components.
To improve inter-individual comparability, we first generated the adult components from the HCP-A cohort to serve as a reference for alignment. We concatenated fMRI data across individuals in HCP-A cohort and performed the SVD analysis on the analytic fMRI signal to extract the temporal and spatial vectors as adult reference. We then used the Procrustes matching to align spatial vectors (size: the number of parcels by the number of components) between each individual from the HCP-D cohort and the reference74,75. The alignment transformation (rotation and reflection) was subsequently applied to the temporal vectors to match with the aligned spatial vectors for each individual. This rotation can reorder the component indices to match the reference while preserving the original spatiotemporal characteristics of each pattern. The resulting reordered component structure in HCP-D individuals is visualized in Figure. S16. Adult resemblance was quantified using a Pearson’s correlation between each individual’s flattened propagation pattern and the group-reference adult pattern.

Occurrence ratio of dynamic state and age effect analysis
To assess the time participants spent (i.e., occurrence ratio) in each of the propagation states defined by the CPCA components, each time point of the cortical fMRI frames was assigned to one of the CPCA components. The temporal strength of the components at every time point was quantified using weighted temporal vector (i.e., UΣ, size: the number of time points by the number of components). The moment-by-moment propagation state was determined using a ‘winner-takes-all’ strategy, where the state with the highest temporal magnitude (Fig. 2A) at each time point was assigned as the dominant state. This simple assignment approach allows for a straightforward temporal interpretation of the dominant state over time. Here, temporal magnitude was calculated as the amplitude of the temporal vector weighted by the singular value of the corresponding propagation pattern at that moment. For the time points where none of the first three components were dominant, we grouped them as “other patterns” state. For each individual, the occurrence ratio was calculated as the proportion of time points during which a specific state was dominant, relative to the total number of time points of fMRI scans. We also replicated the above analysis using a restricted subsample of the HCP-A cohort (age < 50 years) as well as the HCP-Young Adult (HCP-YA, age 21–35 years, n=892) as an alternative reference. Similar results were observed (Figure. S3). Considering the brain state at a given time point could be a combination of multiple propagation patterns instead of a single ‘winner-takes-all’ pattern. In addition to the ‘winner-takes-all’ strategy, we also used a weighted approach that quantifies the relative contribution (i.e., proportion) of each dynamic pattern by comparing their respective magnitudes at a given point. We then averaged these proportions across all time points to estimate the mean occurrence contribution of each pattern for each individual (Figure. S17).
To compare the occurrence ratio across the developmental stages from childhood to adulthood, participants from the HCP-D were first grouped into three developmental stages within dataset children (8 ≤ age < 13), adolescents (13 ≤ age < 18), and early adulthood (18 ≤ age < 22). The late adulthood stage was measured using HCP-A dataset. ANOVA was used to compare the difference across four age groups, followed by sequential pairwise comparisons between age stages using t-tests. To account for multiple comparisons, we applied the Benjamini-Hochberg procedure to control the false discovery rate in the Python package ‘statsmodels’76.

Dynamic directionality analysis
To assess the age effect on the directionality of dynamic propagations along each CPCA component, we counted the time (i.e., occurrence ratio) that participants spent in each of the 32 phase bins (2π/32). Next, we calculated the Spearman’s correlation between age and occurrence ratio for each phase bin. For Pattern S-A, propagations during the 0 to π phase captured spatial transitions from visual and somatomotor to default mode regions (i.e., bottom-up), while the π to 2π phase captured a reversed, top-down propagation. We applied the same approach to assess the age effect on occurrence ratios within each propagation phase bin for Patterns TP-D and SM-V.

Brain-wide prediction
In this study, the brain-wide prediction was conducted using partial least squares regression (PLSR) with stratified five-fold cross-validation to ensure a balanced distribution of the dependent variable77,78. In PLSR, we z-scored the input features for normalization and set the hyperparameter n_components=3 in the ‘scikit-learn’ Python package79. To examine the extent to which age can be represented by all three dynamic patterns in youth, we flattened spatiotemporal patterns of three components (32 spatial maps by 3 components) and concatenated them as input features. Model performance was evaluated by calculating the mean square error and Spearman’s correlation between actual and predicted age in the HCP-D dataset.
A similar approach was used for brain-behavior prediction of cognitive scores. To identify which propagation stage (phase segment) best predicts cognitive performance, we performed separate predictions for each component and each propagation stages (0 to π/2, π/2 to π, π to 3π/2, and 3π/2 to 2π), also using PLSR with five-fold cross-validation. The prediction was repeated 1000 times to ensure robustness of the prediction score. To assess statistical significance, a null distribution was generated using a permutation test (1000 iterations) by randomly shuffling cognitive scores. The mean observed prediction performance was then compared against this null distribution to determine statistical significance.

Functional gradients of functional connectivity
To obtain the low-dimensional representation of the static functional connectome, we calculated the first ten functional gradients of the functional connectivity for each individual. Previous studies have demonstrated that gradients derived from the principal component analysis (PCA) provided higher reliability and prediction power than the traditional diffusion embedding method57. Accordingly, we used PCA to extract the functional gradients here. Specifically, we computed the adjacency matrix by applying Pearson’s correlation to the row-wise thresholded (top 10%) functional connectivity matrix using the BrainSpace toolbox80. Procrustes alignment (rotation and reflection) was then applied to match the individual gradients with the reference gradient derived from the averaged functional connectivity matrix of the HCP-A cohort. Finally, we compared the EVR of the functional gradients with that of the propagation components using Spearman’s correlation.

Reliability analysis
To examine the test-retest reliability of the spatiotemporal patterns, we employed a multivariate non-parametric method, discriminability81. For the HCP-D cohort, each subject had four rs-fMRI scans with two different phase encoding directions, we first concatenated two scans with the same phase encoding direction, and then evaluated the reliability between two resulting concatenated subsets. For the NKI-RS cohorts, the test-retest reliability analysis was assessed by comparing pairs of group-level spatiotemporal patterns from rs-fMRI scans acquired at TR=645ms, TR=1400ms, and TR=2500ms.

Replication analysis in NKI-RS
To evaluate the generalizability of our findings, we replicated the above analyses, including dynamic spatiotemporal pattern analysis, occurrence ratio, and directionality analysis using the NKI-RS dataset. The dataset was collected across three runs during the same session, each with different acquisition parameters (TR=645ms, 3mm isotropic, 9.7 min; TR=1400ms, 2mm isotropic, 9.4 min; TR=2500ms, 3mm isotropic, 5.0 min). We calculated the spatiotemporal patterns for each run separately and replicated the findings in Figure S4-7. To increase data quantity robustness estimation per individual, we also merged three runs together by averaging spatiotemporal patterns and occurrence ratio weighted by the scan duration per run for each individual.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Byeon and colleagues investigated the maturation of spatiotemporal propagations from childhood to early adulthood. They\nused the state-of-the-art complex PCA method (Ref 10) to capture the time-lag dynamic propagations of functional brain\nnetworks. The authors showed the maturation of three propagation patterns: sensory-association (S-A), ‘task-positive’ to\ndefault networks (TP-D), and somatomotor-visual (SM-V) networks, and their predictive power for cognitive scores. The\nresults were replicated in two independent cohorts, highlighting the robustness and generalizability of these findings.\nI really enjoyed reading this paper. The propagation model is complicated, but the authors provided a clear description and a\nrigorous analysis. The sensitivity analysis was also well conducted with different parcellations and the low-motion\nsubdataset. The two independent datasets are important for replication analysis. These results provide insights into how\nnetwork dynamics in the human brain mature from childhood to early adulthood and how this maturation process supports\ncognitive development. Here, I only have a few minor comments that I hope are useful in further clarifying the conception,\nmethods, and results.\nIntroduction. The authors pose the scientific question of how the dynamic propagation emerges and shapes the functional\norganization during childhood and early adulthood. If the authors wanted to address the emergence of network dynamics,\nthey might need to include data from very early life, such as infancy. Static functional connectivity studies have demonstrated\na substantial reorganization of functional brain networks from infancy to childhood (Gilmore et al NRN 2018; Cao et al TINS\n2017), in particular the default-mode network and frontoparietal networks that are closely related to the propagation of the S-\nA pattern. I am not urging the authors to add infant data, but they may need to rephrase the relevant descriptions throughout\nthe manuscript. It would be also encouraging to add a small discussion about the emergence, development, and\ndegeneration of the patterns.\nAnother related point is that the authors hypothesize that dynamic propagation shapes the functional organization of brain\nnetworks, such as functional gradients. This is an interesting but strong claim. Although they have provided results from a\nregression analysis in Figure 3, I would argue that the functional organization such as gradient architecture is also likely to\nprovide a backbone to shape the network dynamics.\nThe authors used a middle-aged adult dataset as a reference (HCP-A, age range:36-64 years). Why did the authors not use\nthe high-quality HCP dataset with a narrower age range as a reference? The age range in HCP-A is a bit wide and there\nmay very possibly be significant age-related decline (e.g., 50-64 years).\nThe authors showed the explained variance ratio of the three patterns in Figure 1 and indicated that these patterns appeared\nin childhood, even at 8 years of age. I am wondering how the explained variance ratios of these three patterns change with\nincreasing age. Is it possible that there is a reversal in the order of these patterns? For example, is the S-A pattern a\ndominant or secondary pattern during childhood? Ref 20 has shown a transition of the S-A gradient from secondary to\ndominant during late childhood. It would be interesting to see the developmental timing of these patterns.\nFigure 2. The authors showed the developmental effects in subgroups. Why did they not show continuous, age-related\nchanges as in Figures 1 and 3?", "reviewer_2": "This work investigates how dynamic propagation patterns of brain activity evolve across development and adulthood.\nUtilizing large-scale fMRI datasets (HCP and NKI), the authors apply a novel complex PCA (cPCA) framework to extract\nspatiotemporal patterns of brain activity. They identified three core propagation patterns—Sensorimotor-Visual (SM-V),\nTask-Positive–Default (TP-D), and Sensorimotor–Association (S-A)—and reveal robust, age-related changes in their\noccurrence and directionality. These changes track well with cognitive development (top down more than bottom up\npropagation), supporting the hypothesis that brain dynamics become more structured and top-down driven over the lifespan.\nThe paper is methodologically innovative, with thorough replication and sensitivity analyses.\nThere are a few issues to be addressed:\n1. There are significant changes in overall hierarchical structure derived from static brain FC from children to adolescents to\nadults (see Dong et al, Nguyen et al), which seem to also occur in the propagations identified by the authors. Procrustes\nalignment, while necessary for cross-group comparison can distort the results. It’d be good to show the raw unaligned\npatterns and their explained variance to show that Procrustes alignment only affects the sign/scale of the data value and\nnothing else.\n2. The dominant state was determined using a winner-take-all approach, which resulted in straightforward analysis and\nfindings, but might miss out the subtlety of brain activity at a certain time point, which is a combination of different\npropagation components identified. Sometimes a state is probably highly dominant, sometimes the difference between them\nmight not be so clear. Given that you have the magnitude of each state at a specific time point, a more nuanced method\nwould be good to derive a more accurate picture. Importantly, how do these relate to age and cognition?\n3. The interpretation of top-down vs. bottom-up propagation is intriguing but may be oversimplified. Since the directionality is\ninferred from relative phase differences, can the authors expand on how confident we can be in attributing cognitive control\nprocesses to such phase shifts?\n4. Do the propagation patterns seen here similar to what have been identified using a dynamic gradient approach (e.g.\nBrown et al, Neuroimage 2022)? It would be good to discuss the implications.\n5. The discussion mainly focuses on previous findings on static FC and gradient. I would be interested to see at least a point\non how the CPCA-derived brain states can be interpreted specifically in comparison to the brain state findings using\ntraditional dynamic methods like HMM.", "reviewer_3": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts."}, "author_response": "We are grateful for your detailed and constructive feedback and suggestions. In response,\nwe have added additional analyses and revised the manuscript. We address each of the\ncomments in details below. The point-by-point responses are highlighted in blue and new\ntext added in the revised manuscript is shown in green.\nReviewer #1 (Remarks to the Author):\nByeon and colleagues investigated the maturation of spatiotemporal propagations from\nchildhood to early adulthood. They used the state-of-the-art complex PCA method (Ref\n10) to capture the time-lag dynamic propagations of functional brain networks. The\nauthors showed the maturation of three propagation patterns: sensory-association (S-A),\n‘task-positive’ to default networks (TP-D), and somatomotor-visual (SM-V) networks, and\ntheir predictive power for cognitive scores. The results were replicated in two independent\ncohorts, highlighting the robustness and generalizability of these findings.\nI really enjoyed reading this paper. The propagation model is complicated, but the authors\nprovided a clear description and a rigorous analysis. The sensitivity analysis was also\nwell conducted with different parcellations and the low-motion subdataset. The two\nindependent datasets are important for replication analysis. These results provide insights\ninto how network dynamics in the human brain mature from childhood to early adulthood\nand how this maturation process supports cognitive development. Here, I only have a few\nminor comments that I hope are useful in further clarifying the conception, methods, and\nresults.\nIntroduction. The authors pose the scientific question of how the dynamic propagation\nemerges and shapes the functional organization during childhood and early adulthood. If\nthe authors wanted to address the emergence of network dynamics, they might need to\ninclude data from very early life, such as infancy. Static functional connectivity studies\nhave demonstrated a substantial reorganization of functional brain networks from infancy\nto childhood (Gilmore et al NRN 2018; Cao et al TINS 2017), in particular the default-\nmode network and frontoparietal networks that are closely related to the propagation of\nthe S-A pattern. I am not urging the authors to add infant data, but they may need to\nrephrase the relevant descriptions throughout the manuscript. It would be also\nencouraging to add a small discussion about the emergence, development, and\ndegeneration of the patterns.\nResponse: We thank the reviewer for the positive comment and for highlighting important\nreferences. We agree that the emergence of network dynamics begins in infancy, with\nsubstantial reorganization occurring during early development. From childhood to early\nadulthood, development involved continued refinement and maturation of such dynamics\nalong the S-A axis. Rather than representing their initial emergence, this period reflects\nthe strengthening of the S-A axis as a dominant mode of functional organization. In the\nrevision, we have rephrased relevant descriptions throughout the manuscript to more\naccurately reflect this developmental process (e.g., by replacing the word \"emerge\" with\n\"mature\" or \"develop\" in relevant sections). As suggested, we have also added a\ndiscussion paragraph as follows:\n“Notably, in the neonate-infancy stage, the spatial configuration of the hierarchical axis\ndoes not strictly follow the canonical somatomotor-to-association, but instead aligns more\nclosely with a posterior-anterior axis37. This axis undergoes gradual spatial reorganization\nthroughout development: the anterior end progressively expands to encompass the\nlateral temporal cortex and posterior cingulate cortex, both key regions of the default\nmode network. The posterior end increasingly strengthens the visual and somatomotor\ncortex, forming the primary somatosensory anchor of the hierarchical S-A axis. By\nchildhood, this hierarchical axis more closely resembles the canonical S-A axis, and by\nadolescence, it emerges as the dominant organizing principle of functional connectivity19.\nIn contrast, the SM-V axis appears to be well-established in infancy but gradually declines\nin prominence, particularly from childhood through adolescence38,39. Interestingly, our\nfindings revealed that the relationship between dynamic propagations and static\nfunctional brain organization during development may extend beyond a simple one-to-\none correspondence (Fig. S15).”\nAnother related point is that the authors hypothesize that dynamic propagation shapes\nthe functional organization of brain networks, such as functional gradients. This is an\ninteresting but strong claim. Although they have provided results from a regression\nanalysis in Figure 3, I would argue that functional organization such as gradient\narchitecture is also likely to provide a backbone to shape the network dynamics.\nResponse: We appreciate the reviewer’s perspective and agree that the hypothesis that\ndynamic propagation shapes functional organization might have been stated too strongly.\nTheir relationship might be more complex and bidirectional. As the reviewer pointed out,\nthe gradient architecture may provide a backbone that constrains how activity propagates\nacross the brain. The moment-to-moment activity may also vary and be influenced by\ninternal (e.g., respiration, vigilance, etc.) and external factors (e.g., learning, stimuli), and\nover time, may contribute to the refinement of static functional connectivity and gradient\narchitecture.\nIn response to this comment, we have revised the relevant text to be more accurate as\nfollows.\nThe Results section: “These findings provide compelling evidence that the time\nparticipants spend in dynamic states plays a critical role in the developmental progression\nin static functional organization from childhood to early adulthood.”\nThe Discussion section: “From one perspective, functional organization axes may serve\nas a backbone that constrains how activity propagates across the brain. At the same time,\nmoment-to-moment neural activity can vary and be modulated by internal (e.g. autonomic\nsignal) and external factors (e.g. learning)34,40. Over time, these dynamic fluctuations may\ncontribute to the refinement of static functional connectivity architecture (i.e. functional\ngradient). Thus, moment-to-moment neural dynamics may play a cortical role in a\nsubstantial portion of age-related differences in static macroscale brain networks during\ndevelopment.”\nThe authors used a middle-aged adult dataset as a reference (HCP-A, age range:36-64\nyears). Why did the authors not use the high-quality HCP dataset with a narrower age\nrange as a reference? The age range in HCP-A is a bit wide and there may very possibly\nbe significant age-related decline (e.g., 50-64 years).\nResponse: We appreciate the reviewer highlighting the potential concern regarding the\nwide age range and possible age-related decline within the HCP-A dataset. To directly\naddress this point, we have replicated our analysis using the restricted subsample of the\nHCP-A cohort (age < 50 years). The updated results are presented in the revised\nSupplementary Figure 3C (lower panel shown below). Results for participants older than\n50 years were also presented in this figure. As the reviewer pointed out, we observed a\nsignificant age-related decline in the occurrence ratio of both Pattern S-A and Pattern\nSM-V in the older middle-age group (age ≥ 50 years), compared to the younger middle-\nage group (age < 50 years), while no significant difference was observed for Pattern TP-\nD. We also replicated the analysis using the HCP-Young Adult (HCP-YA, age 21-35\nyears, n=892) as an alternative reference. The results from this analysis are included in\nSupplementary Figure S3 (upper panel shown below).\nNotably, although changes from younger (age < 50 years) to older (age ≥ 50 years)\nmiddle-aged adult subsamples were observed, the age-related changes from early\nadulthood to middle adulthood are consistent and reproducible using both HCP-A (ages\n36-50) and HCP-YA (ages 21-35) as middle-aged adult groups.\nFigure S3. Propagation patterns and state occurrence across Human Connectome Project Cohorts. (C) State occurrence\nratios were analyzed across various age stages. The upper panel depicts these results, with the HCP-YA cohort serving as an adult\nreference. The lower panel displays findings for the two distinct HCP-A sub-cohorts (age < 50 and age >= 50), where the HCP-A\n(age < 50) group is utilized as an adult reference cohort.\nAdditionally, Reviewer 2 raised a concern that the ‘winner-take-all’ approach capturing\nthe occurrence of dynamic patterns might miss out the subtlety of brain activity, since the\nbrain state at a given time point could be a combination of multiple propagation\ncomponents instead of a single ‘winner’ pattern. Therefore, we replicated our analysis\nusing the proportion approach, in addition to the ‘winner-take-all’ strategy. Specifically,\nwe compared the strength (i.e. magnitude) of the dynamic patterns to calculate their\nproportion at a given time point. The results (revised Supplementary Figure S17B, shown\nbelow as well) revealed similar age-related changes as the main result, Figure 2C: an\nincreased occurrence in magnitude of both Pattern S-A and Pattern TP-D across all age\ngroups from childhood to middle-aged adulthood, and a decrease in Pattern SM-V from\nadolescence to early adulthood. In the revision, all these updated results were now\nincluded to provide a comprehensive perspective.\nWe have added the relevant text in the Method section as follows:\n“We also replicated the above analysis using a restricted subsample of the HCP-A cohort\n(age < 50 years) as well as the HCP-Young Adult (HCP-YA, age 21-35 years, n=892) as\nan alternative reference. Similar results were observed (Fig. S3). Considering the brain\nstate at a given time point could be a combination of multiple propagation patterns instead\nof a single ‘winner-takes-all’ pattern, we also used a weighted approach that quantifies\nthe relative contribution (i.e. proportion) of each dynamic pattern by comparing their\nrespective magnitudes at a given point. We then averaged these proportions across all\ntime points to estimate the mean occurrence contribution of each pattern for each\nindividual (Fig. S17).”\nFigure S17. Occurrence contribution of dynamic patterns across age groups. (B) Mean proportion of pattern magnitude across\ndevelopmental stages: children, adolescents, early adulthood in HCP-D cohorts, and adult reference (HCP-A, age < 50). Age effect\nbetween consecutive development stages was tested using t-tests with significance denoted by asterisks (*: p < 0.05).\nThe authors showed the explained variance ratio of the three patterns in Figure 1 and\nindicated that these patterns appeared in childhood, even at 8 years of age. I am\nwondering how the explained variance ratios of these three patterns change with\nincreasing age. Is it possible that there is a reversal in the order of these patterns? For\nexample, is the S-A pattern a dominant or secondary pattern during childhood? Ref 20\nhas shown a transition of the S-A gradient from secondary to dominant during late\nchildhood. It would be interesting to see the developmental timing of these patterns.\nResponse: We appreciate the reviewer’s comments. As suggested, we conducted a\nsimilar analysis following the approach in Xia et al., 2022. Specifically, we plotted the\nexplained variance ratio (EVR) of each pattern against age and calculated its age effect\n(revised Supplementary Figure 14A, shown below as well). As the reviewer noted, the\nEVR from Pattern S-A is higher than both Pattern TP-D and Pattern SM-V even at age 8.\nTo directly assess whether Pattern S-A remains the dominant pattern during childhood,\nwe calculated the within-subject ratio of EVR between Pattern S-A and Pattern TP-D.\nAs shown in revised Supplementary Figure S14B (also included below, left panel), the ratio\nis greater than 1 (blue dashed line) in the vast majority of participants from age 8,\nindicating that Pattern S-A remains the dominant pattern from childhood through\nadulthood. We performed a similar ratio analysis for Pattern S-A versus Pattern SM-V,\nwith consistent results: Pattern S-A explains the most variance across all ages, and the\nratio increases significantly during childhood (figure below, middle panel). When\ncomparing the EVR between Pattern TP-D and Pattern SM-V, Pattern TP-D is generally\nrecognized as the secondary component, although in some individuals it may fall below\nSM-V (figure below, right panel). From childhood to adulthood, Pattern TP-D shows a\nsignificant age-related increase in its likelihood of being the second dominant component.\nWe have now included these results in the revision (Supplementary Figure S14).\nFigure S14. Age-related change in explained variance and state occurrence of propagation patterns. (A) Age-related\nexplained variance ratio (EVR) of dynamic propagation patterns. This panel displays the explained variance ratio EVR for each of\nthe three primary dynamic propagation patterns as a function of age. Each data point represents an individual's EV ratio. The\noverlaid regression lines and their corresponding correlation values demonstrate how the EV ratio for each pattern changes with\nage, with statistically significant relationships (p < 0.05) marked by an asterisk (*).\nFigure S14. Age-related change in explained variance and state occurrence of propagation patterns. (B) Age-related\nexplained variance proportion between propagation patterns. This panel presents the comparison of explained variance proportions\nbetween specific pairs of dynamic patterns: Pattern S-A relative to Pattern TP-D (left panel), Pattern S-A relative to Pattern SM-V\n(middle panel), and Pattern TP-D relative to Pattern SM-V (right panel). The dashed blue line signifies a proportion of 1, indicating\nequal explained variance between the two patterns being compared. While the proportion between Pattern S-A and Pattern TP-D\ndoes not show a statistically significant age effect, the remaining two pairs (S-A/SM-V and TP-D/SM-V) exhibit significant age\neffects. This analysis addresses potential transitions or reversals in pattern dominance with age.\nFigure 2. The authors showed the developmental effects in subgroups. Why did they not\nshow continuous, age-related changes as in Figures 1 and 3?\nResponse: We thank the reviewer’s comment. Our initial rationale for using subgroups\nin Figure 2 was to identify potential stage-specific differences associated with key\nneurodevelopmental milestones. Here, to address the reviewer’s point, we have added\nthe continuous, age-related changes in the revised Supplementary Figure S14C (shown\nbelow as well).\nFigure S14. Age-related change in explained variance and state occurrence of propagation patterns. (C) Age-Related State\nOccurrence Ratio of Each Propagation Pattern. This panel depicts the continuous age-related changes in the state occurrence ratio\nfor each individual propagation pattern (Pattern S-A, Pattern TP-D, and Pattern SM-V). The scatter plots, alongside their regression\nlines and correlation values, illustrate how the proportion of time individuals spend in each dominant propagation state evolves\nacross the development from children to early adulthood.\nReviewer #2 (Remarks to the Author):\nThis work investigates how dynamic propagation patterns of brain activity evolve across\ndevelopment and adulthood. Utilizing large-scale fMRI datasets (HCP and NKI), the\nauthors apply a novel complex PCA (cPCA) framework to extract spatiotemporal patterns\nof brain activity. They identified three core propagation patterns—Sensorimotor-Visual\n(SM-V), Task-Positive–Default (TP-D), and Sensorimotor–Association (S-A)—and\nrevealed robust, age-related changes in their occurrence and directionality. These\nchanges track well with cognitive development (top down more than bottom up\npropagation), supporting the hypothesis that brain dynamics become more structured and\ntop-down driven over the lifespan. The paper is methodologically innovative, with\nthorough replication and sensitivity analyses.\nThere are a few issues to be addressed:\n1. There are significant changes in overall hierarchical structure derived from static brain\nFC from children to adolescents to adults (see Dong et al, Nguyen et al), which seem to\nalso occur in the propagations identified by the authors. Procrustes alignment, while\nnecessary for cross-group comparison, can distort the results. It’d be good to show the\nraw unaligned patterns and their explained variance to show that Procrustes alignment\nonly affects the sign/scale of the data value and nothing else.\nResponse: We agree that Procrustes alignment can potentially distort certain aspects of\nthe data, especially when examining individual-level variation. Specifically, Procrustes\nalignment involves reflection, rotation, and scaling. Reflection and rotation adjust the sign\nand order of components, while scaling can stretch or compress the component\nmagnitude. In our study, we disabled scaling (i.e., scale factor was fixed to 1), restricting\nthe alignment to rotation and reflection only. This ensures that the scale of the component\nis preserved for each individual. Notably, in eigenvector decomposition, the sign of an\neigenvector is arbitrary; only the axis matters. Hence, reflection was permitted. Rotation\nwas used to match individual components to the group-level reference patterns (S-A, TP-\nD, SM-V). Importantly, we preserved the explained variance of the raw component in this\nprocess. In other words, for each individual, the explained variance reported for Pattern\nS-A (or TP-D, SM-V) corresponds to the explained variance of its raw component from\ncomplex PCA, independent of whether it originally appeared as the first, second, or third\ncomponent.\nTo directly address the reviewer’s comments, we visualized the first three raw, unaligned\ncomponents at each developmental stage (revised Supplementary Figure S16, panel A\nshown below), as well as the raw component order for each individual (panel B shown\nbelow). As shown in Figure S16 below, the spatial pattern of the raw components Patterns\n1-3 is highly consistent across age groups, aside from sign, which is arbitrary. Across\nindividuals, 93.87% individuals show Pattern S-A as their first raw component; 61.51%\nindividuals have Pattern TP-D identified as their second component; and 46.32% of\nPattern SM-V as their third.\nFigure S16. (A) The first three raw patterns at the phase of maximal spatial variance for each age group. Participants from the HCP-\nD dataset were grouped into 1-year age bins from 8 to 21 years old. (B) The raw order of Pattern S-A, Pattern TP-D and Pattern SM-\nV for each individual. Notably, 93.87% of individuals exhibit Pattern S–A as the first component, 61.51% show Pattern TP–D as the\nsecond, and 46.32% show Pattern SM–V as the third.\nIn addition, we examined the developmental effects of the raw explained variance ratio\n(EVR) of each pattern (revised Figure S14A, shown below). Consistent with prior work on\nstatic gradients (e.g., Dong et al.), the explained variance of Pattern S-A increases with\nage, while that of Pattern SM-V decreases. To directly assess whether Pattern S-A\nremains the dominant pattern during childhood, we calculated the within-subject ratio of\nEVR between Pattern S-A and Pattern TP-D (figure below, panel B). The ratio is greater\nthan 1 (blue dashed line) in the vast majority of participants from age 8, indicating that\nPattern S-A remains the dominant pattern from childhood through adulthood. We\nperformed a similar ratio analysis for Pattern S-A versus Pattern SM-V, with consistent\nresults: Pattern S-A explains the most variance across all ages, and the ratio increases\nsignificantly during childhood (Figure S14 below, panel B, middle). When comparing the\nEVR between Pattern TP-D and Pattern SM-V, Pattern TP-D is generally recognized as\nthe secondary component, although in some individuals it may fall below SM-V (Figure\nS14 below, panel B, right). From childhood to adulthood, Pattern TP-D shows a significant\nage-related increase in its likelihood of being the second dominant component. In the\nrevision, we have included all of these results (Figure S14, Figure S16).\nFigure S14. Age-related change in explained variance and state occurrence of propagation patterns. (A) Age-related\nexplained variance ratio (EVR) of dynamic propagation patterns. This panel displays the explained variance ratio EVR for each of\nthe three primary dynamic propagation patterns as a function of age. Each data point represents an individual's EV ratio. The\noverlaid regression lines and their corresponding correlation values demonstrate how the EV ratio for each pattern changes with\nage, with statistically significant relationships (p < 0.05) marked by an asterisk (*).(B) Age-related explained variance proportion\nbetween propagation patterns. This panel presents the comparison of explained variance proportions between specific pairs of\ndynamic patterns: Pattern S-A relative to Pattern TP-D (left panel), Pattern S-A relative to Pattern SM-V (middle panel), and Pattern\nTP-D relative to Pattern SM-V (right panel). The dashed blue line signifies a proportion of 1, indicating equal explained variance\nbetween the two patterns being compared. While the proportion between Pattern S-A and Pattern TP-D does not show a statistically\nsignificant age effect, the remaining two pairs (S-A/SM-V and TP-D/SM-V) exhibit significant age effects. This analysis addresses\npotential transitions or reversals in pattern dominance with age.\n2. The dominant state was determined using a winner-take-all approach, which resulted\nin straightforward analysis and findings, but might miss out the subtlety of brain activity at\na certain time point, which is a combination of different propagation components\nidentified. Sometimes a state is probably highly dominant, sometimes the difference\nbetween them might not be so clear. Given that you have the magnitude of each state at\na specific time point, a more nuanced method would be good to derive a more accurate\npicture. Importantly, how do these relate to age and cognition?\nResponse: We appreciate the reviewer's comment and agree that the brain state could\nbe a combination of multiple propagation components instead of a single ‘winner’ pattern.\nAs suggested, we replicated our analysis using the magnitude of the propagation patterns\n(i.e. the contribution percentage of each Pattern at a given time point). The results\n(revised Figure S17, shown below) revealed similar age-related changes as Figure 2 in\nthe main manuscript: an increased occurrence in magnitude of both Pattern S-A and\nPattern TP-D across all age groups from childhood to middle-aged adulthood, and a\ndecrease in Pattern SM-V from adolescence to early adulthood. As the reviewer noted,\nthe magnitude-based approach appears to be more sensitive in detecting age-related\nchange in Pattern TP-D, particularly from childhood to early adulthood. However, neither\n‘winner-take-all’ nor magnitude ratio alone yielded a significant prediction of the total\ncognitive score. The successful prediction was only achieved when the spatial pattern of\npropagation was included in the model (Figure 4). In the revision, we have now included\nthese updated results.\nFigure S17. Occurrence contribution of dynamic patterns across age groups. (B) Mean proportion of pattern magnitude across\ndevelopmental stages: children, adolescents, early adulthood in HCP-D cohorts, and adult reference (HCP-A, age < 50). Age effect\nbetween consecutive development stages was tested using t-tests with significance denoted by asterisks (*: p < 0.05).\n3. The interpretation of top-down vs. bottom-up propagation is intriguing but may be\noversimplified. Since the directionality is inferred from relative phase differences, can the\nauthors expand on how confident we can be in attributing cognitive control processes to\nsuch phase shifts?\nResponse: We thank the reviewer for this insightful comment. In our study, we employed\ncomplex PCA, which captures both zero-lag and time-lagged spatiotemporal\npropagations of brain activity. The phase information allows us to infer directional\npropagation across the cortex, for example, along the S-A axis from higher-order\nassociation areas and later in lower-order sensory regions. While this method does not\ndirectly measure causality, the observed directionality is theoretically grounded and has\nbeen shown to be consistent with prior studies (Gu et al, Cereb. Cortex 2021; Bolt et al,\nNat. Neurosci. 2022). In our study the brain-behavior prediction analysis demonstrates\nthat the top-down dynamic propagations are significantly associated with cognitive\nperformance, more than bottom-up propagations - supporting the functional relevance of\nthe directionality of these propagations.\nThat said, we agree that establishing causal attribution to cognitive control is limited with\nresting-state fMRI data and requires task-based validation. For example, recent work by\nPines et al., (2023) used optical flow analysis to investigate propagation directionality\nduring a Go/NoGo task compared to resting-state fMRI (Pines et al, Neuron 2023). They\nfound that task demands evoke an increased prominence of top-down propagations.\nThis supports the idea that top-down propagation may underlie cognitive control.\nIn the revision, we have added the following text in the discussion with a more cautious\ntone and explicitly state the interpretative limitations.\n“In the present study, CPCA captures both zero-lag and time-lagged spatiotemporal\npropagations of the brain activity, which allow us to infer relative directionality of cortical\npropagation. While this does not provide direct evidence of causality to cognitive control,\nour brain-behavior prediction supports the idea that top-down propagations were more\nstrongly associated with cognitive performance than bottom-up propagations. Prior study\nhas also demonstrated increased top-down propagations during cognitive tasks\ncompared to the resting state23.”\n4. Do the propagation patterns seen here similar to what have been identified using a\ndynamic gradient approach (e.g. Brown et al, Neuroimage 2022)? It would be good to\ndiscuss the implications.\nResponse: Thank you for the comment and reference on the dynamic gradient approach.\nBrown et al. (2022) used gradients to construct a low-dimensional space and tracked\ntemporal changes of fMRI BOLD activity in this space. In our study, we characterize\nsimilar dynamic changes in a low-dimensional space, as spatiotemporal patterns using\ncomplex PCA (CPCA). Importantly, CPCA identifies both spatial and temporal changes\nsimultaneously, and it characterizes the temporal changes in the phase domain without\nrequiring a fixed time window. This offers greater flexibility to examine the individual\nvariability, particularly well-suited for developmental data. In the revision, we have added\nthe discussion as suggested, as follows.\n“Compared to traditional methods, such as the sliding-window approach or quasi-periodic\npattern analysis with a predefined window size8, CPCA demonstrates superior precision\nin detecting the individual dynamic patterns because it doesn't need a fixed arbitrary\nwindowing across individuals from different developmental stages. Additionally, unlike\nmethods that decouple the temporal dynamic changes in a predefined low-dimensional\nspace (e.g., functional gradients)24,54, CPCA extracts spatiotemporal patterns by\nsimultaneously identifying spatial topographies and their corresponding temporal phase\nrelationships. This joint estimation allows CPCA to capture both spatial and temporal\ndynamic changes within a unified framework.”\n5. The discussion mainly focuses on previous findings on static FC and gradient. I would\nbe interested to see at least a point on how the CPCA-derived brain states can be\ninterpreted specifically in comparison to the brain state findings using traditional dynamic\nmethods like HMM.\nResponse: We appreciate the reviewer’s comment. As the reviewer suggested, we\nadded the following text in the discussion in our revision.\n“Traditional approaches, such as Hidden Markov Models, characterize brain dynamic\nchanges as transitions between a set of latent states. However, each of these brain states\nis distinct zero-lag spatial pattern. In contrast, CPCA enables the identification of both\nzero-lag and time-lagged directed propagation patterns. Importantly, prior work by Bolt et\nal. showed that CPCA patterns are closely aligned with those derived from methods\nincluding ICA, PCA, and Hidden Markov Models10, suggesting CPCA captures\nrepresentative brain functional dynamics. Overall, these combined advantages make\nCPCA particularly flexible and suitable for studying functional dynamics, especially in\nneurodevelopmental samples, where spatiotemporal variability is often more\nheterogeneous10,55.”\nReviewer #2 (Remarks on code availability):\nThe code looks good to me.\nResponse: Thank you for your positive feedback.\nReviewer #3 (Remarks to the Author):\nI co-reviewed this manuscript with one of the reviewers who provided the listed reports.\nThis is part of the Nature Communications initiative to facilitate training in peer review and\nto provide appropriate recognition for Early Career Researchers who co-review\nmanuscripts.\nResponse: We are grateful for your time and effort in co-reviewing our manuscript.\nResponse:\nWe thank the reviewers for their comments and suggestions. In this revision, we have\naddressed editorial requirements and production guidelines, including: 1) enhancing\nfigure documentation and statistical reporting, 2) providing complete source data and\ncode availability through GitHub, and 3) ensuring full compliance with reporting and\nethics standards. We enclosed the Author Checklist for your review."}, {"round": 2, "reviewer_comments": {"reviewer_1": "The authors have nicely addressed my previous concerns and I dont have further comments. Congratulations.", "reviewer_2": "All comments have been addressed. I thank the authors for the comprehensive experiments and discussion.", "reviewer_3": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts."}, "author_response": ""}]

Decision letter