Data acquisition and pre-processing

Data were acquired on a 3T Philips Achieva with a dedicated neonatal imaging system including a neonatal 32-channel phased-array head coil. Fifteen minutes of high temporal and spatial resolution rs-fMRI data were acquired using a multislice gradient-echo EPI sequence with multiband excitation (TE = 38 ms; TR = 392 ms; MB factor = 9×; 2.15 mm isotropic, 2300 volumes). In addition, single-band EPI reference (sbref) scans were also acquired with bandwidth-matched readout, along with additional spin-echo EPI acquisitions with 4×AP and 4×PA phase-encoding directions. To correct susceptibility distortion in rs-fMRI data, field maps were also obtained from an interleaved (dual TE) spoiled gradient-echo sequence (TR = 10 ms; TE1 = 4.6 ms; TE2 = 6.9 ms; flip angle = 10°; 3 mm isotropic in-plane resolution, 6 mm slice thickness). High-resolution T1- and T2-weighted anatomical imaging were also acquired in the same scan session, with a spatial resolution of 0.8 mm isotropic. For T1w image: TR = 4795 ms and the field of view (FOV) = 145 × 122 × 100 mm. For T2w image: TR = 12 000 ms, TE = 156 ms and the FOV = 145 × 122 × 100 mm.

The dHCP rs-fMRI data were pre-processed by dHCP group using the project’s in-house pipeline optimized for neonatal imaging. See SI and Fitzgibbon et al.⁵⁶ for full details. In order to reduce signal artefacts related to head motion, the cardiorespiratory fluctuations and multiband acquisition, the 24 extended rigid-body motion parameters together with single-subject independent component analysis (ICA) noise components were regressed out. To further reduce the effect of motion on functional connectivity (FC) measures, motion-outlier volumes were identified, and a scrubbing procedure was applied to retain a continuous sub-sample of the data (∼70%) with the lowest motion for each participant. The subjects who still had a high level of motion after scrubbing procedure were excluded from further analyses. We discarded the first 5 volumes to allow for adaptation to the environment and equilibrium of the MR signal at first. Then, motion outliers were identified from the remaining 2295 volumes. Volumes with the root mean square intensity difference between successive volumes higher than 1.5 interquartile range above the 75th centile, after motion and distortion correction, were considered motion outliers. Then, a continuous sub-sample of 1600 volumes with the minimum number of motion outliers was retained for each subject. Subjects with more than 160 motion-outlier volumes (10% of the cropped dataset) in the continuous subset were labelled ‘high level of motion’ and excluded entirely. Thus, 8 preterm neonates scanned before TEA, 14 preterm neonates scanned at TEA and 61 full-term neonates were excluded. In addition, we performed a temporal low-pass filter (0.08 Hz low-pass cutoff) on the pre-processed dHCP rs-fMRI to conduct FC analyses, as previous studies found that oscillations were primarily detected within grey matter in 0.01–0.08 Hz.^57,58Supplementary Fig. 1A provides a schematic of the processing steps for dHCP fMRI data.

Data were acquired on a customized 3T Siemens ‘Connectome Skyra’ with a 32-channel head coil. Resting state images were collected using gradient-echo EPI sequence: TR = 720 ms; TE = 33.1 ms; flip angle = 52°; FOV = 208 × 180 mm (RO × PE), slice thickness = 2 mm, 72 slices, 2.0 mm isotropic voxels, 1200 volumes per run. rs-fMRI data were pre-processed by HCP group. See SI and Van Essen et al.⁵⁵ for full details.