Animal studies of neurodevelopmental plasticity have shown that intrinsic brain activity evolves from high amplitude and globally synchronized to suppressed and sparse as plasticity declines and the cortex matures. Leveraging functional MRI data from 1033 individuals (8-23 years), we reveal that this stereotyped refinement of intrinsic activity occurs during human development and provides evidence for a cortical gradient of neurodevelopmental plasticity during childhood and adolescence. Specifically, we demonstrate that declines in the amplitude of intrinsic activity are initiated heterochronously across regions, coupled to the maturation of a plasticity-restricting structural feature, and temporally staggered along a hierarchical sensorimotor-association axis from ages 8 to 18. Youth from disadvantaged environments exhibit reduced intrinsic activity in regions further up the sensorimotor-association axis, suggestive of a reduced potential for plasticity in late-maturing cortices. Our results uncover a hierarchical axis of neurodevelopment and offer insight into the temporal sequence of protracted neurodevelopmental plasticity in humans. Continued discovery of temporal axes of development across human’s multi-decade maturational course will provide evidence of how plasticity is distributed across brain regions at different developmental stages. Such insights into the temporal patterning of plasticity may help to guide interventions in youth that align with each child’s neurotemporal context.
Valerie J. Sydnor
Theodore D. Satterthwaite
Bart Larsen, Azeez Adebimpe, Maxwell A. Bertolero, Matthew Cieslak, Sydney Covitz, Yong Fan, Raquel E. Gur, Ruben C. Gur, David R. Roalf, Russell T. Shinohara, Dani S. Bassett, Theodore D. Satterthwaite
RBC PNC-Health Exclude (primary) and LTN (sensitivity)
CBF/: parcel-wise cerebral blood flow maps for each participant, generated with ASLPrep code/: directory with the spatiotemp_dev_plasticity github repo FluctuationAmplitude/PNC/: vertex-wise and parcel-wise fluctuation amplitude maps for each participant, generated with xcp-d and connectome workbench FluctuationAmplitude/GAMRESULTS/: gam model outputs (effect sizes, p-values, fitted values, smooth estimates, smooth characteristics, derivatives) Maps/: surface parcellation files and SNR masks (Maps/parcellations/) and S-A axis github repo (Maps/S-A_ArchetypalAxis/) Myelin/: myelin development maps including the age effect size map (r2), the age of maximal growth map (age of max slope), and the annualized rate of change map (annualized roc) from Baum et al., 2021 sample_info/: sample demographics, factor scores, rbcid-bblid key, and final project participant list (PNC_FinalSample_N1033.csv) software/: project software Structural/: freesurfer output for each participant Timeseries/: vertex-wise and parcel-wise fully processed resting fMRI timeseries data for each participant, generated with fmriprep and xcp-d
The entire analytic workflow implemented in this project is described in the following sections and links to the corresponding github code are provided. This workflow includes quantification of regional fluctuation amplitude, PNC sample selection, fitting of generalized additive models (GAMs), and characterization of relationships between fluctuation amplitude, age, environmental variability, and the sensorimotor-association axis. Scripts were implemented in the order outlined below.
fmriprep 20.2.3 was run with the following parameters:
$ singularity run pennlinc-containers/.datalad/environments/fmriprep-20-2-3/image inputs/data prep participant --output-spaces MNI152NLin6Asym:res-2 --participant-label "$subid" --force-bbr --cifti-output 91k -v -v
xcp-d 0.0.4 was run with the following parameters:
$ singularity run pennlinc-containers/.datalad/environments/xcp-abcd-0-0-4/image inputs/data/fmriprep xcp participant --despike --lower-bpf 0.01 --upper-bpf 0.08 --participant_label $subid -p 36P -f 10 –cifti
Vertex-wise fluctuation amplitude maps were then parcellated with /fluctuation_amplitude/parcellate.Rmd to quantify mean fluctuation amplitude in each cortical region. Regions were defined with the HCP multimodal atlas (i.e. Glasser 360, primary analyses) and with the Schaefer 400 atlas (sensitivity analysis).
Fluctuation amplitude analyses were only conducted in brain regions that reliably exhibited high signal to noise ratio (SNR) in PNC functional MRI data. The vertex level SNR map generated in Cui et al., 2020, Neuron was parcellated with Glasser 360 and Schaefer 400 atlases with the script /fluctuation_amplitude/SNR_mask.Rmd for use in study analyses. Regions wherein >= 25% of vertices had attenuated signal were excluded from analyses.
The final study sample was constructed with /sample_construction/finalsample.Rmd. The final sample was generated from the 1374 PNC participants with dominant group ses-PNC1_task-rest_acq-singleband scans (non-variant CuBIDS acquisitions). The following exclusions were applied to generate the final sample of N = 1033 participants:
Health exclude: 120 participants with medical problems that could impact brain function or incidentally-encountered brain structure abnormalities were excluded from the sample
T1 QA exclude: 23 participants with T1-weighted scans that failed visual QC were excluded from the sample
fMRI motion exclude: 179 participants with a mean relative RMS motion value > 0.2 during the resting state fMRI scan were excluded from the sample
Fluctuation amplitude outlier exclude: from the remaining 1052 participants, 19 individuals that had outlier (+- 4 SD from the mean) fluctuation amplitude data in more than 5% of Glasser 360 parcels were excluded from the sample
To characterize age-dependent changes in spontaneous activity fluctuations across the developing cortex as well as associations between fluctuation amplitude and environmental factors, generalized additive models were fit in each cortical region. GAM models and associated statistics, fitted values, smooths, and derivatives were quantified with the set of functions included in /gam_models/GAM_functions.R. This script includes the following functions:
Model Fitting: Age-Dependent Changes in Regional Fluctuation Amplitude
Developmental models were fit with gam_models/fitGAMs_fluctuationamplitude_age.R, which calls the functions described above. Age-focused GAMs were implemented for all regions included in the Glasser 360 and Schaefer 400 atlases, using the final study sample of N = 1033 generated during the sample construction process.
Model Fitting: Developmental Environment-Dependent Variation in Regional Fluctuation Amplitude
Models examining associations between fluctuation amplitude and environmental and family characteristics were run via gam_models/fitGAMs_fluctuationamplitude_environment.R on the N = 1033 study sample.
Model results were examined and studied within our hierarchical neurodevelopmental plasticity framework in the markdown file /developmental_effects/hierarchical_development.Rmd. This markdown executes the following:
A rendered html of hierarchical_development.Rmd can be viewed here!
The robustness of our developmental findings was confirmed in a series of sensitivity analyses. Sensitiviy analysis GAMs were fit with sensitivity_analyses/fitGAMs_sensitivityanalyses.R and results were examined in sensitivity_analyses/sensitivityresults.Rmd. The following sensitivity analyses were performed: