Resting-state functional connectivity is a widely used approach to study the functional brain network organization during early brain development. However, the estimation of functional connectivity networks in individual infants has been rather elusive due to the unique challenges involved with functional magnetic resonance imaging (fMRI) data from young populations. Here, we use fMRI data from the developing Human Connectome Project (dHCP) database to characterize individual variability in a large cohort of term-born infants (N = 289) using a novel data-driven Bayesian framework. To enhance alignment across individuals, the analysis was conducted exclusively on the cortical surface, employing surface-based registration guided by age-matched neonatal atlases. Using 10 minutes of resting-state fMRI data, we successfully estimated subject-level maps for eight brain networks along with individual functional parcellation maps that revealed differences between subjects. We also found a significant relationship between age and mean connectivity strength in all brain regions, including previously unreported findings in higher-order networks. These results illustrate the advantages of surface-based methods and Bayesian statistical approaches in uncovering individual variability within very young populations.

The human brain undergoes significant changes during late gestation and early infancy. Characterizing the functional brain organization during early postnatal ages is clinically relevant to elucidate how different genetic factors and environmental hazards may impact subsequent development. Resting-state functional connectivity (RSFC) is a widely used approach to delineate the functional brain organization in the perinatal period since subjects are imaged at rest without requiring the performance of any specific task. Typically, RSFC methods look for temporal coherence of spontaneous fluctuations of functional magnetic resonance imaging (fMRI) signals to define resting-state networks (RSNs), also known as the functional connectome. Seminal studies based on RSFC methods have shown that the functional connectome during early development is affected by premature birth (Doria et al., 2010; Smyser et al., 2010), maternal stress (DeSocio, 2018; Sandman et al., 2012), developmental dyslexia (Cao et al., 2014; Yu et al., 2022), and prenatal drug exposure (Merhar et al., 2021), among other factors. While these studies have provided important insights into the typical and atypical functional brain organization in early life, the vast majority have focused on the estimation of RSNs at the group level (i.e., derived from multi-subject data). In recent years, it has been shown that the functional connectome is associated with later cognitive and behavioral outcomes (He et al., 2018). Therefore, individual characterizations are critical for making subject-based predictions of clinical relevance.

The robust estimation of RSNs in individual infants has been rather elusive due to the inherently low signal-to-noise ratio (SNR) of the resting-state fMRI data and the unique sources of noise (e.g., idiosyncratic head motion, faster cardiac and respiration rates) in this population (Smyser & Neil, 2015). Traditional RSFC analysis approaches like group independent component analysis (gICA) produce clean estimations of RSNs when combining data from a large number of infants (Eyre et al., 2021; Gao et al., 2015) but fail to capture individual differences. In contrast, estimations of subject-level RSNs obtained from methods such as dual regression have low power to mitigate noise from single-subject data. Precision functional mapping is an emerging technique that offers an effective strategy to increase the SNR in functional connectivity analysis within individuals, maximizing the power to detect individual differences in RSNs (Gordon et al., 2017; Gratton et al., 2020, 2022; Laumann et al., 2015). However, to achieve an acceptable level of reliability at the subject level, extensive amounts of imaging data need to be collected, which is rarely practical or cost-effective when considering infant and clinical populations. Bayesian approaches offer a powerful complement to longer scans by leveraging shared information across subjects from a representative population, which reduces noise while enabling individual differences to be expressed without requiring extended scans (Kong et al., 2019; Mejia et al., 2020).

In this study, we explore hierarchical Bayesian modeling in combination with surface-based analysis to improve individual estimations of RSNs in infants during the first weeks of life. The proposed Bayesian framework (Mejia et al., 2020) requires population-derived priors or templates, so we leverage the information contained in the large open-access dataset from the developing Human Connectome Project (dHCP) (Edwards et al., 2022). A critical assumption of the Bayesian framework is that all the subjects are anatomically co-registered to a common atlas space. Thus, to optimize the anatomical alignment across individuals with different gestational ages, the analysis is entirely done on the cortical surface using registration techniques guided by age-matched infant atlases. While surface-based registration approaches have become increasingly popular in studies involving older populations (Glasser, Smith, et al., 2016), their application to early brain development has been limited due to the challenges involved in the construction of accurate cortical surfaces and age-appropriate atlases. Surface-based registration approaches have become increasingly popular in studies of functional connectivity among adult populations (Glasser, Smith, et al., 2016). However, their application to early brain development is just now gaining attention. Recently, Wang et al. (2023) used local functional connectivity gradient maps to improve functional alignment across subjects after cortical folding-driven registration, enabling them to generate cortical parcellation maps in a cohort of infants/toddlers. Another study by Hu et al. (2022) employed the same folding-driven registration technique to demonstrate the existence of functional networks in neonates, as well as their longitudinal stability throughout the first two years of life.

Our analysis of a cohort of 289 term-born neonates (age at scan: 37.4–44.8 weeks) shows that the surface-based hierarchical Bayesian approach can produce clean estimates of subject-level cortical RSNs using only 10 minutes of individual fMRI data. The statistical framework enables the computation of subject-level t-statistic maps which, in turn, allows for functional parcellations at the subject level. As a result, we observed individual topographical differences hidden by group-level averages. Furthermore, a positive relationship between age and subject-level connectivity strength was revealed for almost all cortical RSNs, confirming the hypothesis that the functional connectome matures with age in infants. Importantly, our work extends beyond existing research by showing that significant maturational changes are not only restricted to primary sensory networks but also present in higher-order networks in the first weeks of postnatal life. These results illustrate the advantages of combining surface-based processing and hierarchical Bayesian approaches to inform individual variability in very young populations, opening the door to precision neuroimaging studies of early brain development with enhanced accuracy and reliability.

2.1 Subjects and data acquisition

MRI data were obtained from the second release of the dHCP database. For this study, we only considered term-born infants (i.e., gestational age (GA) 37 weeks) with radiological scores lesser than three, indicating no lesions of clinical or analytical significance (see dHCP release notes1 for details). Following these criteria, 305 term-born infants (age at birth: 37.1–42.3 weeks gestational age, GA) scanned shortly after birth (age at scan: 37.4–44.8 weeks postmenstrual age, PMA) were considered for further analysis.

All scans were obtained with a 3T Philips Achieva using a dedicated neonatal head coil at Evelina Newborn Imaging Centre, St. Thomas Hospital, London, UK. Both T1-weighted (TR = 4795 ms; TE = 8.7 ms) and T2-weighted (TR = 12 s; TE = 156 ms) structural scans were obtained with a multi-slice Turbo Spin Echo (TSE) sequence, with in-plane resolution 0.8 × 0.8 mm2 and 1.6 mm slices overlapped by 0.8 mm. Two stacks of images were taken per weighting, sagittal and axial, which were integrated to obtain T2w volumes with an image resolution of 0.8 mm isotropic. Blood oxygen level-dependent (BOLD) scans were obtained with a multi-slice gradient-echo echo planar imaging (EPI) sequence (TE = 38 ms; TR = 392 ms, multiband factor = 9; flip angle = 34°) with an image resolution of 2.15 mm isotropic. A resting-state BOLD fMRI acquisition of 2300 time points (15 minutes) was obtained for each infant.

2.2 MRI preprocessing

2.2.1 Cortical surfaces

After brain extraction using a modified version of FSL BET for unmyelinated brains, tissue-specific masks corresponding to gray matter, low-, and high-intensity white matter were segmented from the brain-extracted T2w images using the Draw-EM tool (Makropoulos et al., 2014). Individual surfaces corresponding to white matter, pial, and midthickness were created using these tissue masks as described in Makropoulos et al. (2018).

For the rest of the analysis, a 40-week symmetrical atlas consisting of 57,700 vertices after excluding the medial wall built from the dHCP cohort (Williams et al., 2023) was used. The original atlas in Bozek et al. (2018) was built based on a spherical registration approach using the Multimodal Surface Matching (MSM) tool for aligning cortical folding (Robinson et al., 2018). For this extended version, left-right vertex correspondence was enforced by co-registering right and left sulcal depth maps during alignment.

We used the MSMSulc approach to register each individual cortical surface to the symmetrical 40-week atlas, concatenating single-week transformations to avoid the accumulation of approximation errors. MSMSulc uses a highly regularized folding-based registration that has been shown to improve functional overlap among subjects (Glasser, Coalson, et al., 2016; Robinson et al., 2014). These transforms were calculated using the scripts from the dHCP GitHub repository2. The spherical registrations from this process were later used when mapping the BOLD fMRI data onto the atlas surface.

2.2.2 Functional data

Resting-state BOLD fMRI data were preprocessed using a modified version of the Human Connectome Project (HCP) (Glasser et al., 2013) pipeline specifically developed for the dHCP (Fitzgibbon et al., 2020). Preprocessing steps were performed on the BOLD volumes in native space and included: (1) field map correction, (2) intra- and inter-volume motion correction, (3) high-pass filtering to remove slow drifts, and (4) spatial ICA using FSL FIX (Griffanti et al., 2014; Salimi-Khorshidi et al., 2014) to regress out structured noise artifacts.

To further reduce the potential effect of motion on functional connectivity measures, we adopted a conservative block censoring approach. Emerging evidence suggests that head motion is correlated with neurophysiological states (e.g., sleep state) (Denisova, 2019). Since aggressive frame censoring may introduce heterogeneities in the BOLD time series associated with different physiological states, we retained a contiguous block of data to minimize this potential issue. Following the criteria proposed in Eyre et al. (2021), frames with DVARS (root-mean-squared BOLD signal intensity change) higher than 1.5 times the interquartile range above the 75th percentile within a session were considered corrupted by motion. A contiguous block of 1600 frames (10 minutes) with the minimum number of motion outliers was retained for each subject. Furthermore, any subject with more than 160 motion-corrupted frames (10%) within the contiguous block was excluded. By following this approach, we were able to retain 274 subjects with 10 minutes of resting-state BOLD fMRI data for further analysis.

The volumetric resting-state BOLD data were mapped onto the individual cortical surfaces using the HCP surface pipeline (Glasser et al., 2013) as implemented in Connectome Workbench. The surface-based BOLD data was then mapped onto the custom 40-week atlas surface using spherical registration (as described in Section 2.2.1). Finally, the individual datasets were spatially smoothed using a geodesic 2D Gaussian kernel (FWHM = 3 mm). In the following analyses, only the cortical data were taken into consideration.

2.2.3 Temporal signal-to-noise ratio

To quantify image quality throughout the cortex and mask out noisy vertices, we computed a measure of temporal signal-to-noise ratio (tSNR) as the ratio between the mean BOLD signal and its standard deviation. A global tSNR map was built in decibels [dB] by averaging the individual tSNR maps across the whole sample as

(1)

where, for each subject i, the SNR is calculated for each vertex v.

As expected from the patterns associated with susceptibility artifacts in BOLD data and previous similar analyses (Fitzgibbon et al., 2020; Thomas Yeo et al., 2011), we observe areas of particularly low tSNR in the inferior and medial temporal lobe, orbitofrontal and insular cortices (Fig. S3).

2.3 Bayesian subject-level ICA

To estimate subject-level RSNs, we use a custom implementation of template ICA (Mejia et al., 2020), an extension of probabilistic ICA (Beckmann & Smith, 2004). In this section, we describe the template-based ICA model, the calculation of population-derived priors, and the strategy for estimating individual posterior mean and variance.

2.3.1 Template-based ICA model

For each subject i, the dimensionally reduced BOLD fMRI data yi with T time-points is modeled at vertex v as

(2)

where Ai(T×Q) is the mixing matrix, si(Q×1) represents the Q independent components (ICs) associated with a different RSN, and ei(T×1) represents any source of Gaussian noise with variance Ω. Like in probabilistic ICA, yi is obtained using singular value decomposition (SVD) of the preprocessed BOLD dataset, and an estimate of the Gaussian noise variance Ω is obtained as the residual variance. Note that the error in equation 2 has a Gaussian distribution since we assume that any source of non-Gaussian noise is eliminated after ICA-FIX on each preprocessed BOLD dataset.

In the template ICA model, we further characterize each IC as

(3)

where μ0(v)(Q×1) is the group-level mean, and δi(v)(Q×1) represents subject-level deviations denoting unique features of spatial topography associated with each RSN. We assume that the δi(v) are normally distributed with covariance matrix Σ0(v)=diag{ν12(v),ν22(v),...,νq2(v)} where each σ0,q2 is the between-subject variance associated to each IC q at each vertex v.

For a subject-specific BOLD dataset yi(v), the proposed problem can be formulated within a hierarchical Bayesian framework where we aim to obtain a vertex-wise posterior distribution for the subject-level IC maps si(v), given a set of population priors or “templates” (mean μ0(v) and between-subject variance Σ0(v)). Note that given the Gaussian prior and likelihood, the vertex-wise posterior spatial maps are normally distributed with a subject-level mean μi and variance σi that have analytical forms. To estimate the posterior mean and variance of the subject-level ICs, along with the model parameters including the mixing matrix and residual variance, a computationally efficient expectation-maximization (EM) algorithm is utilized (see Section 2.3.3 and Figure 1 for a description of the estimation procedure).

Fig. 1.

Estimation of subject-level IC maps. To estimate the posterior mean and variance of the subject-level ICs, along with the model parameters including the mixing matrix and residual variance, a computationally efficient expectation-maximization (EM) algorithm is used. An initial estimation of the model parameters (Ωi) is obtained from dimensionality reduction (yi) of the BOLD timeseries, while the mixing matrix (Ai(0)) is obtained from a rough dual regression estimate of the individual IC maps. The empirical population prior parameters, namely, the mean (μ0(v)) and the between-subject variance σ0(v) of the group-level IC maps are calculated from a uniform sample of the cohort. After the EM algorithm converges, an estimate of the subject-level mean (μ^i,q(v)) and standard deviation (σ^i,q(v)) for each IC is obtained. A binary mask of significant engagement for each IC is obtained from a Bayesian hypothesis test after correcting for multiple comparisons on the individual mean and variance IC maps.

Fig. 1.

Estimation of subject-level IC maps. To estimate the posterior mean and variance of the subject-level ICs, along with the model parameters including the mixing matrix and residual variance, a computationally efficient expectation-maximization (EM) algorithm is used. An initial estimation of the model parameters (Ωi) is obtained from dimensionality reduction (yi) of the BOLD timeseries, while the mixing matrix (Ai(0)) is obtained from a rough dual regression estimate of the individual IC maps. The empirical population prior parameters, namely, the mean (μ0(v)) and the between-subject variance σ0(v) of the group-level IC maps are calculated from a uniform sample of the cohort. After the EM algorithm converges, an estimate of the subject-level mean (μ^i,q(v)) and standard deviation (σ^i,q(v)) for each IC is obtained. A binary mask of significant engagement for each IC is obtained from a Bayesian hypothesis test after correcting for multiple comparisons on the individual mean and variance IC maps.

Close modal

2.3.2 Population priors

The Bayesian framework relies on population-derived priors, which consist of estimates of the mean and between-subject variance of a set of ICs. To first define the set of ICs, we perform group ICA (Beckmann et al., 2005) on a subset of 24 subjects (scanned between 43.5 and 45 weeks PMA) using FSL MELODIC (Beckmann & Smith, 2004). The number of ICs is set at 20, excluding subcortical regions, to achieve a balance between robustness and similarity to previous analysis on infants (Doria et al., 2010; Eyre et al., 2021; Rajasilta et al., 2020; Toulmin et al., 2015). Eight ICs were associated with neurologically relevant cortical networks, namely, medial motor, lateral motor, auditory, somatosensory, primary visual, motor association, visual association, and default mode network (DMN) (1). Examples of the time courses and power spectra associated with signal and nuisance components are included in Figure S6. To obtain noisy pseudo test-retest point estimates of individual IC maps that can be used to obtain population means and between-subject variance estimates, dual regression (Nickerson et al., 2017) is applied to a subset of subjects (N = 35) uniformly scanned between 37 and 43 weeks PMA. To calculate the between-subject variability, each individual BOLD time series is divided into two halves (or pseudo-sessions) prior to dual regression. The population mean μ0(v) for each IC q is estimated as

(4)

where sij(q,v) is the dual regression IC estimate for each subject i and each session j.

The between-subject variance is obtained by decomposing the total variance into within- and between-subject components. The total variance for each individual IC at each vertex was estimated as the average of the variance obtained for each pseudo-session,

(5)

while the within-subject variance can be estimated based on the variance of the difference between the individual ICs from each pseudo-session

(6)

Finally, from equations 5 and 6, the between-subject variance is estimated as the difference between the total and within-subject variances

(7)

where the σ0q are the diagonal elements of the covariance matrix Σ0 in equation 2, representing the variance of the empirical population prior. The resulting empirical population priors are shown in Figure S4.

2.3.3 Estimation of individual posterior mean and variance

A schematic of the estimation procedure is shown in Figure 1. The EM algorithm requires an initial estimate of the parameters Ω and Ai. Ai(0) was obtained from the dual regression procedure described above. Gaussian noise variance σ0(0) was estimated from the residual variance of the dimensionality reduction of the BOLD data.

These two parameters are updated at every iteration of the EM algorithm by maximizing the expected log-likelihood. Alternately, conditional on the current estimate of the parameters, the posterior distribution of the subject-level ICs si is updated to obtain the posterior moments required for the parameter maximum likelihood estimations. The algorithm runs until convergence of the parameter estimators, and the posterior mean μ^i and variance σ^i of s^i(v,y,θ) are obtained. A schematic representation of the entire pipeline is shown in Figure S1.

3.1 Individual patterns of functional connectivity

Using only 10 minutes of resting-state fMRI data, we identified eight cortical RSNs (represented by independent components, ICs) in individual subjects from a cohort of 239 full-term infants (note that the 35 subjects considered in the estimation of the empirical population priors were excluded from further analysis). Overall, the surface-based Bayesian approach produced cleaner estimates of subject-level IC maps than dual regression, as illustrated by five representative IC maps obtained in a single infant (Fig. 2). Notably, the Bayesian approach removed clusters of spurious activations present in the dual regression maps (see, for example, the temporal cortex in the somatosensory IC maps of Fig. 2). Other improvements include the emergence of more defined clusters in distributed networks, as shown in the prefrontal cluster of the default mode IC maps. A complete set of subject-level IC maps obtained for both methods is included in Figure S2.

Fig. 2.

Subject-level maps obtained from 10 minutes of resting-state fMRI data. (A) Five representative IC maps obtained from dual regression for a term-born infant scanned at 42.6 weeks PMA. (B) Posterior mean IC maps derived from the Bayesian model produced cleaner maps than the dual regression approach. White arrows highlight some areas of notable improvement in cluster convexity (e.g., default mode network) or reduction of spurious engagement (e.g., somatosensory network and inferior parietal cluster of the DMN). All maps are projected onto the inflated 40-week atlas.

Fig. 2.

Subject-level maps obtained from 10 minutes of resting-state fMRI data. (A) Five representative IC maps obtained from dual regression for a term-born infant scanned at 42.6 weeks PMA. (B) Posterior mean IC maps derived from the Bayesian model produced cleaner maps than the dual regression approach. White arrows highlight some areas of notable improvement in cluster convexity (e.g., default mode network) or reduction of spurious engagement (e.g., somatosensory network and inferior parietal cluster of the DMN). All maps are projected onto the inflated 40-week atlas.

Close modal

To perform comparisons across subjects, we computed t-statistic maps from the subject-specific posterior mean and standard error maps estimated by the surface-based Bayesian model. Figure 3 shows the unthresholded t-statistic maps for five networks and three term-born infants scanned at different postmenstrual ages within a period of four weeks. While a general spatial topography is preserved for each specific network, individual differences are evident across subjects. For example, the temporal cluster of the default mode IC map of subject A is more defined than in subjects B and C. Similarly, the prefrontal cluster of the same IC map of subject B shows a stronger level of engagement in comparison with subjects A and C.

Fig. 3.

Subject-level maps across different subjects. Unthresholded t-statistic maps were computed from the mean and standard error maps derived from the Bayesian model estimation. The five ICs shown in Figure 2 are displayed for three term-born infants scanned at 40.6 weeks PMA, 41.7 weeks PMA, and 43.9 weeks PMA, respectively. Although a general spatial topography is preserved within each specific network, there are evident variations across individuals (as indicated by the white arrows).

Fig. 3.

Subject-level maps across different subjects. Unthresholded t-statistic maps were computed from the mean and standard error maps derived from the Bayesian model estimation. The five ICs shown in Figure 2 are displayed for three term-born infants scanned at 40.6 weeks PMA, 41.7 weeks PMA, and 43.9 weeks PMA, respectively. Although a general spatial topography is preserved within each specific network, there are evident variations across individuals (as indicated by the white arrows).

Close modal

As an alternative way to summarize the individual differences observed in Figure 3, we used a winner-takes-all strategy to obtain subject-level functional parcellations (Fig. 4) in which a network label was assigned to each vertex based on the highest t-score at that location (see Fig. S5 for a comparison between the individual IC maps and their associated WTA parcellations). The Bayesian approach constitutes a compromise between robust group-level RSNs provided by the population template and subject-level RSNs associated with the individual data. This can be appreciated by looking at the functional parcellation maps obtained for the group ICA maps and the subject-level parcellations for three different infants (Fig. 4). Even though a general resemblance to the group-level parcellation is preserved, remarkable individual differences in shape, size, and location of the individual parcels are evident across subjects. See, for example, the significantly reduced temporal cluster of the default mode IC map in the right hemisphere of Subject B in comparison with the other individuals.

Fig. 4.

Group vs. individual parcellations. Functional parcellations were obtained using a winner-takes-all approach, dividing the cortex into eight distinct brain networks. The group parcellation was obtained from the Z-score maps derived from the group ICA analysis on 24 subjects. The individual parcellations were obtained from the subject-level t-maps derived from the Bayesian estimation of the subjects shown in Figure 3. White arrows highlight notable topographical differences between subjects. The results are displayed only in regions where the temporal signal-to-noise ratio (tSNR) is greater than 17 dB (see Fig. S3 for a whole depiction of the global tSNR computed for this cohort).

Fig. 4.

Group vs. individual parcellations. Functional parcellations were obtained using a winner-takes-all approach, dividing the cortex into eight distinct brain networks. The group parcellation was obtained from the Z-score maps derived from the group ICA analysis on 24 subjects. The individual parcellations were obtained from the subject-level t-maps derived from the Bayesian estimation of the subjects shown in Figure 3. White arrows highlight notable topographical differences between subjects. The results are displayed only in regions where the temporal signal-to-noise ratio (tSNR) is greater than 17 dB (see Fig. S3 for a whole depiction of the global tSNR computed for this cohort).

Close modal

3.2 Effect of age at scan

To explore whether the observed differences are related to maturational factors, we computed the Spearman’s rank correlation coefficient between age at scan and individual connectivity strength for each network, after controlling for sex and motion. Statistical significance was defined as monotonic trends with a p-value lower than 0.05, as calculated using the AS 89 algorithm (Best & Roberts, 1975) implemented in the ‘stats’ R-package. Individual connectivity strength was defined as the mean t-score within the areas of significant engagement or IC mask, for each of the subject-level ICs. Areas of significant engagement were estimated using a Bayesian hypothesis test based on each pair of posterior mean and standard error IC maps after Bonferroni correction for multiple comparisons (as shown in the pipelines of Fig. 1 and Fig. S1). This test is based on the Bayesian credible intervals of the posterior distribution, which is Gaussian in this case (Kruschke, 2011). To account for variable vertex area across the cortical surface, the individual connectivity strength was also weighed by the midthickness vertex area within each IC mask. Within-network connectivity significantly increased with age at scan (37.4–44.8 weeks PMA) for all cortical ICs, including medial motor, lateral motor, somatosensory, auditory, primary visual, motor association, default mode, and visual association networks (Fig. 5). As in similar studies (e.g., Eyre et al., 2021), we present the uncorrected p-values for multiple comparisons across ICs. The full results, including corrected p-values for all tested ICs, are available in Table S1.

Fig. 5.

Effect of age at scan. The relationship between age at scan and individual connectivity strength for each network was assessed by a Spearman’s rank correlation test after controlling for sex and motion. All networks exhibit significant effects with age (p<0.05), as shown by Spearman’s r and uncorrected p-values. A dashed blue line showing a linear trend is also included for illustration purposes.

Fig. 5.

Effect of age at scan. The relationship between age at scan and individual connectivity strength for each network was assessed by a Spearman’s rank correlation test after controlling for sex and motion. All networks exhibit significant effects with age (p<0.05), as shown by Spearman’s r and uncorrected p-values. A dashed blue line showing a linear trend is also included for illustration purposes.

Close modal

3.3 Inter-subject variability

To quantify the spatial variability across individuals, we estimated a frequency map from the subject-level parcellations. Figure 6 shows the percentage of subjects that share the same label or parcel at each cortical vertex. Given that multiple labels may contribute to the same vertex in different subjects, the frequency map only shows the percentage of subjects associated with the dominant label at each location. In addition to all primary networks, the visual association networks and the posterior nodes of the default mode network appear clearly defined in the frequency map and exhibit a high spatial overlap. In contrast, the clusters associated with the prefrontal node of the default mode network and motor association networks exhibit lower spatial overlap across individuals.

Fig. 6.

Frequency map of the entire cohort. Each vertex represents the percentage of subjects that share the same parcel or network. Networks with higher spatial overlap such as the primary sensory RSNs are represented by warm colors. Note that boundaries between networks and areas with low SNR show higher variability (i.e., low spatial overlap) across subjects.

Fig. 6.

Frequency map of the entire cohort. Each vertex represents the percentage of subjects that share the same parcel or network. Networks with higher spatial overlap such as the primary sensory RSNs are represented by warm colors. Note that boundaries between networks and areas with low SNR show higher variability (i.e., low spatial overlap) across subjects.

Close modal

This study extends beyond existing research by characterizing RSNs in term-born neonates at the subject level, leveraging several technical improvements to elucidate individual differences and developmental trajectories within the first weeks of life. Although the existence of resting-state networks (RSNs) is well documented in adults, several neonatal studies have demonstrated the presence of recognizable, albeit sometimes incomplete, RSNs at the group level (Doria et al., 2010; Eyre et al., 2021; Gao et al., 2015; Smyser et al., 2010). Our study presents a more complete picture at the individual level of the functional connectome at the time of birth. To achieve these results, we utilized surface-based methods and age-matched atlases to mitigate co-registration errors and partial volume effects. In addition, we adopted a Bayesian approach that offers a powerful statistical framework to obtain clean estimates of subject-level RSNs using a limited amount of individual data but borrowing strength from population-derived priors. We hypothesize that the cleaner subject-level estimates of RSN maps compared to those obtained with dual regression (Fig. 2) led to stronger associations between age and functional connectivity strength. An in-depth comparison between the Bayesian method and dual regression can be found in the original methods article (Mejia et al., 2020).

Using these methodological innovations, we were able to identify topographical differences across subjects hidden by group averages. This is in line with studies done in adults that show that areas of common activation obscure individual connectivity features (Gordon et al., 2017; Hermosillo et al., 2024; Michon et al., 2022). This is especially important because individual parcellations derived from individual connectivity estimations are susceptible to suffer from co-registration misalignment between subjects (Bijsterbosch et al., 2019). Indeed, inter-subject differences are evident across the individual t-maps (Fig. 3) and functional parcellations (Fig. 4) obtained for different infants. Interestingly, the individual parcellations also show topographical differences with respect to the group-level parcellation. As an example, note the individual differences in the temporal and prefrontal clusters of the default mode network; while the group presents a single representation (Fig. 4 Group), different subjects present variable cluster shape and sizes (Fig. 3, Subjects A-C). It is notable that despite these differences, the results show adult-like higher-order RSNs maps already present at birth, albeit sometimes in sparser, precursory form (Fig. 4).

To quantitatively assess these individual differences, we studied the relationship between subject-level RSNs and age at scan. As in previous reports, significant maturation is observed in all primary sensory RSNs, as measured by the correlation between network strength and age at scan (Fig. 5) (Doria et al., 2010; Eyre et al., 2021; Gao et al., 2015; Hu et al., 2022; Smyser et al., 2010). Furthermore, contrary to prior studies, we note a clear developmental trajectory within a five-week window for most of the independent components (ICs) linked to higher-order RSNs, except for the prefrontal and posterior cingulate cortices. For the first time to our knowledge, we also observe a clear developmental trajectory within a five-week window for most of the independent components (ICs) linked to higher-order RSNs, except for the prefrontal and posterior cingulate cortices.

The frequency maps depicting the spatial overlap of networks across subjects reveal a pattern that progresses from posterior to anterior regions (Fig. 6). In particular, we observe high levels of spatial overlap in the primary sensory networks, along with the higher-order networks located in the posterior regions. In contrast, the prefrontal cluster of the DMN displays low spatial overlap which can be attributed to developmental factors, as supported by histological findings (Keunen et al., 2017). Notably, there are other regions characterized by low spatial overlap, namely, the boundaries between the different RSNs as well as the insular and inferior temporal cortices. The lack of consistency in the delineation of boundaries between the RSNs reflects the topographical variability that our study aims to capture. By contrast, the limited overlap observed in the insula and inferior temporal cortices may be associated with elevated levels of noise in the BOLD signals, as demonstrated by the tSNR map obtained for this specific cohort (Fig. S3).

4.1 Subject-level inference produces more accurate descriptions of brain organization

We hypothesize that variability in the spatial distribution of RSNs across individuals could be the primary factor that hinders the identification of population-level trends and renders growth curves challenging to obtain, especially for distributed higher-order RSNs. In contrast to prior studies, the use of individual regions of interest (ROIs) enabled us to achieve improved outcomes in terms of estimating within-network connectivity strength at the individual level. This is consistent with previous reports showing group-average ROIs that lead to inflated or under-estimated group differences (Levi et al., 2023).

The trend toward individual inference instead of group analysis has emerged as a recent focus in precision neuroimaging. This approach facilitates the acquisition of relevant individual insights in the clinic and also plays a critical role in aggregating individual data without diluting individual features (Gordon et al., 2017; Laumann et al., 2015). To achieve acceptable levels of reliability at the individual level using naive statistical methods, extensive amounts of imaging data need to be collected, which is rarely practical in infants who are typically scanned during natural sleep. As an alternative, we adopt a recently proposed empirical Bayesian model (Mejia et al., 2020) that leverages the growing volume of neonatal data available in existing databases such as the dHCP database. The adopted model produces more accurate subject-level RSN maps by shrinking to the empirical population prior in subject-specific areas of low SNR while maintaining the individual differences. It is worth noting that Bayesian shrinkage pulls the subject-level estimates toward the group average, potentially underestimating individual differences in connectivity strengths (Fig. 5). Therefore, the computed p-values for the age-connectivity strength relationship may be conservative estimates of the true effects.

4.2 Advanced techniques are crucial to maintain precision mapping

Concurrently, we used additional technical enhancements to mitigate any biases that could potentially skew our results toward adult populations. To eliminate any distortions that may arise from mapping neonatal brains onto an adult-based atlas, we utilized a symmetrical atlas (Bozek et al., 2018; Williams et al., 2023) derived from the dHCP cohort. This atlas served as a reference surface onto which the individual BOLD timeseries were projected.

Moreover, the adoption of surface-based analysis effectively reduces the potential influence of partial volume effects encountered in volumetric analysis (Coalson et al., 2018). By applying spatial smoothing on cortical surfaces, contamination between functional regions on opposite sides of a sulcus is mitigated (Brodoehl et al., 2020) because the filtering kernel is applied over 2D geodesic distances, which are more neurobiologically relevant, instead of 3D volumetric distances (Glasser et al., 2013).

Additionally, alignment of the cerebral cortex to the atlas by surface registration algorithms simplifies the 3D problem of volumetric registration to 2D, resulting in improved and more robust co-registration (Robinson et al., 2014). This is especially important in this work because individual parcellations derived from individual connectivity estimations are susceptible to suffer from co-registration misalignment between subjects (Bijsterbosch et al., 2019). In turn, behavioral and age differences have been associated with spatial functional topography, rather than dynamic differences within individual RSNs (Bijsterbosch et al., 2018; Bozek et al., 2018; Kong et al., 2019), meaning that improved alignment between subjects is crucial to improve the accuracy of behavior or age models.

4.3 Methodological considerations

The selection of the number of independent components in the group ICA was based on previous work with this cohort (Eyre et al., 2021; Molloy & Saygin, 2022; Nielsen et al., 2022) as well as adult parcellations (Gordon et al., 2016; Thomas Yeo et al., 2011). While this collection of studies includes both hard and soft parcellations ranging from 7 to 30 components, the exact number that strikes a balance between functional homogeneity and interpretability remains an open debate (Gordon et al., 2016). For consistency with the latest studies on this dataset, where 30 components were used to describe volumetric RSNs (Eyre et al., 2021) that include subcortical structures, we opted for 20 components for our surface-based analysis excluding subcortical regions. As a result, we find that both group ICA and the individual parcellations align with biologically relevant areas previously associated with increased neuronal activity in adult task-based experiments (Nickerson, 2018; Smith et al., 2009).

Another critical methodological consideration in this study was the choice of the co-registration strategy to align individual cortical surfaces. Aligning brain features across subjects is known to significantly impact functional comparisons of the brain, and this issue is further magnified during the perinatal period when the brain undergoes significant growth. Even when alignment is satisfactory, individual differences have been shown to obscure important functional areas in group averages (Gordon et al., 2017, 2023). When we first started this study, we employed the MSM registration method (Robinson et al., 2018) that optimizes cortical folding-driven alignment. However, we observed that certain clusters of higher-order networks were not strongly present in the group ICA results but appeared to emerge in numerous individual ICs. Subsequently, we modified our alignment strategy by increasing the penalization for surface distortion using the ‘MSMSulc’ method (Williams et al., 2023), which enhances functional alignment. As a result, we observed that the group ICA results included a complete and strong default mode network, illustrating the relative sensitivity to alignment strategies across subjects in these types of analyses.

As an alternative strategy, a recent study conducted on a similar neonatal population has explored multimodal alignment methods that complement cortical folding-driven alignment with functional connectivity gradients to refine inter-subject registration (Wang et al., 2023). This study echoes previous efforts in adults (Glasser, Smith, et al., 2016). Applying these methods to our current study presents certain challenges. First, it is uncertain whether the sharpness of the individual functional connectivity gradients is sufficient to drive further alignments due to the relatively smaller number of scans of the dHCP datasets. Additionally, the dHCP acquisition’s multiband factor of 9 results in non-homogeneously lower SNR, further complicating this method (Risk et al., 2021). Due to these complexities, we deemed multimodal registration of the individuals in this dataset to be the scope of future research.

Head movement modulated by arousal states may influence the resting-state functional connectivity patterns in infants (Mitra et al., 2017). Indeed, it is plausible that certain infants may exhibit increased alertness or be in a state of active sleep, leading to increased head movement compared to other infants who are in a quiet sleep state (Denisova, 2019). Applying an aggressive frame censoring approach that selectively retains BOLD data segments with minimal motion can introduce a potential systematic bias toward subjects or data segments associated with a particular state. Thus, to mitigate any potential bias influenced by the arousal state of the infants, we implemented a conservative block censoring approach. This approach guarantees an equal number of consecutive frames for each subject, ensuring that the segments remain qualitatively and quantitatively comparable across the entire cohort.

In this study, we employed a surface-based hierarchical Bayesian framework to leverage shared information across subjects and produced cleaner estimates of individual functional connectivity maps in neonates when compared to classic RSFC analyses. Improved detection of individual differences was evidenced by a maturational path in almost all brain regions. Importantly, our work extends beyond existing research by showing that significant maturational changes are not only restricted to primary sensory networks but also present in higher-order networks in the first weeks of postnatal life. These results illustrate the potential advantages of combining surface-based processing and template-based approaches to inform individual variability in very young populations, opening the door to precision neuroimaging studies of early brain development with enhanced accuracy and reliability.

The neonatal data in this study are part of the second release of the developing Human Connectome Project and are available to download (https://www.developingconnectome.org). The cortical atlas is also publicly available (https://brain-development.org/brain-atlases/atlases-from-the-dhcp-project/cortical-surface-atlas-bozek/). The in-house code used in this study includes shell scripts, python scripts, and R implementations of the templateICAr (https://github.com/mandymejia/templateICAr) and ciftiTools (Pham et al., 2022) libraries. A version of the code, edited for clarity, is available and maintained at https://github.com/FerradalLab/babyBayes.

S.L.F. developed the original idea and research question. S.L.F. and D.D. designed the methodology. A.F.M. and D.D.P. developed the original methods and tools. D.D. designed, coded, and ran the analysis. D.D. and S.L.F. interpreted the results and wrote the manuscript. All authors contributed to editing the manuscript.

The authors declare no competing financial or non-financial interests.

Data were provided by the developing Human Connectome Project, KCL-Imperial-Oxford Consortium funded by the European Research Council under the European Union Seventh Framework Programme (FP/2007-2013) / ERC Grant Agreement no. [319456]. We are grateful to the families who generously supported this work. The non-dHCP pipeline analyses were performed on the Indiana University HPC Cluster (https://uits.iu.edu/). This study was supported by the Indiana Clinical and Translational Sciences Institute (CTSI), funded in part by grant #UL1TR002529 from the National Institutes of Health (NIH). Additional support was provided by NIH grants R01AG083919 and R01EB027119.

Supplementary material for this article is available with the online version here: https://doi.org/10.1162/imag_a_00504.

Beckmann
,
C.
, &
Smith
,
S.
(
2004
).
Probabilistic independent component analysis for functional magnetic resonance imaging
.
IEEE Transactions on Medical Imaging
,
23
(
2
),
137
152
. https://doi.org/10.1109/TMI.2003.822821
Beckmann
,
C. F.
,
DeLuca
,
M.
,
Devlin
,
J. T.
, &
Smith
,
S. M.
(
2005
).
Investigations into resting-state connectivity using independent component analysis
.
Philosophical Transactions of the Royal Society B: Biological Sciences
,
360
(
1457
),
1001
1013
. https://doi.org/10.1098/rstb.2005.1634
Best
,
D. J.
, &
Roberts
,
D. E.
(
1975
).
Algorithm AS 89: The upper tail probabilities of Spearman’s rho
.
Journal of the Royal Statistical Society. Series C (Applied Statistics)
,
24
(
3
), pp. 377–379. https://doi.org/10.2307/2347111
Bijsterbosch
,
J. D.
,
Beckmann
,
C. F.
,
Woolrich
,
M. W.
,
Smith
,
S. M.
, &
Harrison
,
S. J.
(
2019
).
The relationship between spatial configuration and functional connectivity of brain regions revisited
.
eLife
,
8
,
e44890
. https://doi.org/10.7554/eLife.44890
Bijsterbosch
,
J. D.
,
Woolrich
,
M. W.
,
Glasser
,
M. F.
,
Robinson
,
E. C.
,
Beckmann
,
C. F.
,
Van Essen
,
D. C.
,
Harrison
,
S. J.
, &
Smith
,
S. M.
(
2018
).
The relationship between spatial configuration and functional connectivity of brain regions
.
eLife
,
7
,
e32992
. https://doi.org/10.7554/eLife.32992
Bozek
,
J.
,
Makropoulos
,
A.
,
Schuh
,
A.
,
Fitzgibbon
,
S.
,
Wright
,
R.
,
Glasser
,
M. F.
,
Coalson
,
T. S.
,
O’Muircheartaigh
,
J.
,
Hutter
,
J.
,
Price
,
A. N.
,
Cordero-Grande
,
L.
,
Teixeira
,
R. P. A.
,
Hughes
,
E.
,
Tusor
,
N.
,
Baruteau
,
K. P.
,
Rutherford
,
M. A.
,
Edwards
,
A. D.
,
Hajnal
,
J. V.
,
Smith
,
S. M.
, …
Robinson
,
E. C
. (
2018
).
Construction of a neonatal cortical surface atlas using Multimodal Surface Matching in the Developing Human Connectome Project
.
NeuroImage
,
179
,
11
29
. https://doi.org/10.1016/j.neuroimage.2018.06.018
Brodoehl
,
S.
,
Gaser
,
C.
,
Dahnke
,
R.
,
Witte
,
O. W.
, &
Klingner
,
C. M.
(
2020
).
Surface-based analysis increases the specificity of cortical activation patterns and connectivity results
.
Scientific Reports
,
10
(
1
),
5737
. https://doi.org/10.1038/s41598-020-62832-z
Cao
,
M.
,
Shu
,
N.
,
Cao
,
Q.
,
Wang
,
Y.
, &
He
,
Y.
(
2014
).
Imaging functional and structural brain connectomics in attention-deficit/hyperactivity disorder
.
Molecular Neurobiology
,
50
(
3
),
1111
1123
. https://doi.org/10.1007/s12035-014-8685-x
Coalson
,
T. S.
,
Van Essen
,
D. C.
, &
Glasser
,
M. F.
(
2018
).
The impact of traditional neuroimaging methods on the spatial localization of cortical areas
.
Proceedings of the National Academy of Sciences
,
115
(
27
),
E6356
E6365
. https://doi.org/10.1073/pnas.1801582115
Denisova
,
K.
(
2019
).
Neurobiology, not artifacts: Challenges and guidelines for imaging the high risk infant
.
NeuroImage
,
185
,
624
640
. https://doi.org/10.1016/j.neuroimage.2018.07.023
DeSocio
,
J. E.
(
2018
).
Epigenetics, maternal prenatal psychosocial stress, and infant mental health
.
Archives of Psychiatric Nursing
,
32
(
6
),
901
906
. https://doi.org/10.1016/j.apnu.2018.09.001
Doria
,
V.
,
Beckmann
,
C. F.
,
Arichi
,
T.
,
Merchant
,
N.
,
Groppo
,
M.
,
Turkheimer
,
F. E.
,
Counsell
,
S. J.
,
Murgasova
,
M.
,
Aljabar
,
P.
,
Nunes
,
R. G.
,
Larkman
,
D. J.
,
Rees
,
G.
, &
Edwards
,
A. D.
(
2010
).
Emergence of resting state networks in the preterm human brain
.
Proceedings of the National Academy of Sciences
,
107
(
46
),
20015
20020
. https://doi.org/10.1073/pnas.1007921107
Edwards
,
A. D.
,
Rueckert
,
D.
,
Smith
,
S. M.
,
Abo Seada
,
S.
,
Alansary
,
A.
,
Almalbis
,
J.
,
Allsop
,
J.
,
Andersson
,
J.
,
Arichi
,
T.
,
Arulkumaran
,
S.
,
Bastiani
,
M.
,
Batalle
,
D.
,
Baxter
,
L.
,
Bozek
,
J.
,
Braithwaite
,
E.
,
Brandon
,
J.
,
Carney
,
O.
,
Chew
,
A.
,
Christiaens
,
D.
, …
Hajnal
,
J. V.
(
2022
).
The Developing Human Connectome Project neonatal data release
.
Frontiers in Neuroscience
,
16
,
886772
. https://doi.org/10.3389/fnins.2022.886772
Eyre
,
M.
,
Fitzgibbon
,
S. P.
,
Ciarrusta
,
J.
,
Cordero-Grande
,
L.
,
Price
,
A. N.
,
Poppe
,
T.
,
Schuh
,
A.
,
Hughes
,
E.
,
O’Keeffe
,
C.
,
Brandon
,
J.
,
Cromb
,
D.
,
Vecchiato
,
K.
,
Andersson
,
J.
,
Duff
,
E. P.
,
Counsell
,
S. J.
,
Smith
,
S. M.
,
Rueckert
,
D.
,
Hajnal
,
J. V.
,
Arichi
,
T.
, …
Edwards
,
A. D.
(
2021
).
The Developing Human Connectome Project: Typical and disrupted perinatal functional connectivity
.
Brain
,
144
(
7
),
2199
2213
. https://doi.org/10.1093/brain/awab118
Fitzgibbon
,
S. P.
,
Harrison
,
S. J.
,
Jenkinson
,
M.
,
Baxter
,
L.
,
Robinson
,
E. C.
,
Bastiani
,
M.
,
Bozek
,
J.
,
Karolis
,
V.
,
Cordero Grande
,
L.
,
Price
,
A. N.
,
Hughes
,
E.
,
Makropoulos
,
A.
,
Passerat-Palmbach
,
J.
,
Schuh
,
A.
,
Gao
,
J.
,
Farahibozorg
,
S.-R.
,
O’Muircheartaigh
,
J.
,
Ciarrusta
,
J.
,
O’Keeffe
,
C.
, …
Andersson
,
J.
(
2020
).
The developing Human Connectome Project (dHCP) automated resting-state functional processing framework for newborn infants
.
NeuroImage
,
223
,
117303
. https://doi.org/10.1016/j.neuroimage.2020.117303
Gao
,
W.
,
Alcauter
,
S.
,
Smith
,
J. K.
,
Gilmore
,
J. H.
, &
Lin
,
W.
(
2015
).
Development of human brain cortical network architecture during infancy
.
Brain Structure and Function
,
220
(
2
),
1173
1186
. https://doi.org/10.1007/s00429-014-0710-3
Glasser
,
M. F.
,
Coalson
,
T. S.
,
Robinson
,
E. C.
,
Hacker
,
C. D.
,
Harwell
,
J.
,
Yacoub
,
E.
,
Ugurbil
,
K.
,
Andersson
,
J.
,
Beckmann
,
C. F.
,
Jenkinson
,
M.
,
Smith
,
S. M.
, &
Van Essen
,
D. C
. (
2016
).
A multi-modal parcellation of human cerebral cortex
.
Nature
,
536
(
7615
),
171
178
. https://doi.org/10.1038/nature18933
Glasser
,
M. F.
,
Smith
,
S. M.
,
Marcus
,
D. S.
,
Andersson
,
J. L. R.
,
Auerbach
,
E. J.
,
Behrens
,
T. E. J.
,
Coalson
,
T. S.
,
Harms
,
M. P.
,
Jenkinson
,
M.
,
Moeller
,
S.
,
Robinson
,
E. C.
,
Sotiropoulos
,
S. N.
,
Xu
,
J.
,
Yacoub
,
E.
,
Ugurbil
,
K.
, &
Van Essen
,
D. C.
(
2016
).
The Human Connectome Project’s neuroimaging approach
.
Nature Neuroscience
,
19
(
9
),
1175
1187
. https://doi.org/10.1038/nn.4361
Glasser
,
M. F.
,
Sotiropoulos
,
S. N.
,
Wilson
,
J. A.
,
Coalson
,
T. S.
,
Fischl
,
B.
,
Andersson
,
J. L.
,
Xu
,
J.
,
Jbabdi
,
S.
,
Webster
,
M.
,
Polimeni
,
J. R.
,
Van Essen
,
D. C.
, &
Jenkinson
,
M.
(
2013
).
The minimal preprocessing pipelines for the Human Connectome Project
.
NeuroImage
,
80
,
105
124
. https://doi.org/10.1016/j.neuroimage.2013.04.127
Gordon
,
E. M.
,
Chauvin
,
R. J.
,
Van
,
A. N.
,
Rajesh
,
A.
,
Nielsen
,
A.
,
Newbold
,
D. J.
,
Lynch
,
C. J.
,
Seider
,
N. A.
,
Krimmel
,
S. R.
,
Scheidter
,
K. M.
,
Monk
,
J.
,
Miller
,
R. L.
,
Metoki
,
A.
,
Montez
,
D. F.
,
Zheng
,
A.
,
Elbau
,
I.
,
Madison
,
T.
,
Nishino
,
T.
,
Myers
,
M. J.
, …
Dosenbach
,
N. U.
(
2023
).
A somato-cognitive action network alternates with effector regions in motor cortex
.
Nature
,
617
,
351
359
. https://doi.org/10.1038/s41586-023-05964-2
Gordon
,
E. M.
,
Laumann
,
T. O.
,
Adeyemo
,
B.
,
Huckins
,
J. F.
,
Kelley
,
W. M.
, &
Petersen
,
S. E.
(
2016
).
Generation and evaluation of a cortical area parcellation from resting-state correlations
.
Cerebral Cortex
,
26
(
1
),
288
303
. https://doi.org/10.1093/cercor/bhu239
Gordon
,
E. M.
,
Laumann
,
T. O.
,
Gilmore
,
A. W.
,
Newbold
,
D. J.
,
Greene
,
D. J.
,
Berg
,
J. J.
,
Ortega
,
M.
,
Hoyt-Drazen
,
C.
,
Gratton
,
C.
,
Sun
,
H.
,
Hampton
,
J. M.
,
Coalson
,
R. S.
,
Nguyen
,
A. L.
,
McDermott
,
K. B.
,
Shimony
,
J. S.
,
Snyder
,
A. Z.
,
Schlaggar
,
B. L.
,
Petersen
,
S. E.
,
Nelson
,
S. M.
, &
Dosenbach
,
N. U.
(
2017
).
Precision functional mapping of individual human brains
.
Neuron
,
95
(
4
),
791
807.e7
. https://doi.org/10.1016/j.neuron.2017.07.011
Gratton
,
C.
,
Kraus
,
B. T.
,
Greene
,
D. J.
,
Gordon
,
E. M.
,
Laumann
,
T. O.
,
Nelson
,
S. M.
,
Dosenbach
,
N. U.
, &
Petersen
,
S. E.
(
2020
).
Defining individual-specific functional neuroanatomy for precision psychiatry
.
Biological Psychiatry
,
88
(
1
),
28
39
. https://doi.org/10.1016/j.biopsych.2019.10.026
Gratton
,
C.
,
Nelson
,
S. M.
, &
Gordon
,
E. M.
(
2022
).
Brain-behavior correlations: Two paths toward reliability
.
Neuron
,
110
(
9
),
1446
1449
. https://doi.org/10.1016/j.neuron.2022.04.018
Griffanti
,
L.
,
Salimi-Khorshidi
,
G.
,
Beckmann
,
C. F.
,
Auerbach
,
E. J.
,
Douaud
,
G.
,
Sexton
,
C. E.
,
Zsoldos
,
E.
,
Ebmeier
,
K. P.
,
Filippini
,
N.
,
Mackay
,
C. E.
,
Moeller
,
S.
,
Xu
,
J.
,
Yacoub
,
E.
,
Baselli
,
G.
,
Ugurbil
,
K.
,
Miller
,
K. L.
, &
Smith
,
S. M.
(
2014
).
ICA-based artefact removal and accelerated fMRI acquisition for improved resting state network imaging
.
NeuroImage
,
95
,
232
247
. https://doi.org/10.1016/j.neuroimage.2014.03.034
He
,
L.
,
Li
,
H.
,
Holland
,
S. K.
,
Yuan
,
W.
,
Altaye
,
M.
, &
Parikh
,
N. A.
(
2018
).
Early prediction of cognitive deficits in very preterm infants using functional connectome data in an artificial neural network framework
.
NeuroImage: Clinical
,
18
,
290
297
. https://doi.org/10.1016/j.nicl.2018.01.032
Hermosillo
,
R. J. M.
,
Moore
,
L. A.
,
Feczko
,
E.
,
Miranda-Domínguez
,
Ó.
,
Pines
,
A.
,
Dworetsky
,
A.
,
Conan
,
G.
,
Mooney
,
M. A.
,
Randolph
,
A.
,
Graham
,
A.
,
Adeyemo
,
B.
,
Earl
,
E.
,
Perrone
,
A.
,
Carrasco
,
C. M.
,
Uriarte-Lopez
,
J.
,
Snider
,
K.
,
Doyle
,
O.
,
Cordova
, M.
, …
Fair
,
D. A.
(
2024
).
A precision functional atlas of personalized network topography and probabilities
.
Nature Neuroscience
,
27
,
1000
1013
. https://doi.org/10.1038/s41593-024-01596-5
Hu
,
D.
,
Wang
,
F.
,
Zhang
,
H.
,
Wu
,
Z.
,
Zhou
,
Z.
,
Li
,
G.
,
Wang
,
L.
,
Lin
,
W.
,
Li
,
G.
, &
UNC/UMN Baby Connectome Project
Consortium
. (
2022
).
Existence of functional connectome fingerprint during infancy and its stability over months
.
The Journal of Neuroscience
,
42
(
3
),
377
389
. https://doi.org/10.1523/JNEUROSCI.0480-21.2021
Keunen
,
K.
,
Counsell
,
S. J.
, &
Benders
,
M. J.
(
2017
).
The emergence of functional architecture during early brain development
.
NeuroImage
,
160
,
2
14
. https://doi.org/10.1016/j.neuroimage.2017.01.047
Kong
,
R.
,
Li
,
J.
,
Orban
,
C.
,
Sabuncu
,
M. R.
,
Liu
,
H.
,
Schaefer
,
A.
,
Sun
,
N.
,
Zuo
,
X.-N.
,
Holmes
,
A. J.
,
Eickhoff
,
S. B.
, &
Yeo
,
B. T. T.
(
2019
).
Spatial topography of individual-specific cortical networks predicts human cognition, personality, and emotion
.
Cerebral Cortex
,
29
(
6
),
2533
2551
. https://doi.org/10.1093/cercor/bhy123
Kruschke
,
J. K.
(
2011
).
Doing Bayesian data analysis: A tutorial with R and BUGS
.
Academic Press OCLC
: ocn653121532. https://doi.org/10.1111/jedm.12029
Laumann
,
T. O.
,
Gordon
,
E. M.
,
Adeyemo
,
B.
,
Snyder
,
A. Z.
,
Joo
,
S. J.
,
Chen
,
M.-Y.
,
Gilmore
,
A. W.
,
McDermott
,
K. B.
,
Nelson
,
S. M.
,
Dosenbach
,
N. U.
,
Schlaggar
,
B. L.
,
Mumford
,
J. A.
,
Poldrack
,
R. A.
, &
Petersen
,
S. E.
(
2015
).
Functional system and areal organization of a highly sampled individual human brain
.
Neuron
,
87
(
3
),
657
670
. https://doi.org/10.1016/j.neuron.2015.06.037
Levi
,
P. T.
,
Chopra
,
S.
,
Pang
,
J. C.
,
Holmes
,
A.
,
Gajwani
,
M.
,
Sassenberg
,
T. A.
,
DeYoung
,
C. G.
, &
Fornito
,
A.
(
2023
).
The effect of using group-averaged or individualized brain parcellations when investigating connectome dysfunction in psychosis
.
Network Neuroscience
,
7
(
4
),
1228
1247
. https://doi.org/10.1162/netn_a_00329
Makropoulos
,
A.
,
Gousias
,
I. S.
,
Ledig
,
C.
,
Aljabar
,
P.
,
Serag
,
A.
,
Hajnal
,
J. V.
,
Edwards
,
A. D.
,
Counsell
,
S. J.
, &
Rueckert
,
D.
(
2014
).
Automatic whole brain MRI segmentation of the developing neonatal brain
.
IEEE Transactions on Medical Imaging
,
33
(
9
),
1818
1831
. https://doi.org/10.1109/TMI.2014.2322280
Makropoulos
,
A.
,
Robinson
,
E. C.
,
Schuh
,
A.
,
Wright
,
R.
,
Fitzgibbon
,
S.
,
Bozek
,
J.
,
Counsell
,
S. J.
,
Steinweg
,
J.
,
Vecchiato
,
K.
,
Passerat-Palmbach
,
J.
,
Lenz
,
G.
,
Mortari
,
F.
,
Tenev
,
T.
,
Duff
,
E. P.
,
Bastiani
,
M.
,
Cordero-Grande
,
L.
,
Hughes
,
E.
,
Tusor
,
N.
,
Tournier
,
J.-D.
, …
Rueckert
,
D.
(
2018
).
The developing human connectome project: A minimal processing pipeline for neonatal cortical surface reconstruction
.
NeuroImage
,
173
,
88
112
. https://doi.org/10.1016/j.neuroimage.2018.01.054
Mejia
,
A. F.
,
Nebel
,
M. B.
,
Wang
,
Y.
,
Caffo
,
B. S.
, &
Guo
,
Y.
(
2020
).
Template independent component analysis: Targeted and reliable estimation of subject-level brain networks using big data population priors
.
Journal of the American Statistical Association
,
115
(
531
),
1151
1177
. https://doi.org/10.1080/01621459.2019.1679638
Merhar
,
S. L.
,
Jiang
,
W.
,
Parikh
,
N. A.
,
Yin
,
W.
,
Zhou
,
Z.
,
Tkach
,
J. A.
,
Wang
,
L.
,
Kline-Fath
,
B. M.
,
He
,
L.
,
Braimah
,
A.
,
Vannest
,
J.
, &
Lin
,
W.
(
2021
).
Effects of prenatal opioid exposure on functional networks in infancy
.
Developmental Cognitive Neuroscience
,
51
,
100996
. https://doi.org/10.1016/j.dcn.2021.100996
Michon
,
K. J.
,
Khammash
,
D.
,
Simmonite
,
M.
,
Hamlin
,
A. M.
, &
Polk
,
T. A.
(
2022
).
Person-specific and precision neuroimaging: Current methods and future directions
.
NeuroImage
,
263
,
119589
. https://doi.org/10.1016/j.neuroimage.2022.119589
Mitra
,
A.
,
Snyder
,
A. Z.
,
Tagliazucchi
,
E.
,
Laufs
,
H.
,
Elison
,
J.
,
Emerson
,
R. W.
,
Shen
,
M. D.
,
Wolff
,
J. J.
,
Botteron
,
K. N.
,
Dager
,
S.
,
Estes
,
A. M.
,
Evans
,
A.
,
Gerig
,
G.
,
Hazlett
,
H. C.
,
Paterson
,
S. J.
,
Schultz
,
R. T.
,
Styner
,
M. A.
,
Zwaigenbaum
,
L.
,
The IBIS Network
,…
Raichle
,
M.
(
2017
).
Resting-state fMRI in sleeping infants more closely resembles adult sleep than adult wakefulness
.
PLoS One
,
12
(
11
),
e0188122
. https://doi.org/10.1371/journal.pone.0188122
Molloy
,
M. F.
, &
Saygin
,
Z. M.
(
2022
).
Individual variability in functional organization of the neonatal brain
.
NeuroImage
,
253
,
119101
. https://doi.org/10.1016/j.neuroimage.2022.119101
Nickerson
,
L. D.
(
2018
).
Replication of resting state-task network correspondence and novel findings on brain network activation during task fMRI in the Human Connectome Project study
.
Scientific Reports
,
8
(
1
),
17543
. https://doi.org/10.1038/s41598-018-35209-6
Nickerson
,
L. D.
,
Smith
,
S. M.
,
Öngür
,
D.
, &
Beckmann
,
C. F.
(
2017
).
Using dual regression to investigate network shape and amplitude in functional connectivity analyses
.
Frontiers in Neuroscience
,
11
,
115
. https://doi.org/10.3389/fnins.2017.00115
Nielsen
,
A. N.
,
Kaplan
,
S.
,
Meyer
,
D.
,
Alexopoulos
,
D.
,
Kenley
,
J. K.
,
Smyser
,
T. A.
,
Wakschlag
,
L. S.
,
Norton
,
E. S.
,
Raghuraman
,
N.
,
Warner
,
B. B.
,
Shimony
,
J. S.
,
Luby
,
J. L.
,
Neil
,
J. J.
,
Petersen
,
S. E.
,
Barch
,
D. M.
,
Rogers
,
C. E.
,
Sylvester
,
C. M.
, &
Smyser
,
C. D.
(
2022
).
Maturation of large-scale brain systems over the first month of life
.
Cerebral Cortex
,
33
(
6
),
2788
2803
. https://doi.org/10.1093/cercor/bhac242
Pham
,
D. D.
,
Muschelli
,
J.
, &
Mejia
,
A. F.
(
2022
).
ciftiTools: A package for reading, writing, visualizing, and manipulating CIFTI files in R
.
NeuroImage
,
250
,
118877
. https://doi.org/10.1016/j.neuroimage.2022.118877
Rajasilta
,
O.
,
Tuulari
,
J. J.
,
Björnsdotter
,
M.
,
Scheinin
,
N. M.
,
Lehtola
,
S. J.
,
Saunavaara
,
J.
,
Häkkinen
,
S.
,
Merisaari
,
H.
,
Parkkola
,
R.
,
Lähdesmäki
,
T.
,
Karlsson
,
L.
, &
Karlsson
,
H.
(
2020
).
Resting-state networks of the neonate brain identified using independent component analysis
.
Developmental Neurobiology
,
80
(
3-4
),
111
125
. https://doi.org/10.1002/dneu.22742
Risk
,
B. B.
,
Murden
,
R. J.
,
Wu
,
J.
,
Nebel
,
M. B.
,
Venkataraman
,
A.
,
Zhang
,
Z.
, &
Qiu
,
D.
(
2021
).
Which multiband factor should you choose for your resting-state fMRI study?
NeuroImage
,
234
,
117965
. https://doi.org/10.1016/j.neuroimage.2021.117965
Robinson
,
E. C.
,
Garcia
,
K.
,
Glasser
,
M. F.
,
Chen
,
Z.
,
Coalson
,
T. S.
,
Makropoulos
,
A.
,
Bozek
,
J.
,
Wright
,
R.
,
Schuh
,
A.
,
Webster
,
M.
,
Hutter
,
J.
,
Price
,
A.
,
Cordero Grande
,
L.
,
Hughes
,
E.
,
Tusor
,
N.
,
Bayly
,
P. V.
,
Van Essen
,
D. C.
,
Smith
,
S. M.
,
Edwards
,
A. D.
, …
Rueckert
,
D.
(
2018
).
Multimodal surface matching with higher-order smoothness constraints
.
NeuroImage
,
167
,
453
465
. https://doi.org/10.1016/j.neuroimage.2017.10.037
Robinson
,
E. C.
,
Jbabdi
,
S.
,
Glasser
,
M. F.
,
Andersson
,
J.
,
Burgess
,
G. C.
,
Harms
,
M. P.
,
Smith
,
S. M.
,
Van Essen
,
D. C.
, &
Jenkinson
,
M.
(
2014
).
MSM: A new flexible framework for multimodal surface matching
.
NeuroImage
,
100
,
414
426
. https://doi.org/10.1016/j.neuroimage.2014.05.069
Salimi-Khorshidi
,
G.
,
Douaud
,
G.
,
Beckmann
,
C. F.
,
Glasser
,
M. F.
,
Griffanti
,
L.
, &
Smith
,
S. M.
(
2014
).
Automatic denoising of functional MRI data: Combining independent component analysis and hierarchical fusion of classifiers
.
NeuroImage
,
90
,
449
468
. https://doi.org/10.1016/j.neuroimage.2013.11.046
Sandman
,
C. A.
,
Davis
,
E. P.
,
Buss
,
C.
, &
Glynn
,
L. M.
(
2012
).
Exposure to prenatal psychobiological stress exerts programming influences on the mother and her fetus
.
Neuroendocrinology
,
95
(
1
),
8
21
. https://doi.org/10.1159/000327017
Smith
,
S. M.
,
Fox
,
P. T.
,
Miller
,
K. L.
,
Glahn
,
D. C.
,
Fox
,
P. M.
,
Mackay
,
C. E.
,
Filippini
,
N.
,
Watkins
,
K. E.
,
Toro
,
R.
,
Laird
,
A. R.
, &
Beckmann
,
C. F.
(
2009
).
Correspondence of the brain’s functional architecture during activation and rest
.
Proceedings of the National Academy of Sciences
,
106
(
31
),
13040
13045
. https://doi.org/10.1073/pnas.0905267106
Smyser
,
C. D.
,
Inder
,
T. E.
,
Shimony
,
J. S.
,
Hill
,
J. E.
,
Degnan
,
A. J.
,
Snyder
,
A. Z.
, &
Neil
,
J. J.
(
2010
).
Longitudinal analysis of neural network development in preterm infants
.
Cerebral Cortex
,
20
(
12
),
2852
2862
. https://doi.org/10.1093/cercor/bhq035
Smyser
,
C. D.
, &
Neil
,
J. J.
(
2015
).
Use of resting-state functional MRI to study brain development and injury in neonates
.
Seminars in Perinatology
,
39
(
2
),
130
140
. https://doi.org/10.1053/j.semperi.2015.01.006
Thomas Yeo
,
B. T.
,
Krienen
,
F. M.
,
Sepulcre
,
J.
,
Sabuncu
,
M. R.
,
Lashkari
,
D.
,
Hollinshead
,
M.
,
Roffman
,
J. L.
,
Smoller
,
J. W.
,
Zöllei
,
L.
,
Polimeni
,
J. R.
,
Fischl
,
B.
,
Liu
,
H.
, &
Buckner
,
R. L.
(
2011
).
The organization of the human cerebral cortex estimated by intrinsic functional connectivity
.
Journal of Neurophysiology
,
106
(
3
),
1125
1165
. https://doi.org/10.1152/jn.00338.2011
Toulmin
,
H.
,
Beckmann
,
C. F.
,
O’Muircheartaigh
,
J.
,
Ball
,
G.
,
Nongena
,
P.
,
Makropoulos
,
A.
,
Ederies
,
A.
,
Counsell
,
S. J.
,
Kennea
,
N.
,
Arichi
,
T.
,
Tusor
,
N.
,
Rutherford
,
M. A.
,
Azzopardi
,
D.
,
Gonzalez-Cinca
,
N.
,
Hajnal
,
J. V.
, &
Edwards
,
A. D.
(
2015
).
Specialization and integration of functional thalamocortical connectivity in the human infant
.
Proceedings of the National Academy of Sciences
,
112
(
20
),
6485
6490
. https://doi.org/10.1073/pnas.1422638112
Wang
,
F.
,
Zhang
,
H.
,
Wu
,
Z.
,
Hu
,
D.
,
Zhou
,
Z.
,
Girault
,
J. B.
,
Wang
,
L.
,
Lin
,
W.
, &
Li
,
G.
(
2023
).
Fine-grained functional parcellation maps of the infant cerebral cortex
.
eLife
,
12
,
e75401
. https://doi.org/10.7554/eLife.75401
Williams
,
L. Z. J.
,
Fitzgibbon
,
S. P.
,
Bozek
,
J.
,
Winkler
,
A. M.
,
Dimitrova
,
R.
,
Poppe
,
T.
,
Schuh
,
A.
,
Makropoulos
,
A.
,
Cupitt
,
J.
,
O’Muircheartaigh
,
J.
,
Duff
,
E. P.
,
Cordero-Grande
,
L.
,
Price
,
A. N.
,
Hajnal
,
J. V.
,
Rueckert
,
D.
,
Smith
,
S. M.
,
Edwards
,
A. D.
, &
Robinson
,
E. C.
(
2023
).
Structural and functional asymmetry of the neonatal cerebral cortex
.
Nature Human Behaviour
,
7
(
6
),
942
955
. https://doi.org/10.1038/s41562-023-01542-8
Yu
,
X.
,
Ferradal
,
S.
,
Dunstan
,
J.
,
Carruthers
,
C.
,
Sanfilippo
,
J.
,
Zuk
,
J.
,
Zöllei
,
L.
,
Gagoski
,
B.
,
Ou
,
Y.
,
Grant
,
P. E.
, &
Gaab
,
N.
(
2022
).
Patterns of neural functional connectivity in infants at familial risk of developmental dyslexia
.
JAMA Network Open
,
5
(
10
),
e2236102
. https://doi.org/10.1001/jamanetworkopen.2022.36102
This is an open-access article distributed under the terms of the Creative Commons Attribution 4.0 International (CC BY 4.0) license, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. For a full description of the license, please visit https://creativecommons.org/licenses/by/4.0/legalcode.

Supplementary data