NeuroImage 229 (2021) 117769

Contents lists available at ScienceDirect

NeuroImage
journal homepage: www.elsevier.com/locate/neuroimage

Development of the brain functional connectome follows

puberty-dependent nonlinear trajectories
Zeus Gracia-Tabuenca, Martha Beatriz Moreno, Fernando A. Barrios, Sarael Alcauter∗
Instituto de Neurobiología, Universidad Nacional Autónoma de México, Blvd Juriquilla 3001, Querétaro 76230, Mexico

a r t i c l e i n f o a b s t r a c t

Keywords: Adolescence is a developmental period that dramatically impacts body and behavior, with pubertal hormones
puberty playing an important role not only in the morphological changes in the body but also in brain structure and
adolescence function. Understanding brain development during adolescence has become a priority in neuroscience because
functional connectivity
it coincides with the onset of many psychiatric and behavioral disorders. However, little is known about how
resting state fMRI
puberty influences the brain functional connectome. In this study, taking a longitudinal human sample of typically
cognitive flexibility
brain network developing children and adolescents (of both sexes), we demonstrate that the development of the brain functional
connectome better fits pubertal status than chronological age. In particular, centrality, segregation, efficiency,
and integration of the brain functional connectome increase after the onset of the pubertal markers. We found
that these effects are stronger in attention and task control networks. Lastly, after controlling for this effect,
we showed that functional connectivity between these networks is related to better performance in cognitive
flexibility. This study points out the importance of considering longitudinal nonlinear trends when exploring
developmental trajectories, and emphasizes the impact of puberty on the functional organization of the brain in
adolescence.

1. Introduction The widespread use of in vivo neuroimaging techniques has shifted

the study of developmental trajectories from focal brain areas to a sys-
Adolescence is a developmental period that dramatically impacts tem perspective (Di Martino et al., 2014; Zuo et al., 2017). Thus, model-
body and behavior. Pubertal hormones play an important role in ing the brain as a set of interconnected elements, or connectome, allows
adrenal, gonadal, and growth axes, having a substantial impact not only the description of its functional organization in terms of network-based
in the morphological changes in the body but also in brain structure properties, such as the functional segregation and integration, and the
and function (Vijayakumar et al., 2018). Therefore, measuring brain efficiency of the brain network (Rubinov and Sporns, 2010; Stam and
development along this stage is complicated given that these changes Van Straaten, 2012). This framework has been extensively applied in
are not only influenced by chronological age (Blakemore et al., 2010). other lifespan studies, showing dramatic developmental changes in rela-
Other factors play a role, such as sex, ethnic, and environmental influ- tion with chronological age (Zuo et al., 2011; Tomasi and Volkow, 2012;
ences, and all of them with a high degree of variability between indi- Alcauter et al., 2014; Betzel et al., 2014; Gao et al., 2015; Gracia-
viduals (Sawyer et al., 2018). Also, adolescence is considered a critical Tabuenca et al., 2018). Particularly in adolescence, previous cross-
period for mental health, when many behavioral and psychiatric disor- sectional studies have shown an increase in functional segregation along
ders first manifest (Kessler et al., 2005; Paus et al., 2008). Additionally, with cognitive specialization (Fair et al., 2009; Satterthwaite et al.,
cognitive abilities as important as task-switching or cognitive flexibility 2013a; Gu et al., 2015), whereas other studies have emphasized an in-
strengthen in adolescence (Somsen, 2007; Hauser et al., 2015), together crease in the functional integration between brain networks with age
with other executive functions, including working memory, inhibition, (Hwang et al., 2012; Marek et al., 2015). Such potentially divergent re-
and attention (Gur et al., 2012). Such cognitive improvement, that usu- sults may be influenced by differences in samples and/or methodologies;
ally starts in late infancy and continues improving until early adulthood, however, the implementation of longitudinal samples and/or nonlinear
may result from the morphological and metabolic changes of the brain, trends will potentially reduce the bias when characterizing neurodevel-
which are expected to improve their efficiency to solve complex prob- opmental trajectories (Zuo et al., 2017; Vijayakumar et al., 2018). Addi-
lems (Baum et al., 2017; Chai et al., 2017). tionally, the motion artifact is another potential bias given that younger
subjects tend to move more in pediatric populations affecting the sam-
ple inferences (Satterthwaite et al., 2013a). Previous findings have de-
scribed mitigating strategies at the edgewise level (Ciric et al., 2017),
∗
Corresponding author.
E-mail address: [email protected] (S. Alcauter).

https://doi.org/10.1016/j.neuroimage.2021.117769
Received 26 October 2020; Received in revised form 18 December 2020; Accepted 5 January 2021
Available online 20 January 2021
1053-8119/© 2021 The Authors. Published by Elsevier Inc. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/)
Z. Gracia-Tabuenca, M.B. Moreno, F.A. Barrios et al. NeuroImage 229 (2021) 117769

sponses include: 1) absence, 2) first signs, 3) evident, and 4) complete.

This scale allows taking into account sex-specific patterns of puberty de-
velopment even at the same chronological age. Children under 10 years
old were set to level 1 of the PDS, as previous studies have shown simi-
lar values at these ages (Hibberd et al., 2015; van Duijvenvoorde et al.,
2019). Missing PDS values corresponding to 8 time-points (4 females)
were estimated with a Generalized Additive Mixed Model (GAMM) in-
cluding Locally Estimated Scatterplot Smoothing (LOESS) terms for the
age-sex interaction.

2.3. MRI sessions

MRI sessions consisted of a resting-state fMRI sequence followed

with a high-resolution T1-weighted scan. 150 whole-brain fMRI vol-
umes were acquired using a gradient recalled T2∗ echo-planar imaging
Fig. 1. Longitudinal sample. Each dot represents a subject at the age of as- sequence (TR/TE = 2000/40 ms, voxel size 4 × 4 × 4 mm3 ). Partici-
sessment. Lines represent longitudinal assessments of a subject. Abbreviations: pants were instructed to keep their eyes closed and not to fall asleep.
female (F), male (M). To facilitate participants to remain awake, this sequence was always
applied in the morning and at the beginning of the MRI session. Fur-
thermore, participants were asked to avoid staying up after their usual
however little is known regarding functional organization measures in bedtime. Structural T1 sequence was obtained using a 3D spoiled gradi-
a developmental context. ent recalled (SPGR) acquisition (TR/TE = 8.1/3.2 ms, flip angle = 12.0,
Furthermore, the pubertal spurt is an important source of variability voxel size 1 × 1 × 1 mm3 ). MR images were acquired with a 3T MR
in adolescence when considering developmental trajectories based on GE750 Discovery scanner (General Electric, Waukesha, WI), using an
chronological age, but it is rarely considered (Blakemore et al., 2010; 8-channel-array head coil. 20 sessions were acquired with a 32-channel
Vijayakumar et al., 2018). A pair of recent studies have taken an impor- coil, therefore a dummy covariate was included in the subsequent anal-
tant step in this direction, by fitting the developmental trajectories of yses.
fronto-subcortical connectivity with the pubertal stage (Spielberg et al.,
2014; van Duijvenvoorde et al., 2019). However, the relationship be- 2.4. Preprocessing
tween the pubertal stage and the developmental of the brain functional
organization still remains largely unexplored. Anatomical T1 images were denoised with non-local means
This work is focused on elucidating the aforementioned discrepan- (Manjón et al., 2010) and N4 bias field correction (Tustison et al.,
cies and unknowns. First, we tested how mitigation strategies of the mo- 2010). fMRI preprocessing was implemented using FMRIB’s Software
tion artifact may impact the functional connectome features. Then, we Libraries (FSL v.5.0.6; Jenkinson et al., 2012; RRID: SCR_002823) in-
tested the development of the functional organization of the adolescent cluding slice timing, head motion correction, brain extraction, inten-
brain, by exploring the longitudinal trajectories with generalized addi- sity normalization, regression of confounding variables, spatial normal-
tive mixed models, and their relation with pubertal stage and cognitive ization, and band-pass temporal filtering (0.01–0.08 Hz). Each fMRI
development. dataset was registered to its corresponding structural image with a rigid-
body transformation, followed by a two-step diffeomorphic SyN regis-
2. Methods tration (Avants et al., 2008; RRID: SCR_004757) to a pediatric brain
atlas (NIHPD 4.5–18.5; Fonov et al., 2011), and then to the MNI-152
2.1. Sample standard template.
Framewise displacement was estimated and summarized by the
The total sample consisted of 98 typically developing subjects (45 means of FD-RMS computed with FSL’s MCFLIRT (Jenkinson et al.,
male, ages from 6.71 to 18.06 years old), of whom 41 and 16 returned 2002), and those volumes with more than 0.25 mm were considered
to a second and a third session, respectively (Fig. 1). The average in- “spikes”, and they were regressed out from the fMRI time series. Given
terval for sessions 1–2 was 5.06 years (s.d. 1.24 years), and for sessions that pediatric samples tend to move more inside the scanner than
2–3 was 2.8 years (s.d. 0.23 years). A general invitation describing the adult populations, we implemented two rigorous confounding regres-
study protocol and inclusion criteria was sent to local schools. Inclusion sion strategies, one including the global signal regression (GSR) and the
criteria included full-term gestation (at least 37 weeks) and enrollment other one without GSR:
in its corresponding academic year. Exclusion criteria consisted of any
• GSR: 36 parameters were regressed out (Satterthwaite et al., 2013a),
neurological or psychiatric disorder identified with a semi-structured in-
including the time series, derivatives, and their corresponding
terview: MINI-KID for minors and MINI for those who were older than
quadratic terms of the six estimated motion parameters, and the
18 years old. For every session, written informed consent was required,
average signal from global and the eroded white matter (WM) and
and in the case of minors, parents were required to sign the consent
cerebrospinal fluid (CSF).
and children to verbally assent. Study protocols followed the Declara-
• no-GSR: 37 parameters were regressed out (Gracia-Tabuenca et al.,
tion of Helsinki and were approved by the Institutional Review Board.
2020), including the time series, derivatives, and their correspond-
Follow-ups occurred after 5 years and the second after 2 years.
ing quadratic terms of the six estimated motion parameters, and the
average signal from the eroded WM and CSF, plus five principal com-
2.2. Pubertal status
ponents of the WM+CSF, this last procedure is termed aCompCor
(Chai et al., 2012; Muschelli et al., 2014).
Participants completed the Pubertal Development Scale (PDS;
Petersen et al., 1988). PDS is a self-report questionnaire that averages Eighteen datasets from seventeen subjects (9/8 females/males; age:
five Likert-like scales to respond about growth spurt in height, pubic 6.97–13.16 years old) with less than four minutes of “spikes”-clean time
hair, and skin change for both sexes; plus facial hair growth and voice series were discarded for further analyses (Satterthwaite et al., 2013a),
change for males, and breast growth and menarche for females. Re- resulting in a final sample of 89 typically developing subjects (39 male,

2
Z. Gracia-Tabuenca, M.B. Moreno, F.A. Barrios et al. NeuroImage 229 (2021) 117769

age: 6.71 – 18.06 years old), of which 37 and 11 returned for a second density; significance was corrected with a false discovery rate (FDR q <
and third session, respectively. 0.05; Benjamini and Hochberg, 1995) for the thirteen functional net-
works defined in the atlas (Power et al., 2011).
2.5. Network features
2.7. Brain-cognition inference
Brain network construction was based on the 264 regions of inter-
est (ROIs) described by Power et al. (2011), which consist of 5-mm ra- A cross-sectional subset of the sample (N = 59; 24 males; age range:
dius spheres, encompassed in 13 different functional networks. This at- 9.67 - 17.98 years old) fulfilled a neuropsychological assessment that in-
las has been widely applied in pediatric studies (Satterthwaite et al., cluded a Card Sorting Task (CST; Berg, 1948; Lázaro et al., 2012). The
2013a, Satterthwaite et al., 2014; Gu et al., 2015; Marek et al., 2015; CST measures executive control, reasoning, and learning, but is particu-
Ciric et al., 2017; Chai et al., 2017; Gracia-Tabuenca et al., 2020). Edges larly sensitive to task shifting or cognitive flexibility. Briefly, a set of 64
were calculated as the functional connectivity between every pair of cards are given one at a time to the participant; each card can be catego-
ROIs, via the Pearson’s cross-correlation from their average time se- rized based on color, shape, or the number of figures. The participant has
ries. Based on these matrices, weighted graph theory metrics were com- to guess the correct category by putting the card behind four templates,
puted to characterize network organization features. The sum of edges and after that is told if the selection was correct or incorrect. After ten
or weighted degree (d), clustering coefficient (C), efficiency (E), and consecutive correct attempts, the category is shifted without informing
characteristic path length (L) were estimated to account for centrality, the participant, thus s/he has to guess again. Several scores are recorded:
segregation, and integration in terms of functional network organiza- the score of corrected attempts, perseverations (consecutive errors in
tion (Rubinov and Sporns, 2010). C was calculated based on the Barrat the same category), and recurrent perseverations (non-consecutive per-
formula (Barrat et al., 2004). L was calculated as the average shortest- severations in a three attempt span). The total score can be interpreted
paths based on the Dijkstra algorithm (Dijkstra, 1959), where pairwise as a measure of cognitive flexibility; meanwhile, perseverations reflect
distances were set as the inverse of the corresponding functional connec- impulsivity and a lack of executive control and task-switching.
tivity edge. E was computed as the average of the inverse of the shortest- The relationship between the task scores and the functional connec-
paths of the connectome. For L, the shortest-path of disconnected nodes tomes was assessed with the Network-Based Statistics approach (NBS;
was set to the maximum observed value of the network in order to avoid Zalesky et al., 2010) using the R-package Network-Based R-statistics
infinite outcomes (Fornito et al., 2010). Barrat and Dijkstra’s algorithms (NBR; Gracia-Tabuenca and Alcauter, 2020; RRID: SCR_019114). NBR
were calculated via the R package ‘igraph’ (Csardi and Nepusz, 2006; implements the NBS framework allowing to implement other models
RRID:SCR_019225). instead of the general linear model, in this case the edgewise relation-
Additionally, we tested the relationship of the head motion (i.e., FD- ship was tested via Generalized Additive Models (GAM) splines. This
RMS) and the brain network properties via linear mixed-effects models algorithm computes the likelihood of each set of connected edges that
(LME) for the two confounding regression strategies, in order to iden- surpass an initial statistical threshold (also referred to as clusters of con-
tify the best preprocessing strategy, as suggested by Ciric et al. (2017). nections, sub-networks, or components). Such likelihood estimation is
Besides, we tested if such relationships are affected by the Euclidean based on their size (sum of the weights of the edges in this case) and
distance of the network’s nodes (Satterthwaite et al., 2013b). how it compares to a null distribution generated by permutations of
the original data. Specifically, the null distribution is the density of the
2.6. Neurodevelopmental trajectories maximum size of the suprathreshold component for each permutation of
the original dataset. In this case, the three standardized cognitive scores
Six models were tested to delineate the neurodevelopmental ef- were fit for each intra- and inter-connection of a connectome of 13 × 13
fects. Two linear mixed-effect models (LME) for the age (1) and functional networks (Power et al., 2011). The model includes the GAM
the age-sex interaction (2), plus four generalized additive mixed spline of the residual degree at the 25% connectivity density. LOESS
models (GAMM) with smooth splines for the age (3), age by sex PDS, sex, head-motion, and coil terms were previously regressed out. A
(4), and with nonparametric Locally Estimated Scatterplot Smooth- null distribution was computed based on 1000 permutations setting p <
ing (LOESS) for the PDS (5) and the PDS-sex interaction (6). Every 0.05 for individual connections.
model included random-effects intercepts, head-coil, and average head-
motion as covariates. All models were implemented in R, using the 2.8. Code accessibility
libraries ‘lme4’ (Bates et al., 2014; RRID:SCR_015654) for the LME,
‘gamm4’ for the GAMM with smooth splines (Wood and Scheipl, 2014; The preprocessed data and the code/software described in the paper
https://CRAN.R-project.org/package=gamm4), and ‘gamlss’ for the is freely available online at https://github.com/BrainMapINB/Pubertal_
GAMM with LOESS (Stasinopoulos and Rigby 2007; https://CRAN.R- Connectome.
project.org/package=gamlss). Basis functions for the smooth splines
were set to four (k = 4) according to van Duijvenvoorde et al. (2019), 3. Results
meanwhile LOESS terms were set as ‘gamlss’ defaults (span = 0.75; de-
gree = 2). Mixed-effects were estimated via Maximum Log-likelihood 3.1. Head-motion mitigation strategies over functional connectivity
in order to compare between models. Model selection was based on
the Akaike Information Criterion (AIC). The AIC evals the goodness of We tested two confounding strategies to minimize the effect of head
fit of the model via the log-likelihood function (lnL) plus a penaliza- motion on the brain connectivity estimates; one included the global sig-
tion of the number of parameters of the model (k), being AIC=2k-2lnL nal regression (GSR) and the other did not (no-GSR). For the functional
(Akaike, 1974). Therefore, the model with lower AIC is the one with a connectivity density, both strategies showed a sharp distribution cen-
better trade-off of adjustment and parsimony. To evaluate the robustness tered on zero. On the other side, graph theory measures showed different
at different connectivity thresholds, network parameters were computed associations with head motion after implementing such strategies, but in
in a range of 1% to 48% of connectivity densities. general for all measures, using GSR resulted in lower associations with
At the 25% connectivity density, the terms of the model with lower head motion, as reflected by their distributions being centered closer to
AIC were tested via the Likelihood Ratio Test (LRT). That is, we com- 0 (Fig. 2A). Additionally, we tested if these confounding strategies were
puted the LRT of the full model and a null model dropping the term of affected by the Euclidean distance between the nodes for the head mo-
interest (Stasinopoulos and Rigby 2007). In addition, network proper- tion at the edgewise level, however, both approaches showed a pattern
ties were estimated for each functional network at the 25% connectivity with a slope near to zero (Fig. 2B).

3
Z. Gracia-Tabuenca, M.B. Moreno, F.A. Barrios et al. NeuroImage 229 (2021) 117769

When comparing the adjustment of the six different models for each
network measure at different connectivity densities, we found that for
every measure, the models including PDS (without the PDS-sex interac-
tion) showed the lowest AIC for the vast majority of connectivity densi-
ties (Fig. 4).
Taking a representative connectivity density of 25%, we can observe
the neurodevelopmental trajectories via GAMM LOESS, for the four net-
work measures the trend along the PDS was nonlinear, with an inflec-
tion point around the level 2 of the pubertal scale (Fig. 5). In particu-
lar, the degree, clustering coefficient, and efficiency showed a nonlinear
pattern, where these measures increase in value after the turning point,
while the path-length showed the inverse pattern (Fig. 5). The PDS term
was significant for each network measure, as well as the random-effect
term for the degree, clustering, and path-length (Table 1).
At the functional network level, neurodevelopmental effects along
the pubertal stage were found mainly in frontal and parietal related net-
works for every graph theory measure (Fig. 6). Particularly, for the de-
gree, stronger effects were found in the default mode (DMN), showing a
nonlinear pattern similar to the whole-brain trend (Table 2; Fig. 7). The
clustering coefficient showed a similar spatial pattern, but with no sig-
nificant effects after correcting for multiple comparisons. For efficiency,
several frontal and parietal related functional networks such as the de-
fault mode (DMN), dorsal attention (DAN), and frontoparietal (FPN),
plus the cerebellar (CBL), subcortical (SUB), and visual (VIS) networks
showed a strong effect along the PDS (Table 2; Fig. 7). While for the
path-length, as expected, the effects were inverse to those of efficiency
Fig. 2. A) Density plots of the linear mixed-effects (LME) t-scores of the rela- (Fig. 7), being the stronger effects in the cerebellar (CBL), default mode
tionship between the head motion (measured by the average FD-RMS) and the (DMN), frontoparietal (FPN), salience (SAL), ventral attention (VAN),
functional connectivity of edgewise and nodal graph theory variables, for both and visual (VIS) networks (Table 2).
preprocessing approaches: with and without Global Signal Regression (GSR).
B) Hexbin plots of the Euclidean distance between pairs of nodes and the as-
sociation of head motion with functional connectivity (LME t-scores), for both
preprocessing strategies: with and without GSR. Colormap expresses the count 3.3. Brain-cognition relationships
of edges with that effect. The fitted line of this relationship is depicted in red.
Regarding the relationship of the functional connectome and cogni-
tive performance, NBS identified a significant set of network interac-
tions that relate with performance in the Card Sorting Task (CST) as-
3.2. Neurodevelopmental trajectories sessment (Fig. 8). A widespread inter- and intra-network connectivity
was associated with higher scores in the CST (pFWE = 0.013), where
Regarding the pubertal scale along chronological age, both sex the cingulo-opercular network (CON), dorsal attention (DAN), default
groups showed a sigmoid shape with faster growth in females than mode (DMN), and fronto-parietal (FPN) networks are those with greater
males. There was a significant effect for the age (LRT = 349.58; number of positive connections. Additionally, two connections of the
DOF = 15.1; p < 2.2e-16), age-sex interaction (LRT = 52.81; sensorimotor-hand network (SMH) showed a concave pattern. No sig-
DOF = 14.99; p = 4.11e-06), and random-effects terms (LRT = 24.88; nificant components were found for the perseveration nor recursive per-
DOF = 10.56; p = 0.0076; Fig. 3). severation scores after the permutation test correction.

Fig. 3. Pubertal Development Scale scores (PDS) along age.

Thin lines represent individual trajectories; thick lines repre-
sent GAMM LOESS curves per sex group (with 95% confidence-
interval shadow). Abbreviations: female (F), male (M).

4
Z. Gracia-Tabuenca, M.B. Moreno, F.A. Barrios et al. NeuroImage 229 (2021) 117769

4. Discussion

By exploring a longitudinal sample and with a strict control of the ef-

fects of head motion, the present study showed that the brain functional
organization along adolescence follows nonlinear puberty-dependent
trajectories, with an inflection point around the beginning of puberty,
when the first signs of body changes appear (level 2 of the PDS). The
frontal and parietal functional networks exhibited the most dramatic
rates of change after the inflection point. Finally, a specific pattern of
network interactions reflected higher performance in a cognitive task
associated with executive functions, mainly cognitive flexibility. Such
a pattern includes the higher integration of attention and task-control
networks.

4.1. Head-motion effects on graph theory

Taking two rigorous confounding strategies to mitigate the effect of

head motion on the network estimates, considering or not considering
the global signal regression (GSR), we found that both strategies prop-
erly account for this effect for the edgewise functional connectivity, even
when considering the Euclidean distance between nodes. When explor-
ing such associations for the graph theory measures, they are closer to
zero for the approach using GSR (Fig. 2). In addition to the confounding
strategies at the subject level, the average head motion for each subject
was included at the group level analyses, to further minimize its effect
on the identified trajectories.
Previous pediatric studies have found similar results of the motion
effects on edgewise functional connectivity, with the 36 parameters and
the aCompCor being two of the most effective strategies amid other
widely used approaches (Ciric et al., 2017; Taymourtash et al., 2019).
However, when considering graph theory measures only one previous
work in adult populations has tested the in-scanner motion effects.
Fig. 4. Akaike Information Criterion (AIC) is depicted for each graph theory Aurich et al. (2015) tested several preprocessing strategies of motion
measure (y-axis) and every connectivity density (x-axis) here explored. To facil- mitigation on graph theory measures, finding that censoring of “spikes”
itate the visualization, the x-axis was split into two parts from 1–10 and 11–48%.
reduces the association with the head-motion estimates, even lower than
Model abbreviations: linear mixed-effects model for age (LM-Age), linear mixed-
GSR. Nevertheless, they did not test more rigorous strategies such as
effects model for age-sex interaction (LM-Age.Sex), GAMM nonlinear model for
age (NL-Age), GAMM nonlinear model for age by sex (NL-Age.Sex), GAMM with derivative, quadratic trends, PCA, etc. and the inclusion of several miti-
a nonparametric term for PDS (NP-Age), GAMM with nonparametric terms for gators in the same strategy, such as we did in our study. Therefore, when
PDS-Sex interaction (NP-Age.Sex). considering a rigorous amount of motion regressors the head motion ef-
fect on graph theory measures is effectively controlled for the degree,
clustering, and efficiency, but less for path-length (Fig. 2A).

Fig. 5. Scatter-plots of the GAMM LOESS-PDS model at the

25% connectivity density for each graph theory residuals (after
regressing out in-scanner motion and head-coil). Degree (top-
left), clustering coefficient (top-right), efficiency (bottom-left),
and characteristic path length (bottom-right) in relation to the
pubertal scale (PDS). Thin lines represent individual trajecto-
ries; thick black lines represent the sample LOESS curve (with
95% confidence-interval shadow).

5
Z. Gracia-Tabuenca, M.B. Moreno, F.A. Barrios et al. NeuroImage 229 (2021) 117769

Table 1
For every graph theory measure computed at the 25% connectivity density, the Likelihood Ratio Test (LRT), its degrees of
freedom (DOF), and significance (p) of the Pubertal Development Scale LOESS (PDS) and random-effects intercepts (RE) terms
within the Generalized Additive Mixed Model (GAMM).

PDS-LRT PDS-DOF PDS-p RE-LRT RE-DOF RE-p

Degree 10.83 0.77 6.17e-04 36.56 15.52 0.0019

Clustering 12.77 3.04 0.0054 12.97 5.24 0.028
Efficiency 19.65 3.27 2.72e-04 0.0029 1 0.96
Path-length 15.33 2.4 7.89e-04 96.08 35.95 2.21e-07

Table 2
Graph theory measures with significant LOESS term for Pu-
bertal Development Scale (PDS) within the additive model
(GAMM) at the functional network level, corrected for FDR q
< 0.05. Values for PDS Likelihood Ratio Test (LRT), its degrees
of freedom (DOF), and p-value (p).

Measure Network PDS-LRT PDS-DOF PDS-p

Degree DMN 24.07 4.53 1.74e-04

Efficiency CBL 10.74 3.27 0.017
Efficiency DMN 23.72 3.21 3.69e-05
Efficiency DAN 9.51 1.47 0.0043
Efficiency FPN 12.07 1.32 9.01e-04
Efficiency SUB 4.59 0.55 0.014
Efficiency VIS 19.96 5.77 0.0024
Path-Length CBL 8.96 2.33 0.016
Path-Length DMN 25.89 3.98 3.25e-05
Path-Length FPN 13.08 3.53 0.0073
Path-Length SAL 2.7 0.068 0.0045
Path-Length VAN 30.62 6.38 4.26e-05
Fig. 6. Brain maps of the Likelihood Ratio Test of the GAMM LOESS of the
Path-Length VIS 35.69 10.61 1.46e-04
Pubertal Developmental Scale (PDS) term for each graph theory measure at the
functional networks level. Mapping was based on ROI corresponding to the con- Network abbreviations: cerebellar (CBL), dorsal attention
sensus area according to Power et al. (2011), using BrainNet Viewer (Xia et al., (DAN), default-mode (DMN), fronto-parietal (FPN), salience
2013; RRID:SCR_009446). (SAL), subcortical (SUB), ventral attention (VAN), visual (VIS).

4.2. Puberty and neurodevelopmental trajectories Notably, turning points coincide with the start of puberty. No pre-
vious work has focused on the functional organization of the brain net-
Pubertal status showed a faster growth rate in females compared to work; however, animal studies have demonstrated how the start of the
males. Although self-reported questionnaires are sensitive to inaccura- hormone events of puberty provoke plastic changes in the structure
cies (Dorn et al., 2006), the pattern observed has been replicated in other and function of the brain (Laube et al., 2020). Pubertal hormones at
populations (Braams et al., 2015; Wierenga et al., 2018; van Duijvenvo- the onset of puberty are related to an increase of the synaptic pruning
orde et al., 2019). Pubertal scale scores showed a faster growth rate in in frontal areas (Drzewiecki et al., 2016), as well as a growth in in-
females around level 2 and 12 years of age compared to males. Also the hibitory neurotransmission in these areas while not in sensory regions
PDS fits better to brain functional organization measures than chrono- (Piekarski et al., 2017). These results evince the changes in the brain
logical age, with inflection points at the first signs of puberty (level 2 connectivity patterns at lower scales following the onset of puberty, par-
of the PDS). Functional centrality, segregation, and efficiency increases ticularly in association areas. In addition, human studies have shown the
after the inflection, while path-length decreases after such a point, being relevance of pubertal status in gray matter (Goddings et al., 2014), white
more prominent in frontal and parietal functional networks in all cases. matter microstructure (Genc et al., 2017), and even fronto-subcortical

Fig. 7. Fitted line with 95% confidence in-

terval of the network measures standardized
residuals (after regressing out in-scanner mo-
tion and head-coil) along Pubertal Develop-
ment Scale (PDS) for those functional net-
works with a significant (FDR q<0.05) Likeli-
hood Ratio Test for the PDS term. Abbrevia-
tions: Cerebellar (CBL), Default Mode (DMN),
Dorsal Attention (DAN), Fronto-parietal (FPN),
Salience (SAL), Subcortical (SUB), Ventral At-
tention (VAN), and Visual (VIS) networks.

6
Z. Gracia-Tabuenca, M.B. Moreno, F.A. Barrios et al. NeuroImage 229 (2021) 117769

control related systems: dorsal attention (DAN), default mode (DMN),

fronto-parietal (FPN), salience (SAL), and ventral attention (VAN). In
particular, those effects were more prominent than those found in the
whole-brain inferences. These results are congruent with previous work,
which showed an increase in within and between task control networks.
Hwang et al. (2012) reported higher functional connectivity in adoles-
cence within frontal and between frontal and parietal networks, while
Marek et al. (2015) found an increase in functional integration in the
same task control networks along this period. Our results combined with
previous studies coincide with the strengthening of cognitive control
brain networks, which strengthens along with the control-related be-
havior (Mennes et al., 2010; Baum et al., 2017), but ours highlight the
Fig. 8. Chord diagram (Z. Gu et al., 2014) of the GAM Network-Based Statistics relevance of puberty onset for the turning point of brain functional orga-
(NBS) of the relationship between the functional connectivity degree and the nization, particularly those that include frontal and parietal networks.
score (z) of the Card Sorting Task (CST). The component was significant after
correction of a 1000 permutation test with a pFWE of 0.013. 4.3. Brain-cognition effects

Even controlling for the pubertal effects, higher connectivity be-

functional connectivity (van Duijvenvoorde et al., 2019). Other neu- tween attention and task control networks was associated with higher
roimaging studies have found lower gray matter density in cortical and scores in a task evaluating cognitive flexibility. In contrast, this pattern
subcortical areas in adolescence with higher pubertal status when con- was inverse when considering the sensorimotor-hand network (SMH).
trolling for chronological age (Peper et al., 2009; Bramen et al., 2011). Previously, Marek et al. (2015) demonstrated that higher functional in-
Similarly, puberty has been related with higher white matter density tegration (participation coefficient) between the cingulo-opercular and
(Perrin et al., 2009; Herting et al., 2012), which can be interpreted as a salience networks was associated with better accuracy on inhibition
maturity advance (Gogtay et al., 2004; Lenroot et al., 2007; Lebel and tasks. This shows that the integration of attention/task-control com-
Beaulieu, 2011; Coupé et al., 2017). Also, the same effects of lower gray pared to primary/unimodal systems not only reflects the brain orga-
matter and greater white matter density have been related to blood lev- nization development, but also that such patterns are very relevant for
els of testosterone (Paus et al., 2010). These results, along with those cognitive control and cognitive flexibility.
presented here, demonstrate how the onset of puberty greatly impacts
brain structure and function. In particular, we showed an increase in 4.4. Limitations
the efficiency of the brain network. Furthermore, the pubertal markers
have implicit sex dimorphisms that explain the functional neurodevelop- A potential limitation in this study was the missing values in the
ment better than chronological age or sex-age effects. Hence, it is highly pubertal scale scores. Such information was not collected in the first
recommendable that developmental studies in adolescence include esti- sampling, which included participants with ages lower than 10 years
mations of the pubertal stage. old. These values were set to score 1, as similar values have been shown
Regarding the whole-brain functional network, the results showed an in comparable samples (Hibberd et al., 2015; van Duijvenvoorde et al.,
increase after the onset of puberty (around 13 years old; Supplementary 2019). In addition, a minor number of values (n = 8) in older participants
Fig. 1) in brain functional centrality, segregation, efficiency, and inte- were estimated via a GAMM approach. However, the adjustment of the
gration (i.e., the inverse of path length). Although no studies have tested GAMM LOESS to the observed data together with the reproducibility of
the connectome-puberty relationship, previous cross-sectional studies similar curve shapes in other samples, lead us to consider the estimated
have shown a widespread linear increase in functional connectivity PDS scores as a reliable measure of pubertal status.
along with chronological age during adolescence (Hwang et al., 2012; Another limitation is the absence of assessment of the menstrual
Satterthwaite et al., 2013a; Marek et al., 2015). However, these studies cycle in females. Although the menstrual cycle has not been associ-
differ in how this increase impacts functional segregation or integra- ated with changes in functional connectivity in fronto-parietal net-
tion. Satterthwaite et al. (2013a) emphasized that the linear increase in works (Hjelmervik et al., 2014), menstrual cycle phases are related
functional connectivity along with age is more prominent within func- to structural and functional changes regarding subcortical structures
tional networks, which implies an increase in functional segregation. On (Lisofsky et al., 2015; Hidalgo-Lopez et al., 2020).
the other hand, Marek et al. (2015) showed a decline in the functional Also, the time series length for the fMRI sequence was five minutes in
connectivity from childhood to adolescence followed by an increase in total, which is a short scan for current standards in terms of data quality
late adolescence, particularly between network connectivity, i.e., higher (e.g., 20 min for the ABCD study; Casey et al., 2018), nevertheless, it
functional integration. The latter fits with our nonlinear results for de- was considered sufficient at the time of the start of the study (Van Dijk
gree and path length; however, they did not show the same pattern re- et al., 2010).
garding functional segregation. In summary, previous studies differ in Finally, the sample used for the exploration of brain-cognition as-
the evolution of the connectomic functional segregation and integration, sociations spans from 9.67–17.98 years old, which underrepresents the
potentially due to differences in the methodological approaches; never- younger part of the whole sample. To avoid this and any other potential
theless, the implementation of nonlinear splines in the present study bias, LOESS-PDS, Sex, FD-RMS, and coil terms were added as confound-
was able to extract the main effects previously reported. Therefore, our ings to ensure that the brain-cognition inference was not affected by
results contribute to converging those earlier findings, and also to em- those sources of variability.
phasize the relevance of taking nonlinear approaches in developmental
contexts. Additionally, we demonstrate the factual gains of including 5. Conclusion
longitudinal data in the additive models, despite the increase of model
parameters the inclusion of the random-effects intercepts substantially The present work shows that the brain functional organization fol-
improves the goodness of fit of the data. lows a nonlinear trend that fits the pubertal status during the adolescent
Concerning the development of specific functional networks, signif- period. In particular, brain functional centrality, segregation, efficiency,
icant effects were found in the degree, efficiency, and path length af- and integration increased after the onset of the pubertal markers, this
ter controlling for multiple testing, mostly involving attention and task effect being more pronounced in regions related to attention and task

7
Z. Gracia-Tabuenca, M.B. Moreno, F.A. Barrios et al. NeuroImage 229 (2021) 117769

control networks. In addition, higher connectivity between those re- Berg, E.A., 1948. A simple objective technique for measuring flexibility in thinking. J.
gions was associated with better cognitive flexibility performance. In Gen. Psychol. 39 (1), 15–22. doi:10.1080/00221309.1948.9918159.
Betzel, R.F., Byrge, L., He, Y., Goñi, J., Zuo, X.N., Sporns, O., 2014. Changes in structural
summary, pubertal status along with longitudinal and nonlinear pat- and functional connectivity among resting-state networks across the human lifespan.
terns contribute to better description of the developmental trajectories Neuroimage 102, 345–357. doi:10.1016/j.neuroimage.2014.07.067.
of the functional connectome, even converging previous discrepancies Blakemore, S.J., Burnett, S., Dahl, R.E., 2010. The role of puberty in the developing ado-
lescent brain. Hum. Brain Mapp. 31 (6), 926–933. doi:10.1002/hbm.21052.
regarding both segregation and integration increase in adolescence. Braams, B.R., van Duijvenvoorde, A.C., Peper, J.S., Crone, E.A., 2015. Longitudinal
changes in adolescent risk-taking: a comprehensive study of neural responses to re-
wards, pubertal development, and risk-taking behavior. J. Neurosci. 35 (18), 7226–
Data and code accessibility
7238. doi:10.1523/JNEUROSCI.4764-14.2015.
Bramen, J.E., Hranilovich, J.A., Dahl, R.E., Forbes, E.E., Chen, J., Toga, A.W., Sowell, E.R.,
The preprocessed data and the code/software described in the paper 2011. Puberty influences medial temporal lobe and cortical gray matter maturation
is freely available online at https://github.com/BrainMapINB/Pubertal_ differently in boys than girls matched for sexual maturity. Cereb. Cortex 21 (3), 636–
646. doi:10.1093/cercor/bhq137.
Connectome. Casey, B.J., Cannonier, T., Conley, M.I., Cohen, A.O., Barch, D.M., Heitzeg, M.M.,
Orr, C.A., 2018. The adolescent brain cognitive development (ABCD)
study: imaging acquisition across 21 sites. Dev. Cogn. Neurosci. 32, 43–54.
Author statement doi:10.1016/j.dcn.2018.03.001.
Chai, L.R., Khambhati, A.N., Ciric, R., Moore, T.M., Gur, R.C., Gur, R.E., Bassett, D.S.,
Zeus Gracia-Tabuenca: Methodology, Software, Validation, Formal 2017. Evolution of brain network dynamics in neurodevelopment. Netw. Neurosci. 1
(1), 14–30. doi:10.1162/NETN_a_00001.
analysis, Investigation, Data Curation, Writing – Original Draft. Martha Chai, X.J., Castañón, A.N., Öngür, D., Whitfield-Gabrieli, S., 2012. Anticorrelations in
Beatriz Moreno: Methodology, Investigation, Data Curation. Fernando resting state networks without global signal regression. Neuroimage 59 (2), 1420–
A. Barrios: Methodology, Resources. Sarael Alcauter: Conceptualization, 1428. doi:10.1016/j.neuroimage.2011.08.048.
Ciric, R., Wolf, D.H., Power, J.D., Roalf, D.R., Baum, G.L., Ruparel, K., Gur, R.C., 2017.
Methodology, Investigation, Resources, Writing - Review & Editing
Benchmarking of participant-level confound regression strategies for the control of
motion artifact in studies of functional connectivity. Neuroimage 154, 174–187.
Acknowledgments doi:10.1016/j.neuroimage.2017.03.020.
Coupé, P., Catheline, G., Lanuza, E., Manjón, J.V.Alzheimer’s Disease Neuroimag-
ing Initiative, 2017. Towards a unified analysis of brain maturation and aging
We are extremely grateful to the children and parents that volun- across the entire lifespan: a MRI analysis. Hum. Brain Mapp. 38 (11), 5501–5518.
tarily participated in this study. We are grateful to Nuri Aranda López, doi:10.1002/hbm.23743.
Csardi, G., Nepusz, T., 2006. The igraph software package for complex network research.
Leonor Casanova Rico, Luis Concha, César Arturo Domínguez Frausto, InterJ. Complex Syst. 1695 (5), 1–9. http://igraph.org.
Deisy Gasca Martínez, Leopoldo González-Santos, Ma. de Lourdes Lara Di Martino, A., Fair, D.A., Kelly, C., Satterthwaite, T.D., Castellanos, F.X., Thomason, M.E.,
Ayala, Juan J. Ortiz, and Erick H. Pasaye for their technical support. Milham, M.P., 2014. Unraveling the miswired connectome: a developmental perspec-
tive. Neuron 83 (6), 1335–1353. doi:10.1016/j.neuron.2014.08.050.
We are also grateful to Rober Mercadillo for his support in the psycho- Dijkstra, E.W., 1959. A note on two problems in connexion with graphs. Numer. Math. 1
logical assessment and Michael C. Jeziorski for editing the manuscript. (1), 269–271. doi:10.1007/BF01386390.
The authors thankfully acknowledge the imaging resources and support Dorn, L.D., Dahl, R.E., Woodward, H.R., Biro, F., 2006. Defining the boundaries of early
adolescence: a user’s guide to assessing pubertal status and pubertal timing in research
provided by the “Laboratorio Nacional de Imagenología por Resonancia
with adolescents. Appl. Dev. Sci. 10 (1), 30–56. doi:10.1207/s1532480xads1001_3.
Magnética”, CONACYT network of national laboratories. Zeus Gracia Drzewiecki, C.M., Willing, J., Juraska, J.M., 2016. Synaptic number changes in the medial
Tabuenca is a doctoral student of “Programa de Doctorado en Ciencias prefrontal cortex across adolescence in male and female rats: a role for pubertal onset.
Synapse 70 (9), 361–368. doi:10.1002/syn.21909.
Biomédicas, Universidad Nacional Autónoma de México” (UNAM) and
Fair, D.A., Cohen, A.L., Power, J.D., Dosenbach, N.U., Church, J.A., Miezin, F.M., Pe-
received a fellowship (330142) from “Consejo Nacional de Ciencia y tersen, S.E., 2009. Functional brain networks develop from a “local to distributed” or-
Tecnología” (CONACYT). CONACYT had no role in study design, data ganization. PLoS Comput. Biol. 5 (5), e1000381. doi:10.1371/journal.pcbi.1000381.
collection, analyses nor writing the manuscript. This research was par- Fonov, V., Evans, A.C., Botteron, K., Almli, C.R., McKinstry, R.C., Collins, D.L.Brain Devel-
opment Cooperative Group, 2011. Unbiased average age-appropriate atlases for pedi-
tially supported by grant UNAM-DGAPA-PAPIIT IN212219 to SA. atric studies. Neuroimage 54 (1), 313–327. doi:10.1016/j.neuroimage.2010.07.033.
Fornito, A., Zalesky, A., Bullmore, E.T., 2010. Network scaling effects in graph analytic
studies of human resting-state FMRI data. Front. Syst. Neurosci. 4, 22. doi:10.3389/fn-
Supplementary materials sys.2010.00022.
Gao, W., Alcauter, S., Smith, J.K., Gilmore, J.H., Lin, W., 2015. Development of human
Supplementary material associated with this article can be found, in brain cortical network architecture during infancy. Brain Struct. Funct. 220 (2), 1173–
1186. doi:10.1007/s00429-014-0710-3.
the online version, at doi:10.1016/j.neuroimage.2021.117769.
Genc, S., Seal, M.L., Dhollander, T., Malpas, C.B., Hazell, P., Silk, T.J., 2017.
White matter alterations at pubertal onset. Neuroimage 156, 286–292.
References doi:10.1016/j.neuroimage.2017.05.017.
Goddings, A.L., Mills, K.L., Clasen, L.S., Giedd, J.N., Viner, R.M., Blakemore, S.J., 2014.
Akaike, H., 1974. A new look at the statistical model identification. IEEE Trans. Autom. The influence of puberty on subcortical brain development. Neuroimage 88, 242–251.
Control 19 (6), 716–723. doi:10.1109/TAC.1974.1100705. doi:10.1016/j.neuroimage.2013.09.073.
Alcauter, S., Lin, W., Smith, J.K., Short, S.J., Goldman, B.D., Reznick, J.S., Gao, W., 2014. Gogtay, N., Giedd, J.N., Lusk, L., Hayashi, K.M., Greenstein, D., Vaituzis, A.C.,
Development of thalamocortical connectivity during infancy and its cognitive corre- Rapoport, J.L., 2004. Dynamic mapping of human cortical development during
lations. J. Neurosci. 34 (27), 9067–9075. doi:10.1523/JNEUROSCI.0796-14.2014. childhood through early adulthood. Proc. Natl. Acad. Sci. 101 (21), 8174–8179.
Aurich, N.K., Alves Filho, J.O., Marques da Silva, A.M., Franco, A.R., 2015. Evaluating the doi:10.1073/pnas.0402680101.
reliability of different preprocessing steps to estimate graph theoretical measures in Gracia-Tabuenca, Z., Moreno, M.B., Barrios, F.A., Alcauter, S., 2018. Hemispheric
resting state fMRI data. Front. Neurosci. 9, 48. doi:10.3389/fnins.2015.00048. asymmetry and homotopy of resting state functional connectivity correlate
Avants, B.B., Epstein, C.L., Grossman, M., Gee, J.C., 2008. Symmetric diffeo- with visuospatial abilities in school-age children. Neuroimage 174, 441–448.
morphic image registration with cross-correlation: evaluating automated label- doi:10.1016/j.neuroimage.2018.03.051.
ing of elderly and neurodegenerative brain. Med. Image Anal. 12 (1), 26–41. Gracia-Tabuenca, Z., Díaz-Patiño, J.C., Arelio, I., Alcauter, S., 2020. Topological
doi:10.1016/j.media.2007.06.004. data analysis reveals robust alterations in the whole-brain and frontal lobe
Barrat, A., Barthelemy, M., Pastor-Satorras, R., Vespignani, A., 2004. The architec- functional connectomes in attention-deficit/hyperactivity disorder. eNeuro 7 (3).
ture of complex weighted networks. Proc. Natl. Acad. Sci. 101 (11), 3747–3752. doi:10.1523/ENEURO.0543-19.2020.
doi:10.1073/pnas.0400087101. Gracia-Tabuenca, Z., Alcauter, S., 2020. NBR: network-based R-statistics for (unbalanced)
Bates, D., Mächler, M., Bolker, B., & Walker, S. (2014). Fitting linear mixed-effects models longitudinal samples. bioRxiv doi:10.1101/2020.11.07.373019.
using lme4. arXiv:1406.5823. Gu, S., Satterthwaite, T.D., Medaglia, J.D., Yang, M., Gur, R.E., Gur, R.C., Bassett, D.S.,
Baum, G.L., Ciric, R., Roalf, D.R., Betzel, R.F., Moore, T.M., Shinohara, R.T., 2015. Emergence of system roles in normative neurodevelopment. Proc. Natl. Acad.
Cook, P.A., 2017. Modular segregation of structural brain networks supports the Sci. 112 (44), 13681–13686. doi:10.1073/pnas.1502829112.
development of executive function in youth. Curr. Biol. 27 (11), 1561–1572. Gu, Z., Gu, L., Eils, R., Schlesner, M., Brors, B., 2014. circlize implements and enhances
doi:10.1016/j.cub.2017.04.051. circular visualization in R. Bioinformatics 30 (19), 2811–2812. doi:10.1093/bioinfor-
Benjamini, Y., Hochberg, Y., 1995. Controlling the false discovery rate: a practical and matics/btu393.
powerful approach to multiple testing. J. R. Stat. Soc. Ser. B Methodol. 57 (1), 289– Gur, R.C., Richard, J., Calkins, M.E., Chiavacci, R., Hansen, J.A., Bilker, W.B., Abou-
300. doi:10.1111/j.2517-6161.1995.tb02031.x. Sleiman, P.M., 2012. Age group and sex differences in performance on a comput-

8
Z. Gracia-Tabuenca, M.B. Moreno, F.A. Barrios et al. NeuroImage 229 (2021) 117769

erized neurocognitive battery in children age 8− 21. Neuropsychology 26 (2), 251. Petersen, A.C., Crockett, L., Richards, M., Boxer, A., 1988. A self-report measure of puber-
doi:10.1037/a0026712. tal status: reliability, validity, and initial norms. J. Youth Adolesc. 17 (2), 117–133.
Hauser, T.U., Iannaccone, R., Walitza, S., Brandeis, D., Brem, S., 2015. Cognitive flexi- doi:10.1007/BF01537962.
bility in adolescence: neural and behavioral mechanisms of reward prediction error Piekarski, D.J., Boivin, J.R., Wilbrecht, L., 2017. Ovarian hormones organize the matura-
processing in adaptive decision making during development. Neuroimage 104, 347– tion of inhibitory neurotransmission in the frontal cortex at puberty onset in female
354. doi:10.1016/j.neuroimage.2014.09.018. mice. Curr. Biol. 27 (12), 1735–1745. doi:10.1016/j.cub.2017.05.027.
Herting, M.M., Maxwell, E.C., Irvine, C., Nagel, B.J., 2012. The impact of sex, puberty, Power, J.D., Cohen, A.L., Nelson, S.M., Wig, G.S., Barnes, K.A., Church, J.A., Petersen, S.E.,
and hormones on white matter microstructure in adolescents. Cereb. Cortex 22 (9), 2011. Functional network organization of the human brain. Neuron 72 (4), 665–678.
1979–1992. doi:10.1093/cercor/bhr246. doi:10.1016/j.neuron.2011.09.006.
Hibberd, E.E., Hackney, A.C., Lane, A.R., Myers, J.B., 2015. Assessing biological ma- Rubinov, M., Sporns, O., 2010. Complex network measures of brain con-
turity: chronological age and the pubertal development scale predict free testos- nectivity: uses and interpretations. Neuroimage 52 (3), 1059–1069.
terone in adolescent males. J. Pediatr. Endocrinol. Metab. 28 (3–4), 381–386. doi:10.1016/j.neuroimage.2009.10.003.
doi:10.1515/jpem-2014-0187. Satterthwaite, T.D., Wolf, D.H., Ruparel, K., Erus, G., Elliott, M.A., Eickhoff, S.B., Hakonar-
Hidalgo-Lopez, E., Mueller, K., Harris, T., Aichhorn, M., Sacher, J., Pletzer, B., son, H., 2013a. Heterogeneous impact of motion on fundamental patterns of devel-
2020. Human menstrual cycle variation in subcortical functional brain opmental changes in functional connectivity during youth. Neuroimage 83, 45–57.
connectivity: a multimodal analysis approach. Brain Struct. Funct. 1–15. doi:10.1016/j.neuroimage.2013.06.045.
doi:10.1007/s00429-019-02019-z. Satterthwaite, T.D., Elliott, M.A., Gerraty, R.T., Ruparel, K., Loughead, J., Calkins, M.E.,
Hjelmervik, H., Hausmann, M., Osnes, B., Westerhausen, R., Specht, K., 2014. Resting Wolf, D.H., 2013b. An improved framework for confound regression and filtering for
states are resting traits–an FMRI study of sex differences and menstrual cycle effects in control of motion artifact in the preprocessing of resting-state functional connectivity
resting state cognitive control networks. PLoS ONE 9 (7), e103492. doi:10.1371/jour- data. Neuroimage 64, 240–256. doi:10.1016/j.neuroimage.2012.08.052.
nal.pone.0103492. Satterthwaite, T.D., Wolf, D.H., Roalf, D.R., Ruparel, K., Erus, G., Vandekar, S., Verma, R.,
Hwang, K., Hallquist, M.N., Luna, B., 2012. The development of hub architecture in the 2014. Linked sex differences in cognition and functional connectivity in youth. Cereb.
human functional brain network. Cereb. Cortex 23 (10), 2380–2393. doi:10.1093/cer- Cortex 25 (9), 2383–2394. doi:10.1093/cercor/bhu036.
cor/bhs227. Sawyer, S.M., Azzopardi, P.S., Wickremarathne, D., Patton, G.C., 2018.
Jenkinson, M., Bannister, P., Brady, M., Smith, S., 2002. Improved optimization for the The age of adolescence. Lancet Child Adolesc. Health 2 (3), 223–228.
robust and accurate linear registration and motion correction of brain images. Neu- doi:10.1016/S2352-4642(18)30022-1.
roimage 17 (2), 825–841. doi:10.1006/nimg.2002.1132. Somsen, R.J., 2007. The development of attention regulation in the Wisconsin Card Sorting
Jenkinson, M., Beckmann, C.F., Behrens, T.E., Woolrich, M.W., Smith, S.M., 2012. Fsl.. Task. Dev. Sci. 10 (5), 664–680. doi:10.1111/j.1467-7687.2007.00613.x.
Neuroimage 62 (2), 782–790. doi:10.1016/j.neuroimage.2011.09.015. Spielberg, J.M., Forbes, E.E., Ladouceur, C.D., Worthman, C.M., Olino, T.M., Ryan, N.D.,
Kessler, R.C., Berglund, P., Demler, O., Jin, R., Merikangas, K.R., Walters, E.E., 2005. Dahl, R.E., 2014. Pubertal testosterone influences threat-related amygdala–
Lifetime prevalence and age-of-onset distributions of DSM-IV disorders in the na- orbitofrontal cortex coupling. Soc. Cogn. Affect. Neurosci. 10 (3), 408–415.
tional comorbidity survey replication. Arch. Gen. Psychiatry 62 (6), 593–602. doi:10.1093/scan/nsu062.
doi:10.1001/archpsyc.62.6.593. Stam, C.V., Van Straaten, E.C.W., 2012. The organization of physiological brain networks.
Lázaro, J., Ostrosky, F., Lozano, A., 2012. Batería Neuropsicológica de Funciones Ejecuti- Clinical neurophysiology 123 (6), 1067–1087. doi:10.1016/j.clinph.2012.01.011.
vas y Lóbulos Frontales-2. Manual Moderno, México. Stasinopoulos, D.M., Rigby, R.A., 2007. Generalized additive models for location scale and
Laube, C., van den Bos, W., Fandakova, Y., 2020. The relationship between pubertal hor- shape (GAMLSS) in R. J. Stat. Softw. 23 (7), 1–46. doi:10.18637/jss.v023.i07.
mones and brain plasticity: implications for cognitive training in adolescence. Dev. Taymourtash, A., Schwartz, E., Nenning, K.H., Sobotka, D., Diogo, M., Kasprian, G.,
Cogn. Neurosci. 42, 100753. doi:10.1016/j.dcn.2020.100753. Langs, G., 2019. Quantifying residual motion artifacts in fetal fMRI data. In: Smart
Lebel, C., Beaulieu, C., 2011. Longitudinal development of human brain wiring Ultrasound Imaging and Perinatal, Preterm and Paediatric Image Analysis. Springer,
continues from childhood into adulthood. J. Neurosci. 31 (30), 10937–10947. Cham, pp. 171–180. doi:10.1007/978-3-030-32875-7_19.
doi:10.1523/JNEUROSCI.5302-10.2011. Tomasi, D., Volkow, N.D., 2012. Aging and functional brain networks. Mol. Psychiatry 17
Lenroot, R.K., Gogtay, N., Greenstein, D.K., Wells, E.M., Wallace, G.L., Clasen, L.S., (5), 549. doi:10.1038/mp.2011.81.
Thompson, P.M., 2007. Sexual dimorphism of brain developmental trajec- Tustison, N.J., Avants, B.B., Cook, P.A., Zheng, Y., Egan, A., Yushkevich, P.A., Gee, J.C.,
tories during childhood and adolescence. Neuroimage 36 (4), 1065–1073. 2010. N4ITK: improved N3 bias correction. IEEE Trans. Med. Imaging 29 (6), 1310.
doi:10.1016/j.neuroimage.2007.03.053. doi:10.1109/TMI.2010.2046908.
Lisofsky, N., Mårtensson, J., Eckert, A., Lindenberger, U., Gallinat, J., Kühn, S., 2015. Hip- Van Dijk, K.R., Hedden, T., Venkataraman, A., Evans, K.C., Lazar, S.W., Buck-
pocampal volume and functional connectivity changes during the female menstrual ner, R.L., 2010. Intrinsic functional connectivity as a tool for human connec-
cycle. Neuroimage 118, 154–162. doi:10.1016/j.neuroimage.2015.06.012. tomics: theory, properties, and optimization. J. Neurophysiol. 103 (1), 297–321.
Manjón, J.V., Coupé, P., Martí-Bonmatí, L., Collins, D.L., Robles, M., 2010. Adaptive non- doi:10.1152/jn.00783.2009.
local means denoising of MR images with spatially varying noise levels. J. Magn. van Duijvenvoorde, A.C.K., Westhoff, B., de Vos, F., Wierenga, L.M., Crone, E.A, 2019.
Reson. Imaging 31 (1), 192–203. doi:10.1002/jmri.22003. A three-wave longitudinal study of subcortical–cortical resting-state connectivity in
Marek, S., Hwang, K., Foran, W., Hallquist, M.N., Luna, B., 2015. The contribution of adolescence: testing age- and puberty-related changes. Hum. Brain Mapp. 40, 3769–
network organization and integration to the development of cognitive control. PLoS 3783. doi:10.1002/hbm.24630.
Biol. 13 (12), e1002328. doi:10.1371/journal.pbio.1002328. Vijayakumar, N., de Macks, Z.O., Shirtcliff, E.A., Pfeifer, J.H., 2018. Puberty and the hu-
Mennes, M., Kelly, C., Zuo, X.N., Di Martino, A., Biswal, B.B., Castellanos, F.X., man brain: insights into adolescent development. Neurosci. Biobehav. Rev. 92, 417–
Milham, M.P., 2010. Inter-individual differences in resting-state functional con- 436. doi:10.1016/j.neubiorev.2018.06.004.
nectivity predict task-induced BOLD activity. Neuroimage 50 (4), 1690–1701. Wierenga, L.M., Bos, M.G., Schreuders, E., vd Kamp, F., Peper, J.S., Tamnes, C.K.,
doi:10.1016/j.neuroimage.2010.01.002. Crone, E.A., 2018. Unraveling age, puberty and testosterone effects on subcorti-
Muschelli, J., Nebel, M.B., Caffo, B.S., Barber, A.D., Pekar, J.J., Mostofsky, S.H., 2014. Re- cal brain development across adolescence. Psychoneuroendocrinology 91, 105–114.
duction of motion-related artifacts in resting state fMRI using aCompCor. Neuroimage doi:10.1016/j.psyneuen.2018.02.034.
96, 22–35. doi:10.1016/j.neuroimage.2014.03.028. Wood, S., & Scheipl, F. (2014). Gamm4: generalized additive mixed models using mgcv
Paus, T., Keshavan, M., Giedd, J.N., 2008. Why do many psychiatric disorders emerge and lme4. R package version 0.2-3.
during adolescence? Nat. Rev. Neurosci. 9 (12), 947. doi:10.1038/nrn2513. Xia, M., Wang, J., He, Y., 2013. BrainNet Viewer: a network visualization tool for human
Paus, T., Nawaz-Khan, I., Leonard, G., Perron, M., Pike, G.B., Pitiot, A., Pausova, Z., 2010. brain connectomics. PLoS ONE 8 (7), e68910. doi:10.1371/journal.pone.0068910.
Sexual dimorphism in the adolescent brain: role of testosterone and androgen receptor Zalesky, A., Fornito, A., Bullmore, E.T., 2010. Network-based statistic: iden-
in global and local volumes of grey and white matter. Horm. Behav. 57 (1), 63–75. tifying differences in brain networks. Neuroimage 53 (4), 1197–1207.
doi:10.1016/j.yhbeh.2009.08.004. doi:10.1016/j.neuroimage.2010.06.041.
Peper, J.S., Schnack, H.G., Brouwer, R.M., Van Baal, G.C.M., Pjetri, E., Szekely, E., Zuo, X.N., Ehmke, R., Mennes, M., Imperati, D., Castellanos, F.X., Sporns, O., Mil-
Boomsma, D.I, 2009. Heritability of regional and global brain structure at the on- ham, M.P., 2011. Network centrality in the human functional connectome. Cereb.
set of puberty: a magnetic resonance imaging study in 9-year-old twin pairs. Hum. Cortex 22 (8), 1862–1875. doi:10.1093/cercor/bhr269.
Brain Mapp. 30 (7), 2184–2196. doi:10.1002/hbm.20660. Zuo, X.N., He, Y., Betzel, R.F., Colcombe, S., Sporns, O., Milham, M.P., 2017. Human
Perrin, J.S., Leonard, G., Perron, M., Pike, G.B., Pitiot, A., Richer, L., Paus, T., 2009. Sex connectomics across the life span. Trends Cogn. Sci. Regul. Ed. 21 (1), 32–45.
differences in the growth of white matter during adolescence. Neuroimage 45 (4), doi:10.1016/j.tics.2016.10.005.
1055–1066. doi:10.1016/j.neuroimage.2009.01.023.