Reinforcement of the Brain's Rich-Club Architecture Following Early Neurodevelopmental Disruption Caused by Very Preterm Birth.

The second half of pregnancy is a crucial period for the development of structural brain connectivity, and an abrupt interruption of the typical processes of development during this phase caused by the very preterm birth (<33 weeks of gestation) is likely to result in long-lasting consequences. We used structural and diffusion imaging data to reconstruct the brain structural connectome in very preterm-born adults. We assessed its rich-club organization and modularity as 2 characteristics reflecting the capacity to support global and local information exchange, respectively. Our results suggest that the establishment of global connectivity patterns is prioritized over peripheral connectivity following early neurodevelopmental disruption. The very preterm brain exhibited a stronger rich-club architecture than the control brain, despite possessing a relative paucity of white matter resources. Using a simulated lesion approach, we also investigated whether putative structural reorganization takes place in the very preterm brain in order to compensate for its anatomical constraints. We found that connections between the basal ganglia and (pre-) motor regions, as well as connections between subcortical regions, assumed an altered role in the structural connectivity of the very preterm brain, and that such alterations had functional implications for information flow, rule learning, and verbal IQ.

Introduction

The second half of pregnancy is a crucial period for the development of brain connectivity. Thalamocortical fibers reach the cortical plate, leading to the formation of the first sensory-driven circuits. The elaboration of these thalamocortical fibers within the cortical plate promotes the development of cortico-cortical connections, and is quickly followed by the first callosal connections of the cortex (Kostovic and Jovanov-Milosevic 2006). During this period, the subplate peaks in size, before being superseded as the dominant structure by the developing cortex, which quadruples in size over the final 12 weeks of gestation (Huppi et al. 1998). These processes lay the foundation for the brain's complex network of white matter connections, known as the structural connectome.

Two particular network characteristics emerge in the brain by the 30th week of gestation that are preserved throughout life (Ball et al. 2014; Fig. 1A). One of them, modularity, quantifies the ease with which the whole-brain network can be divided into distinct subnetworks, or “modules” (Chen et al. 2008; Rubinov and Sporns 2011). Brain regions that group into a particular module tend to be interconnected with each other and sparsely connected with areas outside the module. The brain's modules do not exist in complete isolation, however. It has been shown that a set of regions, known as a “rich-club,” provides for intramodular connectivity and is thought to play a central role in maintaining network integrity (van den Heuvel et al. 2012). The name “rich-club” originates from the observation that, in social networks, certain (often wealthy) people have a high number of connections, and are also very strongly connected with each other. Similarly, rich-club regions in the brain have a high number of connections and, critically, are also more connected with each other than would be expected by chance. The latter feature effectively implies a redundancy of connections that link segregated modules. Taken together, modularity and rich-club organization characterizes a relatively sparse but, nevertheless, efficient connectivity structure at both local and global levels (Markov et al. 2013).

Figure 1.

Graph-theoretical concepts utilized in this study. (A) Modularity and rich club. Red circles refer to nodes that are likely candidates for the rich club. Blue circles refer to peripheral nodes. Green circles refer to nodes that, despite being well connected, are expected to contribute less to rich-club index than red nodes, as a large proportion of their connections are shared with peripheral nodes (blue circles). (B) Simulated “lesions.” The black boxes represent connectivity matrices, where the intensity of a point indicates the strength of connection between 2 nodes. Red markers in connectivity matrices show the entries that would be assigned zero values in order to simulate a “lesion” to a network; left—simulation of the “lesion” to a node, with all its connection assigned zero values; right—simulation of the “lesion” to a tract.

Individuals who are born very preterm (<33 weeks of gestation) experience a sudden interruption of the typical processes of development during a time in which critical maturational events occur, which affects the subsequent formation of white matter structures (Ajayi-Obe et al. 2000; Inder et al. 2003; Nagy et al. 2009; Nosarti et al. 2014). Brain injuries that are common following very preterm birth, such as subtle diffuse white matter injury and/or acute periventricular hemorrhage and periventricular leukomalacia, can lead to profound white matter alterations and may affect the connectivity of diverse brain regions, due to the periventricular migratory path of many of the brain's connections (Volpe 2009). Connectivity of subcortical gray matter regions, such as the thalamus and striatum, appears to be particularly vulnerable, due to their proximity to the periventricular foci of brain injury (Luciana 2003; Ball et al. 2013). In the long term, the very preterm brain may be required to reorganize its structural connectome around its anatomical constraints, suggesting that very preterm birth can be utilized as a neurodevelopmental model to study the principles that govern the establishment of the brain's connectivity patterns. While a few studies to date have investigated structural and functional brain connectivity in very preterm samples in childhood and beyond (Ball et al. 2013, 2014; Fischi-Gomez et al. 2015; White et al. 2014), it remains to be ascertained whether early white matter damage exerts a long-term impact on the expression of whole-brain network characteristics, including modularity and rich-club organization.

In this study, we examined structural network architecture in preterm infants grown to adulthood. Recent findings suggest that atypical neurodevelopment may be associated with compromised global connectivity. It has been shown that alterations to the structural connectome are implicated in a range of developmental disorders (Crossley et al. 2014). For example, individuals with a diagnosis of schizophrenia display a decrease in biologically expensive (Markov et al. 2013) rich-club organization compared with controls (van den Heuvel et al. 2013). Two hypotheses could be formulated with regard to the very preterm brain. First, the costs of establishing efficient rich-club organization may be difficult to meet following early white matter injury. According to this hypothesis, the very preterm brain would display a decrease in rich-club organization compared with controls, and such a decreased expression of global connectivity might underlie the greater prevalence of psychiatric disorders described in very preterm samples (Nosarti et al. 2012). The alternative hypothesis would be that the (suspected) importance of global connectivity may trigger mechanisms that prioritize the formation of a rich-club architecture, possibly at the expense of other network characteristics such as modularity. This would imply that the characteristics of rich-club organization in the brain of adults who were born very preterm should at least be similar to those seen in controls. A prioritization of the rich club in the face of anatomical constraints could also entail the need for structural network reorganization, whereby certain brain regions may assume an altered role in whole-brain network architecture. The present study intended to explore these 2 alternatives.

Materials and Methods

Participants

Individuals who were born very preterm (before 33 weeks of gestation) were recruited from consecutive cohorts that were born and admitted to the Neonatal Unit of University College Hospital, London (UCLH) between 1979 and 1984. All participants had neonatal ultrasonographic scans and were enrolled for follow-up at 5 developmental phases: infancy (1 year); young and middle childhood (4 and 8 years; Stewart et al. 1989; Roth et al. 1994), and adolescence (14–15 and 18–19 years; Cole et al. 2015; Nam et al. 2015). The sample studied here included all very preterm individuals (n=51) who were part of a further follow-up study who had been assessed from February 2012 to January 2014. Exclusion criteria were any history of neurological complications including meningitis, head injury, and cerebral infections. A sample of 60 age-matched controls was also studied. Controls were recruited from advertisements in the local press and universities. Inclusion criteria were full-term birth (38–42 weeks of gestation). Exclusion criteria were the same used for the very preterm sample, with the addition of any history of birth complications (e.g., endotracheal mechanical ventilation).

Sample statistics for all study participants are summarized in Table 1.

Table 1

Socio-demographic, neonatal, and anthropometric characteristics of the sample

	Very preterm	Controls
Number of subjects	51	60
Males/females	26/25	22/38	χ2 (1) = 1.72, P = 0.19
Age at assessment	29.26 (2.13)	28.88 (3.43)	t(109) = .67, P = NS
Birth weight (g)	1284 (339)	n/a
Gestational age (weeks)	29.21 (2.25)	n/a
Lateral ventricular volume (% of the intracranial volume)a	1.05 (0.55–2.01)	0.82 (0.49–1.29)	t(109) = 2.61, P < 0.01
Socio-economic categories:
I	17.6%	28.3%
II	49.0%	30.0%
III	15.7%	8.3%
IV	0.0%	1.7%
V	2.0%	0.0%
Students	2.0%	23.3%
Unemployed/out of work	13.7%	6.7%
Missing	0.0%	1.7%
Neonatal US classification:
Normal	26	n/a
Uncomplicated PVHb	7	n/a
PVH + DIL	17	n/a
Missing data	1	n/a

Note: SES categories of participants are in accordance with Standard Occupational Classification 1980 (SOC1980).

PVH, periventricular hemorrhage; PVH + DIL, periventricular hemorrhage with ventricular dilatation [see Materials and Methods and Nosarti et al. (2011)].

aStatistics were calculated for the log-transformed data. Mean and confidence intervals (1 SD) were obtained after transforming back to the original scale.

bDue to a small number of participants in this group, no post hoc analyses were performed.

Ethical approval for the study was obtained from the Institute of Psychiatry, King's College London, Ethics Committee (Research). Written informed consent for the assessment, including MRI, was obtained from all participants.

Neonatal Variables

For the very preterm group, the following neonatal risk factors were considered: (1) Length of gestation in weeks; (2) birth weight in grams, and (3) severity of perinatal brain injury, based on results of ultrasonography collected daily for the first 4 days of life, at 1 week, and weekly until discharge from the hospital. A linear-array, real-time ultrasound (US) scanner fitted with a 5- or 7-MHz probe (ADR Model 2130, Ktretztechnik, UK) was used, and the images were stored on a videotape. Diagnosis was reached by consensus between at least 2 neuroradiologists or neonatologists. Hemorrhages limited to the germinal matrix, and those extending to the brain parenchyma or ventricular system were grouped together as periventricular hemorrhage (PVH) (Stewart et al. 1983) and their grade was classified according to Papile et al. (1978). Ventricular dilatation was defined as a clear dilatation of one or both lateral ventricles with cerebrospinal fluid, although not sufficient to meet the conditions for a diagnosis of hydrocephalus; the maximum ventricular width was measured in the coronal plane at the level of the foramen of Monro (Levene and Starte 1981). Ultrasonographic results were summarized as (1) normal (no-PVH), (2) uncomplicated periventricular hemorrhage (grade I-II), without ventricular dilatation (PVH), and (3) periventricular hemorrhage (grade III-IV) with ventricular dilatation (PVH + DIL) (Nosarti et al. 2011). There were no cases of periventricular leukomalacia. In addition, in the analysis of predictors of structural connectivity alterations, we also used ventricular volume (as a ratio of intracranial volume) at current assessment, as it is likely to be related to the extent of early periventricular injury.

MRI Acquisition Parameters

Participants were scanned using a GE Signa HDx 3.0-T MR scanner (General Electric, USA) with an 8-channel head coil. T1-weighted images were acquired using an Enhanced Fast Gradient Echo 3-Dimensional (efgre3D) sequence with the following parameters: repetition time (TR) = 7.2 ms, echo time (TE) = 2.8 ms, flip angle = 20°, which allowed for reconstruction of 1.1 mm3 voxels. For the diffusion-weighted imaging, 60 contiguous near-axial slices were acquired with no gap and the following parameters: rostro-caudal phase encoding, b-value = 1300 s/mm2, TE = 105 ms, voxel size = 2.4 mm3, 32 diffusion-weighted directions, and 4 nondiffusion-weighted volumes, using a spin-echo EPI sequence. Peripheral cardiac gating was applied, with an effective TR of 20/30 R–R intervals.

Parcellation Procedure

Cortical reconstruction was performed using FreeSurfer (Fischl 2012). The cerebral cortex and subcortical regions were then automatically parcellated using geometric information and inference from a hand-dissected training data set (Fischl et al. 2004). Cortical gray matter regions were then dilated by one voxel in order to include the gray/white matter boundary. Subcortical regions were eroded by one voxel on all sides in order to avoid contamination by large neighboring white matter tracts, such as the internal and external capsules.

Diffusion-Weighted MRI Processing

Diffusion-weighted and b0 images were processed as follows. Brain masks were created using FSL's BET (Smith 2002). Motion and eddy-current correction was performed on the masked diffusion data using ExploreDTI (Leemans et al. 2009). A spherical deconvolution approach was chosen to allow for the estimation of multiple fiber directions within a single voxel (Tournier et al. 2004). This was calculated using a damped version of the Richardson-Lucy algorithm (Dell'Acqua et al. 2010). We optimized parameters in order to find the best possible balance between resolving multiple fiber directions and creating minimal spurious (false-positive) fiber orientation distribution (FOD) components. After carefully altering the parameters and testing on multiple subjects, we used a regularization threshold of η = 0.02; a fiber response function = 2; 300 algorithm iterations and regularization parameter v = 20. Fiber orientation estimates were taken from the orientation of the peaks of the FOD profile. We applied an “absolute” (equal to 4 times the amplitude of a spherical FOD obtained from a gray matter voxel) and a “relative” threshold (equal to 7% of the amplitude of the maximum amplitude of the FOD at that voxel) at each voxel to remove the general noise floor and surviving noise local maxima, respectively.

Tractography Algorithm

Whole-brain tractography was performed in native diffusion space using each FOD peak (that survived thresholding) as a seed. Fiber orientation streamlines were propagated using Euler integration with a step-size of 1 mm. Propagation stopped if the track reached an FOD peak that did not survive thresholding, or if the angle between 2 successive steps exceeded 60°. Spherical deconvolution, fiber orientation estimation, and tractography were performed using the in-house software written in Matlab 7.8 (http://www.mathworks.co.uk/products/matlab/) (Dell'Acqua et al. 2013).

Coregistration of Diffusion and Cortical and Subcortical Data

FSL flirt (Jenkinson et al. 2002) was used to register the each individual's diffusion maps to the subject-specific T1-weighted space that contained the cortical and subcortical regions of interest. This transformation matrix was saved, inverted, and used to move the cortical and subcortical regions of interest into each individual's diffusion space.

Connectivity Measure and Graph-Theoretical Analysis

The number of streamlines connecting 2 regions was used as a measure of connection density (van den Heuvel and Sporns 2011). Connectivity matrices were constructed following bias correction and thresholding procedures described in detail in Supplementary Material. The analysis was run in Matlab R2013b using the BCT toolbox (http://www.brain-connectivity-toolbox.net/) and the MIT network analysis toolbox (http://strategic.mit.edu/downloads.php?page=matlab_networks) supplemented with in-house written routines. To enable statistical evaluation of network properties, for each individual network, a sample of 1000 degree-preserving networks was generated using rewiring algorithm proposed in Maslov and Sneppen (2002). Twenty thousand rewiring steps were made.

Principal Measures

Two measures were calculated for each individual's weighted network: (1) rich-club index and (2) modularity index. A raw rich-club index was calculated using the formula proposed in Opsahl et al. (2008). The unnormalized measure proposed that there is somewhat difficult to interpret as it is not independent from the rich-club degree K, or, more precisely, it is negatively related to the number of nodes included in the rich-club pool (see Supplementary Fig. 1). The conventional correction for this dependency involves normalizing raw rich-club indices by rich-club indices derived from randomly generated networks (Opsahl et al. 2008; Allstott et al. 2014). However, this correction may be undesirable in between-group analyses, because rich-club indices derived from randomly generated networks may also differ between groups (see Results below and Discussion). To circumvent this issue, we used a rich-club index obtained by multiplying the raw rich-club index by the ratio of the number of nodes included into the rich-club pool (varying as a function of rich-club degree K) and the number of all nodes in the network. The interpretation of this corrected index is essentially the same as the interpretation of the rich-club index normalized with respect to random networks: A value above 1 indicates an increase in rich-club network architecture. We analyzed rich-club indices calculated for a range of K-values between 6 and 20. Rich-club indices below this range showed very little intersubject variability, explained by inclusion of 70% of nodes into the rich-club pool at the low threshold value. At the high threshold value, the number of nodes included in the rich-club pool was <10% for the majority of subjects (Ball et al. 2014), but could also be as little as 5% for a particular individual. Modularity index was calculated as the average difference between present weights within the modules and ones expected by chance (Newman 2004). The Louvain method for community detection was used to identify modules in the network (Blondel et al. 2008). As certain heuristics are embedded in the algorithm, the modularity index was estimated as an average of 1000 algorithm iterations.

Simulated “Lesioning” Analysis

We investigated 2 types of “lesioning” on modularity and rich-club organization (Irimia and Van Horn 2014): (1) “lesioning,” an entire node and (2) “lesioning,” an individual connection (Fig. 1B). We calculated the normalized differences in rich-club index and, separately, modularity, between the “lesioned” and “non-lesioned” networks for each subject and for each “lesioned” entity, and used their principal component representations (hereafter, “features”) as an input to a linear classifier of group membership (Irimia and Van Horn 2014). To determine a set of features that could reliably differentiate the connectome profiles of the 2 groups (very preterm vs. controls), we combined penalized (LASSO) logistic regression (Tibshirani 1996) with an optimal model selection approach. A detailed description and graphical representation of model fitting pipeline is provided in Supplementary Fig. 2.

Results

Demographic, Neonatal, and Anthropometric Measures

No significant differences were found between the control and very preterm groups in sex (χ2 (1) = 1.72, P = 0.19) and age (t(109) = 0.67, NS). In terms of socio-economic status (SES), the only statistically significant between-group difference at a corrected level of significance was a greater proportion of students in our control sample (P = 0.008). Twenty-six very preterm participants had normal neonatal US results (no-PVH), 7 participants had uncomplicated perinatal brain injury (PVH), and 15 participants had brain injury with complications (PVH + DIL). The PVH group was excluded from post hoc subgroup analyses due to the small number of participants comprising it. The very preterm group as a whole showed an increased volume of the lateral ventricles (t(109) = 2.61, P < 0.01), which was mostly due to the measures in the PVH + DIL group (Fig. 2).

Figure 2.

Ventricular size as a proportion of intracranial volume.

The current very preterm participants did not significantly differ in birth weight or gestational age (both t < 1), and parental SES at birth (P = 0.18) from individuals included in a previous assessment of the same very preterm cohort at 15 years of age (Nosarti et al. 2008).

Topology and Weight Distributions

We first assessed 2 basic characteristics that, from a graph-theoretical analysis perspective, quantify the available resources for network construction. The first one, node-degree distribution (i.e., the distribution of the number of node connections), showed very little between-subject variation (Fig. 3A). Even though some very preterm participants had a smaller number of reconstructed connections compared with controls, the median values did not statistically differ between the groups. In contrast, the second characteristic we assessed, the shape of the consensus distribution of connection weights, suggested a shift toward smaller values in the very preterm group (Fig. 3B). We will refer to this finding as a “relative paucity of white matter resources” in the very preterm group hereafter.

Figure 3.

Perturbation analysis. (A) Average distribution of the degree nodes in the real networks. (B) Consensus distribution (smoothed) for the connection weights in the real networks of each group. (C) Factorizing the process of network randomization into 2 operations. (D) Top row: Average baseline rich-club indices (±SD of a group) as a result of topological, simultaneous topological + weights, and weights-only perturbations. The data are scaled by the grand average across all subjects. Bottom row: same for modularity indices. Asterisks—significant group differences (Bonferroni-corrected).

To evaluate differences between the control and very preterm groups in connectivity parameters, we performed an analysis that aimed to investigate whether group topological and weighting characteristics imposed constraints on a network's capacity to create a rich-club and modular organization. For this purpose, we studied random networks that preserved the degrees and weights of the networks observed in real brains. An alternative view on this kind of analysis is that it assesses network robustness—its resilience to perturbation to specific network characteristics, that is, node degree and/or weight distribution. A recent report demonstrated that the process of generating a random network while preserving the degree (or connectedness) of each node of the original network can be factored out into 2 operations (Allstott et al. 2014; Fig. 3C). First, the pattern of connections between nodes can be permuted in order to create a new unweighted connectivity fingerprint that preserves the same number of connections for the node in question (topological randomness). The second operation does not alter the topology (or existence) of connections, but perturbs weight assignment within a network. Consequently, the baseline for topological perturbations depends solely on the distribution of node degrees, whereas the baseline for weight distribution depends on the existent topology plus the weight distribution.

Statistical analyses were performed using rich-club and modularity indices obtained through averaging across all random networks for each subject separately. Following topological perturbation, rich-club and modularity indices did not differ significantly between the groups (all P > 0.14, uncorrected, Wilcoxon rank-sum test, Fig. 3D), as could be expected when 2 groups are topologically similar. However, following weight perturbation, the very preterm group showed significantly smaller rich-club indices and significantly larger modularity indices for random networks compared with controls (Fig. 3D; Bonferroni-corrected). In summary, the analysis demonstrates that the very preterm brain possesses a limited capacity for constructing global connectivity patterns as a result a relative paucity of the white matter resources in a range of medium-strength connections.

Parameters of Weight Distribution and Neonatal Variables

We performed a correlational analysis in order to investigate the association between weight distribution parameters and neonatal variables. Modularity indices obtained from randomized networks were used as a univariate substitute in this analysis due to their causal dependence on connection weight distribution.

Weight distribution parameters in the group of very preterm individuals who sustained PVH with ventricular dilatation at the time of birth (PVH + DIL) did not differ from those in the very preterm group with normal neonatal US (no-PVH; P = 0.22). Other factors, such as birth weight and gestational age, also did not significantly correlate with baseline modularity indices (Spearman's r = −0.004 and r = −0.084, NS, respectively).

Group Differences in Rich-Club Organization

The above results show that the baseline rich-club-index is dependent on the connection weight distribution. This dependence occurs because the relative paucity of white matter resources implies a smaller probability of assigning at random a large weight to a rich-club connection. In a simulated network this results in a lower rich-club index, and concomitantly, in an increase in modularity. Potentially, this pattern could also be observed in real networks.

Contrary to these expectations, the analysis of real networks has shown that rich-club indices were either significantly greater in the very preterm group at an uncorrected level or numerically higher for a range of rich-club degrees (Fig. 4A). This suggests that rich-club connections have been prioritized in the very preterm brain, as they have been assigned disproportionally greater connection weights. We also observed that the function relating the rich-club indices to the rich-club degree demonstrates a characteristic trend in all studied groups. The function increases at a slow pace for relatively small rich-club degrees and rapidly accelerates for the most highly connected regions, or “hubs,” with this acceleration representing a signature of their high interconnectedness. We will refer to the 2 ranges of rich-club degree as lower-tier and higher-tier rich-club organization. To provide an appropriate statistical cutoff criterion for the analysis of group differences, we note that Bonferroni-corrections for independent tests are unsuitable for the present data due to the presence of significant correlations between rich-club indices for different rich-club degrees. Consequently, we performed a principle component analysis on the rich-club indices for the lower-tier and higher-tier range separately (a range of 6–13 and 13–20, respectively), and used individual loadings on the first component (accounting for 80.5 and 89.5% of variance, respectively) to test for group differences. A t-test analysis showed that the very preterm group significantly differed from controls in lower-tier rich-club architecture (t(109) = 2.34, P = 0.043, corrected), but there was no significant difference in higher-tier rich-club architecture (P = 0.17, uncorrected). There was no significant difference between the PVH + DIL and the no-PVH groups in any comparison of rich-club indices.

Figure 4.

Rich-club indices. (A) Raw anatomical index (±SD). (B) The ratio (±SD) between raw anatomical and averaged simulated baseline indices (following topology + weight perturbation). Values above 1 indicate a rich-club regime.

The contrast between potential resources and an actual state can be further demonstrated by a ratio between raw rich-club indices and those produced by degree-preserving random networks. Consistent with previous studies (van den Heuvel et al. 2012; Irimia and Van Horn 2014), this ratio was greater than 1 for a range of rich-club degrees (Fig. 4B), indicating the presence of a greater-than-expected-by-chance rich-club organization. However, the very preterm group demonstrated systematically greater values than controls, due to a greater difference between baseline and raw anatomical rich-club indices. The difference was particularly strong for lower-tier rich-club architecture (t(109) = 4.46, P < 0.001, higher tier: t(109) = 2.59, P = 0.022; all corrected). A somewhat different perspective on the above observations is provided by the analysis of the proportional cumulative sum of weights allocated to rich-club connections. The difference between the very preterm group and controls was particularly strong for lower-tier rich-club architecture (difference in loadings on the first principal component: t(109) = 4.14, P < 0.001, corrected) and was moderate for higher-tier rich-club architecture (t(109) = 2.41, P = 0.035, corrected).

Group Differences in Modularity

Similar to the rich-club ratio, the ratio between the empirical and randomized network modularity indices was greater than 1 (mean across groups = 1.31, SD = 0.09), which suggests the presence of a modular structure. However, the raw modularity index was similar between the 2 groups (P = 0.21).

Results of Simulated “Lesioning”

The above analysis of global network properties revealed that, despite an apparent paucity of white matter resources, the very preterm brain demonstrates a higher level of rich-club organization than the control brain, suggesting putative structural reorganization in response to a shorter gestational age. The concept of structural reorganization is usually associated with a process, whereby certain network components assume a differential role in network architecture. The lack of topological connectivity differences between the groups argues in favor of this possibility, and against the possibility of “novel” connections in the very preterm group. To investigate this possibility, we used a simulated “lesion” approach that allows for a statistical evaluation of the importance of a network component in supporting specific network characteristics.

Analyses showed that changes in modularity indices, as a result of lesioning a node or a connection, did not reliably improve linear classification of group membership. Similarly, changes in rich-club organization, as a result of lesioning a node, were also not significant predictors of group membership, suggesting that no individual brain region plays a significant role in determining group differences. On the other hand, a single principal component resulting from lesioning individual connections significantly improved classifier performance (P = 0.003). This principal component was observed for a relatively small rich-club degree (node degree = 7), suggesting that the structural alterations preferentially affect lower-tier rich-club architecture. The set of connections that had the greatest coefficients for this component were the connections between the basal ganglia and (pre-)motor cortex in the left hemisphere (Fig. 5) and the connections between subcortical gray matter structures in the right hemisphere. The complete list of connections included those of the left precentral gyrus with the left putamen and pallidum, of the left pallidum with the left superior frontal gyrus (including the premotor area), of the right hippocampus with the right thalamus, and of the right putamen with the right caudate. The strength of the “lesioning” effect was lateralized. “Lesioning” connections in the left hemisphere resulted in a greater decrease in rich-club organization in the very preterm brain compared with the control brain than lesioning the connections in the right hemisphere. Overall, only when subject loadings on the predictive component were used for group classification, classifier accuracy in leave-one-out cross-validation was 74.8%, with controls classified correctly in 80.0% of cases. The composition of the very preterm group appeared to be somewhat more heterogeneous, with classification accuracy of 68.6%. To verify that the observed results were due to the role played by these connections in global connectivity architecture, and not due to group differences in connection density, we repeated the above analysis using connection densities as predictors. This analysis failed to identify a reliable predictor of group classification.

Figure 5.

Basal ganglia connections with precentral and superior frontal gyri in the left hemisphere, showing an altered role in very preterm brain network. Table presents the complete list of connections associated with the principal component that improved classifier performance.

Connectivity Alterations and Neonatal Variables

To investigate in detail the source of heterogeneity of the very preterm group with respect to the alterations detected using a “lesion” approach, we used subject loading on the principal component as a measure of the degree of connectivity alterations and investigated its association with birth weight, gestational age, and severity of neonatal brain injury (Fig. 6A). There was a significant correlation between classifier confidence and birth weight (r = −0.31, P = 0.03), indicating that the smaller the birth weight the greater the degree of alteration. Gestational age at birth appeared to be a less reliable predictor of classifier confidence (r = 0.20, P = 0.15). Finally, it was easier to distinguish the PVH + DIL group than the no-PVH group from controls (t(39) = 3.15, P = 0.003), and when considered alone, the classification of individuals with PVH + DIL was 88.2% accurate.

Figure 6.

Correlations of altered connectivity (principal component values—PC) with neonatal, anatomical, and cognitive variables. (A) Classifier confidence in relation to severity of neonatal US classification, birth weight, and ventricular volume at assessment in the very preterm group. ICV, intracranial volume. (B) Relation of alterations in connectivity and cognitive outcome. Circled data points mark very preterm participants. For display purposes only, scales for behavioral performance are in accordance with the most suitable transformations applied to correct for nonnormality of the distributions (Table 2). The significance of correlations was tested using nonparametric statistics.

Table 2

Behavioral data: descriptive statistics, group contrasts, and correlations with the predictive principal component

	Controls			Very preterm			Control vs. Very preterm		Correlation with subject loading on PC
	Mean	−/+ SD		Mean	−/+ SD		z-value	P-value	r-value	P-value
IQ verbal (82)*	114.8	98.3	126.6	107.2	91.1	120.2	−2.18	0.030	−0.27	0.016
IQ performance (83)*	117.8	105.5	126.8	107.9	91.7	120.7	−2.82	0.005	−0.15	0.183
Phonological Verbal Fluency (104), N	14.7	11.1	18.4	13.2	9.7	16.7	−1.77	0.078	−0.16	0.10
Semantic Verbal Fluency (103), N	24.7	18.7	30.7	22.3	16.2	28.5	−1.59	0.111	0.00	0.97
Trail making time (95)**, sec	59.1	42.2	82.8	80.5	51.5	125.8	−3.28	0.001	0.05	0.64
Hayling (Part A) Initiation (100)**, sec	4.3	0.9	13.7	4.5	1.0	14.2	−0.18	0.86	0.33	0.001
Hayling (Part B) Inhibition (101)***, sec	19.0	3.9	45.5	30.1	8.8	64.2	−2.30	0.02	−0.19	0.06
Intra-Extra Dimensional Shift (92)**, errors	12.1	6.4	23.0	19.0	8.7	41.6	−3.05	0.002	−0.22	0.033
Paired Associate Learning (91)**, errors	4.2	1.2	11.7	6.95	3.09	14.45	−2.38	0.017	0.00	0.99

Note: For z-statistics, negative values indicate worse performance in the very preterm group. The mean and SD for most of the measurements were calculated for transformed data (in order to correct for nonnormality) and then rescaled backwards onto original raw scale. Transformations: *, arcsine; **, log; ***, square-root. In parenthesis—the number of participants for whom the data were available. Nonparametric statistics were used for hypothesis testing. Shaded cells highlight significant statistics.

It should be noted that neonatal variables are correlated in the population we studied (White et al. 2014). To identify the best predictors of alterations in global connectivity, we performed a model comparison analysis where we used all combinations of (1) (the log of) lateral ventricular volume at current assessment (taken as a proportion of the intracranial brain volume), (2) gestational age, (3) birth weight, and (4) neonatal US diagnosis as linear predictors for the subject loading on the predictive principal component. We found that the model that included ventricular volume as the only predictor outperformed all other models on the basis of the Bayesian Information Criterion and explained 44.7% of the variance. A separate analysis also showed that ventricular volume was more strongly correlated with classifier confidence in the very preterm group than in controls (r = 0.73 vs. 0.36, respectively, Fisher's Z = 2.81, P = 0.005). We did not observe a significant effect of laterality of ventricular dilatation on classifier confidence, with predicted variance being very similar for the volumes of the left and right ventricles (r2 = 0.43 and 0.42, respectively). These results could be explained by a high correlation between ventricular volumes in the 2 hemispheres (r = 0.92).

Association Between Network Characteristics and Cognitive Outcomes

Finally, we performed correlational analyses in order to investigate the association between weight distribution (summarized by the modularity of surrogate random networks) and the classifier confidence described above on the one side, and a range of neuropsychological measures on the other (see Supplementary Material for the description of cognitive tasks). Consistent with previous findings, very preterm adults performed worse than controls on a number of neuropsychological tests (Eryigit Madzwamuse et al. 2014; Nosarti et al. 2014; Solsnes et al. 2014; Table 2). We found that, across groups, classifier confidence correlated with scores on verbal IQ (assessed using the WASI), the HSCT, Part A, which measures the speed of response initiation and the number of errors in the IEDS task, which assesses rule learning (Table 2 and Fig. 6B). Notably, the latter is considered as a normative task for measuring frontostriatal deficits (Downes et al. 1989). For verbal IQ and rule learning, the correlation between classifier confidence and performance was negative, that is, the greater the alterations resulting from simulated “lesions” to basal ganglia connections, the worse the cognitive performance, whereas for speed of response initiation, the correlation was positive. Scores on memory (paired-associates learning), and various aspects of executive functions, namely, response inhibition (HSCT, Part B), conceptual tracking (Trail Making Test), and cognitive flexibility (Trail Making Test, Verbal Fluency), did not show reliable correlations with classifier confidence. Baseline modularity indices, summarizing connection weight distribution, significantly correlated only with scores on the HSCT, Part A (r = 0.233, P = 0.019), although this appeared to be due to the collinearity of the baseline modularity index with classifier confidence (r = 0.47, P < 0.001). When a partial correlation between baseline modularity indices and scores on the HSCT, Part A, controlling for classifier confidence was performed, results were no longer significant (r = −0.06, NS). Conversely, a partial correlation controlling for baseline modularity indices barely affected the association between classifier confidence and performance on the HSCT, Part A (r = −0.24, P = 0.013).

Discussion

Very preterm birth has been associated with a cascade of neurodevelopmental alterations in infancy, childhood, and adolescence (Nosarti et al. 2008; Nagy et al. 2009; Ball et al. 2013; Groeschel et al. 2014), which are likely to exert a long-term impact on the brain's structural connectivity. However, the principles governing the processes associated with the formation of the structural connectome remain unclear. Based on graph-theoretical analysis, the present study provides novel evidence in support of the idea that global structural connectivity may be prioritized over peripheral connectivity during neurodevelopment. We showed that the very preterm brain disproportionately assigns a larger share of white matter resources to its rich club. Our analysis of randomly perturbed networks demonstrated that the shift towards lower values in a connection weight distribution in the very preterm brain (the evidence for relative paucity of white matter resources) resulted in a smaller potential for expressing global network structure (and conversely for modular structure). However, while the actual modularity characteristics of the very preterm and control brain were similar, rich-club organization was quantitatively greater in the very preterm group (Ball et al. 2014). Furthermore, the patterns of white matter resource distribution do not appear to be affected by neonatal US classification, birth weight and gestational age, indicating that they may represent a feature characterising a general “preterm phenotype”.

The findings of a long-term effect of prematurity on rich-club architecture contrast the findings in preterm-born infants showing unaltered rich-club architecture when compared with the term-equivalent controls (Ball et al. 2014). This suggests a possibility that anatomical constraints of the very preterm brain give rise to putative compensatory processes at later developmental stages. To test this hypothesis, we used a simulated “lesion” approach (Alstott et al. 2009; Irimia and Van Horn 2014) that allows one to evaluate the importance of individual network components in supporting network characteristics. At a theoretical level, simulated “lesioning” explores “network robustness,” defined as the ability to preserve a structure intact (unchanged or minimally altered) as a result of perturbation (Holmgren 2007). For example, the presence of highly connected nodes, or “hubs,” has been shown to decrease a network's vulnerability to insult that may occur in a randomly targeted network component (Albert et al. 2000). Consequently, differential susceptibility to a “lesion” can be considered as a proxy for the price paid by the very preterm brain for a reorganization around its anatomical constraints.

We found that the very preterm and control groups could be differentiated on the basis of changes in rich-club organization as a result of “lesioning” connections of subcortical gray matter regions, predominantly belonging to subdivisions of the basal ganglia. These findings suggest an altered role of basal ganglia connections in supporting a global exchange of information and in preserving structural network integrity in the very preterm brain. Notably, the connectivity alterations observed in our “lesion” analyses could be related to neonatal risk factors in very preterm individuals. These alterations were more pronounced in those very preterm-born individuals who suffered moderate-to-severe PVH and/or who had a low birth weight. Increased ventricular volume, which can be regarded as a long-term consequence of early brain injury, was the most reliable predictor of these alterations in a model comparison analysis.

The findings of an altered role of basal ganglia in very preterm whole-brain connectivity provide a new perspective on earlier findings showing structural and functional alterations in cortico-basal ganglia-thalamocortical loops associated with very preterm birth (Ball et al. 2012, 2013; Meng et al. 2015). Notably, direct connections between pallidum and neocortex, which are highly reconstructable in diffusion-weighted imaging (Draganski et al. 2008), are not considered to be part of the cortico-basal ganglia-thalamocortical loops. Several studies ascribe the influence of pallidum on cortical functions to GABAergic modulation of thalamic activity (Alexander et al. 1986, 1990; Joel and Weiner 1994). However, animal studies provide some evidence that can help to reconcile our findings with previously proposed models. As internal and external pallidum nuclei are surrounded by layers of cholinergic neurons, projecting directly to the neocortex (Parent et al. 1981; Parent and De Bellefeuille 1982; Parent and Hazrati 1995), it is likely that a restructuring of basal ganglia loops following neonatal brain injury may affect the role of these cholinergic modulators of cortical activity in order to compensate for early disruption of thalamocortical connectivity that has been observed in very preterm infants (Ball et al. 2013). This interpretation remains somewhat speculative given that we do not find a strong evidence for this disruption in adult very preterm individuals. Future studies could address this issue by attempting to relate alterations in basal ganglia structural connectivity and alterations in functional networks reported in a number of very preterm studies (Bauml et al. 2015; White et al. 2014).

The basal ganglia play a crucial role in numerous functions, including motor learning (Smith et al. 2000; Lehericy et al. 2005), information routing (Stocco et al. 2010), and procedural and stimulus–response association learning (Knowlton et al. 1996). Their involvement in language, possibly reflecting an information-routing component and application of grammatical rules, is also documented (Damasio et al. 1982). We previously observed functional alterations in frontostriatal networks in very preterm-born adolescents during completion of a response inhibition task (Nosarti et al. 2006). Alterations to the function and structure of frontostriatal connections have been found in dopaminergic disorders (Owen 2004; Howes and Kapur 2009), which suggest that impaired regulation of dopaminergic signaling via frontostriatal circuits may contribute to cognitive deficits and the increased risk of developing psychotic and motor disorders in very preterm individuals (Pitcher et al. 2011; Nosarti et al. 2012). Here, we showed that alterations in basal ganglia connectivity correlate with some aspects of cognitive outcome, such as verbal IQ, speeded completion of sentences on the HSCT, and rule learning, as measured by the IEDS task. On the other hand, executive and memory functions did not show reliable associations with changes in basal ganglia connectivity. These results suggest that the observed alterations were related to a particular cognitive component of IQ, associated with an efficient flow of information and rule learning, rather than inhibition or planning.

Taken together, our results suggest that the observed shift in connectivity distribution, indicative of the paucity of medium-strength connections, does not in itself present an obstacle for a typical neurodevelopment. However, the presence of neonatal brain injury appears to be a critical factor, which can hinder typical neurodevelopment and result in specific cognitive deficits.

Implications for Graph-Theoretical Analysis

Our results challenge 2 conventional approaches to the analysis of the brain's rich-club architecture. The first convention is to use normalized rich-club indices in the study of group differences. This metric represents a ratio whose value depends both on a nominator and a denominator. As demonstrated in our analysis, the interpretation of results would differ depending on which of the 2 values is affected: Group differences in nominator (i.e., raw rich-club index) reflect differences in connectivity patterns, whereas differences in denominator (i.e., the normalizing value) reflect differences in resources that are available to support these connectivity patterns. Here, we showed that the denominator is more critical for differentiating between very preterm and control groups. The second convention is to associate rich club with “hubs,” that is, highly connected regions. This convention, however, is not mathematically necessitated, as a rich-club index can be calculated for an arbitrary rich-club degree. More generally, rich-club organization can be analyzed as a tiered quality that reflects the degree of network integrity at different connectivity levels. Even though the expression of rich-club property is particularly strong at the hub connectivity level, it is by no means negligible at the lower-tier level, where in fact the most pronounced between-group differences in rich-club architecture were found.

The above considerations are also relevant to the contrast between the significant between-group differences in rich-club organization with nonsignificant findings in modularity. In other words, the evidence that the rich-club is prioritized over peripheral connectivity (as implied by the fact that the rich-club index represents a ratio, whose magnitude is inversely related to the strength of nonrich-club connections) does not suggest that this occurs at the expense of modularity. To explain these findings, it should be noted that the tradeoff between modularity and rich-club organization is not straightforward, that is, an increase in one characteristic does not always imply a decrease in the other. A study by Zamora-Lopez et al. (2010) showed that the rich club, even at its highest degree, does not exist independently from modular organization but rather it is formed by a partial overlap of modules. Strengthening a part of rich-club connectivity, which is formed by the partial overlap of modular organization, may lead to a decrease of modularity; strengthening the supramodular connectivity may lead to the reinforcement of the rich-club module proper, that is, increase in modularity. The reason for this is that modularity index reflects how well nodes are integrated within “different” modules, whereas rich-club index quantifies how well a subset of nodes, selected on the basis of their high connectivity degree, is integrated into “one” common module. Whether the reinforcement of a particular rich-club connection leads to a decrease in modularity would depend on specifics of an actual network topology. Overall, the present results suggest that modularity characteristics appear to be more robust to the effects of very preterm birth, particularly given that no gross topological alterations could be detected in the majority of very preterm subjects.

Limitations

There are several limitations to our study. First, evidence for structural brain reorganization was obtained using quantitative analysis of in silico lesions. While this tool is becoming increasingly popular (Alstott et al. 2009; Crossley et al. 2014; Irimia and Van Horn 2014), future research should establish whether these findings could be related to other measures, for example, to the expression of dopamine signaling or measures of functional connectivity. Secondly, in this study, we used neonatal US classification as our primary perinatal clinical factor, as this information had been systematically collected in our sample. We acknowledge, however, that it would have been preferable to include information about the underlying causes of brain alterations detected via neonatal US (i.e., PVH and ventricular dilatation), such as hypoxic and ischemic events, which are regarded as the primary causes of damage to cortical pathways (Volpe 2009). Furthermore, clinical classification based on ultrasonographic scans does not allow for fine-tuned quantitative analyses and, in the absence of MRI data at the time of birth, we resorted to using ventricular volumes at current assessment as a parametric approximation of the long-term consequences of neonatal brain injury. With the current data, we are therefore unable to ascertain whether the primary cause of ventricular dilation is general white matter atrophy due to the lesion of specific periventricular pathways, or subsequent white matter loss or both.

Another limitation of this study regards the representativeness of our very preterm sample who had above average IQ scores (although these were lower than those in our control group) and is over-represented with nonmanual professions. This is due to the fact that, at birth, subjects from our cohort lived within the catchment area of a major teaching Hospital in central London, which includes several neighborhoods, which are considered “wealthy.” In addition, the very preterm group was not fully matched with the control group, with a greater proportion of student among the controls. The latter factor, however, is unlikely to be important for interpretation of results, given that majority of those students will move to nonmanual professions after completing their studies. Overall, the findings from a high-functioning very preterm group indicate that the observed alterations in rich-club index cannot be attributed to cognitive impairments but rather suggest a propitious structural reorganization (Baggio et al. 2015).

Finally, on a technical note, connectome analysis requires a set of arbitrary decisions about thresholding and bias correction, with standards yet to be defined (van den Heuvel and Sporns 2011; van den Heuvel et al. 2012, 2013). In our study, we used a simple diagnostic tool, presented in Supplementary Material, to determine study-specific thresholds which we believe may offer advantages over using values informed by previous studies, as they take into consideration differences in scanning sequences and tractography algorithms between our current dataset and others. The criteria derived using this method determined a sparser connectivity pattern than that found in many whole-brain connectivity studies (van den Heuvel and Sporns 2011; van den Heuvel et al. 2012, 2013). Further work is, however, required to test for a more general applicability of this method.

Conclusion

To conclude, our study indicates that very preterm birth may be used as a neurodevelopmental model for investigating compensatory mechanisms associated with early brain injury. Statistically, the very preterm group introduces a variance that may be rare in the general population. Our data suggest that the greater contributors to this variance are very preterm-born individuals who sustained severe neonatal brain injury. The uniqueness of the very preterm model, in comparison with other influential models, such as Parkinson and Huntington's diseases, is that, rather than just reflecting the consequences of a degenerative condition, it provides an insight into the long-term neuroplastic processes occurring in the developing brain. Furthermore, such a model increases our understanding of the limitations of neuroplastic processes with respect to behavioral outcomes from a developmental perspective.

Supplementary Material

Supplementary material can be found at: http://www.cercor.oxfordjournals.org/.

Funding

The study was supported by Medical Research Council, UK (ref. MR/K004867/1). Funding to pay the Open Access publication charges for this article was provided by Research Councils UK.