Differential cortical microstructural maturation in the preterm human brain with diffusion kurtosis and tensor imaging.

During the third trimester, the human brain undergoes rapid cellular and molecular processes that reshape the structural architecture of the cerebral cortex. Knowledge of cortical differentiation obtained predominantly from histological studies is limited in localized and small cortical regions. How cortical microstructure is differentiated across cortical regions in this critical period is unknown. In this study, the cortical microstructural architecture across the entire cortex was delineated with non-Gaussian diffusion kurtosis imaging as well as conventional diffusion tensor imaging of 89 preterm neonates aged 31-42 postmenstrual weeks. The temporal changes of cortical mean kurtosis (MK) or fractional anisotropy (FA) were heterogeneous across the cortical regions. Cortical MK decreases were observed throughout the studied age period, while cortical FA decrease reached its plateau around 37 weeks. More rapid decreases in MK were found in the primary visual region, while faster FA declines were observed in the prefrontal cortex. We found that distinctive cortical microstructural changes were coupled with microstructural maturation of associated white matter tracts. Both cortical MK and FA measurements predicted the postmenstrual age of preterm infants accurately. This study revealed a differential 4D spatiotemporal cytoarchitectural signature inferred by non-Gaussian diffusion barriers inside the cortical plate during the third trimester. The cytoarchitectural processes, including dendritic arborization and neuronal density decreases, were inferred by regional cortical FA and MK measurements. The presented findings suggest that cortical MK and FA measurements could be used as effective imaging markers for cortical microstructural changes in typical and potentially atypical brain development.

During the third trimester, complicated yet precisely regularized cellular and molecular processes including dendritic arborization (1, 2), axonal growth (3), synapse formation (4), and neuronal apoptosis (5, 6) take place in the cerebral cortex. The radial glia scaffold uniquely exists in the early developing brain and serves as a neuronal migration pathway (2, 7, 8). Dendritic arborization disrupts the radial glia scaffold, composed by glial cells in a well-defined radially organized cytoarchitecture (8, 9), and causes dramatic microstructural changes in the cortical plate. In parallel, cell death or apoptosis (5, 6) reduces the neuronal density in the cortical plate. Based on histological studies (e.g., ref. 4), spatial and functional gradients of the cytoarchitectural processes across the cortex were suggested. Besides the maturational processes in the cortex, white matter (WM) associated to certain cortical regions also develops rapidly and differentially through the growth and myelination of axonal tracts (10, 11). Cortical microstructural maturation plays a vital role in indicating neuronal circuit formation and was investigated primarily through histological studies (12, 13). However, most histological studies so far have been limited to localized and discrete cortical regions. Quantifying microstructural changes across the entire cerebral cortex in this early developmental stage could offer fresh insights into maturational signature of cortical microstructure, shedding light on differentiated emergence of certain brain functions.

Cortical microstructural development across the cerebral cortex has been investigated using diffusion tensor imaging (DTI), a noninvasive magnetic resonance imaging (MRI) technique measuring Gaussian diffusion properties of the water in the brain tissue (14). Fractional anisotropy (FA) (15), a DTI-derived metric, characterizes the degree of anisotropy of water diffusion and has been applied to infer the cortical microstructure. High FA values were usually observed in the immature cortex (e.g., refs. 16–23) characterized by a more organized radial glia scaffold. Decreases in cortical FA during the development may result from the increasing synapse formation, myelination of intracortical axons, and dendritic arborization that disrupt the highly organized radial architecture in the cortical plate. Such FA decreases have been documented in a number of studies on developing human (16–28) and animal cerebral cortices (29–35). With the advantage of whole-brain coverage, these differentiated regional cortical FA decrease patterns could be used to infer the spatiotemporal processes of dendritic arborization and potentially brain circuit emergence.

An intrinsic limitation of DTI is that the diffusion tensor model approximates biological water diffusion with a Gaussian distribution (14). Non-Gaussian diffusion is believed to arise from diffusion barriers such as axonal sheaths, cellular membranes, and organelles (36). The deviation from Gaussian behavior of water diffusion becomes more prominent in the tissues with more complex or heterogeneous microenvironments, such as the cerebral cortex. Diffusion kurtosis imaging (DKI) characterizing fourth-order tensor is an extension of conventional DTI representing second-order tensor. DKI is capable of quantifying the non-Gaussian properties of water diffusion (36, 37) with multishell diffusion MRI (dMRI) acquisition. Mean kurtosis (MK), a DKI-derived metric, quantifies tissue microstructure complexity and is not limited to an anisotropic environment. Several studies have shown the promise of DKI in contributing to a better understanding of the microstructure in gray matter or WM from infancy to adulthood (38–41). MK was demonstrated to continue to change even after FA has reached its plateau in some WM and gray matter nuclei regions in the developing brain of children (38), suggesting that MK measurements may offer microstructural information complementary to FA measurements in the developing cortex. Surveying microstructural changes of the entire cerebral cortex with combined DKI and DTI measurements may offer comprehensive insight into differential cortical microstructural changes during the early development age in the third trimester. With its strength of noninvasive and comprehensive coverage, the differential disruption of the radial glia scaffold inferred by FA measures and reduction of diffusion barrier inferred by MK measures could be delineated across the cerebral cortical regions in typical brain development. Cortical FA and MK may also serve as clinical biomarkers by quantifying the alteration of FA or MK change curves in atypical development.

In this study, we hypothesized that both cortical MK and FA change patterns are differential across the cortical regions with cortical MK and FA change curves distinguished from and complementary to each other. Multishell high-resolution (1.5 × 1.5 × 1.6 mm3) dMRI with b values of 1,000 and 1,600 s/mm2 were acquired from 89 preterm neonates with ages of 31–42 postmenstrual weeks (PMW) at scan. After diffusion tensor and kurtosis fitting, “core” cortical FA and MK values were obtained on the cortical skeleton to alleviate partial volume effects in the measurements. Spatiotemporal maps of cortical FA and MK measurements and cortical FA and MK changing curves of all regions were delineated. A data-driven clustering approach based on nonnegative matrix factorization (NMF) was adopted to identify the cortical MK or FA regions with similar changing patterns, respectively. Furthermore, the microstructural changes of WM tracts traced from the cortical MK clusters were revealed to be synchronized with MK changes at these cortical clusters. Support vector regression (SVR) analysis revealed high accuracy of predicted postmenstrual ages with cortical MK and FA measurements. A hypothetic schematic depiction was presented to demonstrate the cellular density change and disruption of radial glia scaffold inferred by the cortical MK and FA changes, respectively. This study comprehensively investigated differential maturation in early developmental brains with DKI and DTI.

Results

Heterogeneous Cortical Microstructural Profile of DKI and DTI Measurements at Individual PMW.

The cortical microstructural profile with MK measurements based on the fourth-order kurtosis tensor is shown in Fig. 1. At each week during 33–40 PMW, cortical MK values are heterogeneous on the lateral (upper row) and medial (lower row) surface of the preterm brain. For the younger neonate brain at 33 PMW, higher MK values are prominent at the lateral occipital and primary visual cortex. For the older neonate brain around term (e.g., 40 PMW), lower MK values are clustered around the superior frontal and primary sensorimotor cortical regions. The cortical FA microstructural profile (Fig. 2) at each week of 32–41 PMW is distinguished from cortical MK microstructural profile at the same age but is also heterogeneous across the cortical regions. For younger neonate brains at 32 PMW, FA values at most lateral cortical regions such as the frontal and temporal regions are relatively high with lower FA values shown only around the primary sensorimotor areas. The age-specific cortical profiles of other DKI-derived metric measurements [e.g., radial and axial kurtosis (RK and AK)] and DTI-derived metric measurements [mean, axial, and radial diffusivity (MD, AD, and RD)] of the preterm brains demonstrate distinctive patterns as shown in in .

Fig. 1.

(A) Spatiotemporal maps of cortical MK measurements of preterm brains at eight consecutive time points from 33 to 40 PMW displayed on a partially inflated cortical surface of a template brain. (B) Linear fitting of MK time courses at different cortical gyrus overlaid onto the template brain in the left hemisphere displayed in lateral (Left) and medial (Right) view. In each scatter plot, the x axis indicates postmenstrual age in PMW, and the y axis indicates MK values. The cortical gyri underlying the scatter plots are listed below: a, inferior frontal; b, middle frontal; c, superior frontal; d, precentral; e, postcentral; f, superior parietal; g, superior occipital; h, middle occipital; i, inferior occipital; k, angular; j, supramarginal; l, superior temporal; m, middle temporal; n, inferior temporal; o, lateral orbitofrontal; p, insula; q, cuneus; r, precuneus; s, cingulate; t, medial orbitofrontal; u, gyrus rectus; v, fusiform; w, lingual.

Fig. 2.

(A) Spatiotemporal maps of cortical FA measurements of preterm brains at 10 consecutive time points from 32 to 41 PMW displayed on a partially inflated cortical surface of a template brain. (B) Biphasic linear fitting of FA time courses at different cortical gyrus overlaid onto the template brain in the left hemisphere displayed in lateral (Left) and medial (Right) view. In each scatter plot, the x axis indicates postmenstrual age in PMW, and the y axis indicates FA values. The cortical gyri underlying the scatter plots are the same as those shown in Fig. 1.

Significant and Differential Cortical Microstructural Changes with DKI and DTI Measurements.

Significant and distinctive age-dependent cortical MK and FA decreases across cortical regions during preterm brain development were found with both region of interest (ROI)-based and data-driven analyses. Statistically significant [P < 0.05, false discovery rate (FDR) corrected] age-dependent decreases of cortical MK in 33–40 PMW were observed over all of the cortical gyri (Fig. 1 and ) with each gyrus as an ROI. These MK time courses are heterogeneous across the cortical regions (Fig. 1 and ). MK declined most rapidly in the occipital gyri (Fig. 1 and ). To delineate the MK decrease pattern across the cortical surface, data-driven clustering analysis with NMF method on MK maps was used and revealed four distinctive clusters. Fig. 3 shows four distinctive cortical MK trajectories (Fig. 3) and corresponding anatomical locations of these four clusters on the cortical surface (Fig. 3). Cluster 3 includes the lateral occipital and temporal cortex and is characterized with the highest initial MK values and the fastest MK decrease. The lowest initial MK values with slowest declining rates were found in cluster 1, consisting of the superior frontal, anterior cingulate, and primary sensorimotor cortices.

Fig. 3.

(A) The MK time courses of four clusters identified by NMF. (B) The respective locations of the four clusters on the cortical surface.

Temporal courses of cortical FA (Fig. 2) are different from those of cortical MK. The age-dependent FA with cortical gyri as ROI was best fitted with a biphasic piecewise linear model with significant and differential FA decreases from 31.7 to 36.4 PMW across the brain regions in the first phase and relatively flat FA measurements from 36.4 to 41.7 PMW in most regions in the second phase (Fig. 2). The averaged inflection point was used for more uniform presentation of the fitted lines in the scatter plots in the Fig. 2, whereas the inflection points from biphasic piecewise fitting vary across the brain regions. Regional variation of inflection points can be appreciated with the Fig. 4. The age-related FA decreases were most significant in the frontal gyri (Fig. 2 and ). Four distinctive FA clusters were identified in Fig. 4 with distinctive cortical FA trajectory in Fig. 4 associated with each cluster. The cortical FA trajectories of the four clusters in Fig. 4 also demonstrate the time courses similar to fitted biphasic linear pattern shown in Fig. 2. The MK and FA decrease slopes, correlation coefficients with age, and P values for the 46 cortical gyri can be found in . The 4D spatiotemporal dynamics of cortical MK, RK, AK, FA, MD, RD, and AD measurements are also presented in Movies S1–S7.

Fig. 4.

(A) The FA time courses of four clusters identified by NMF. (B) The respective locations of the four clusters on the cortical surface.

Cortical MK and FA measurements offer different () yet complementary information on the microstructural time courses at individual regions across the cerebral cortex (Figs. 1–4). Two prominent examples of distinctive yet complementary cortical MK and FA time courses at the precentral and postcentral gyrus are shown in Fig. 5, respectively. Significant age-dependent MK decreases (P < 0.05 in ) were found at both primary sensorimotor regions; however, cortical FA remained relatively flat with no significant decreases (P > 0.05 in ) in the first phase at these cortical areas (Fig. 5). Fig. 5 suggests that MK is sensitive to certain underlying cortical microstructural changes that cannot be captured by FA measurements.

Fig. 5.

The significant MK (Left y axis) and nonsignificant first-phase FA (Right y axis) changes at the precentral (primary motor) (A) and postcentral (primary sensory) (B) gyrus. The red triangles are MK measurement, and the blue diamonds are FA measurement of different subjects. Monophasic MK and biphasic linear FA regression lines with correlational coefficient (r), correlation significance (p), and slopes of the fitted linear trend lines are shown. The change rate (slope) is high for MK and low for the first phase of FA in A and B.

Prediction of the Postmenstrual Ages with Cortical DKI and DTI Measurements.

SVR analysis was used to test whether cortical MK and FA could be used to predict the real postmenstrual age. Fig. 6 demonstrates that both cortical MK (Fig. 6; r = 0.63) and FA (Fig. 6; r = 0.92) could predict the postmenstrual age accurately. Both cortical MK and FA measurements are sensitive to brain development and can serve as accurate predictors of postmenstrual age (P < 0.0001), with FA measurements performing slightly better.

Fig. 6.

Age prediction of the preterm neonates with cortical MK (A) and FA (B) measurements.

Synchronized Microstructural Maturation of Cortical Cluster and Associated White Matter.

To examine whether the regional cortical microstructural changes and microstructural changes of associated WM tracts are synchronized, we measured change rates of FA of the WM tracts traced from each of the four MK clusters shown in Fig. 3. Four different WM tracts, identified as cingulum in cingulate gyrus, forceps major, inferior longitudinal fasciculus, and uncinate fasciculus (Fig. 7), were traced from the four MK clusters serving as the seed ROI for DTI-based tractography (see for further details). FA of the WM tracts associated to all four clusters changes significantly with age. As shown in Fig. 7, significant differences between FA change rates of any two of these four WM tracts were found, except between FA change rates of WM tracts associated with cluster 2 and those associated with cluster 3. With cortical MK change rates also plotted in Fig. 7, it can be observed that the cortical MK change rates tend to be synchronized with the associated WM FA change rates. For example, the lowest cortical MK change rate at cluster 1 is coupled with the lowest FA change rate of associated WM tract. In addition, the highest cortical MK change rate at cluster 3 is coupled with the second highest FA change rate of the associated WM tract.

Fig. 7.

Synchronized change rates of cortical MK and associated WM FA. (A) Three-dimensionally reconstructed WM tracts traced from the four clusters shown in Fig. 3. These tracts were identified as cingulum in cingulate gyrus (purple; traced from cluster 1), forceps major (blue; traced from cluster 2), inferior longitudinal fasciculus (green; traced from cluster 3), and uncinate fasciculus (red; traced from cluster 4). (B) Synchronized cortical MK change rates (blue triangles) and associated WM FA change rates (black circle; shown in mean and SD) of the four clusters. Significant differences between FA change rates of any two of these four WM tracts were found, except between FA change rate of WM tract associated with clusters 2 and that associated with cluster 3. *P < 0.05; **P < 0.001.

Discussion

During preterm development from 31 to 42 PMW, differential cortical microstructural profiles and maturational courses were revealed across the entire cortex with MK and FA measurements from high-resolution multishell dMRI. Significant age-related declines of cortical MK were observed for all cortical regions during the studied period with the most rapid changes found in the occipital regions. This study revealed a differential 4D spatiotemporal cytoarchitectural signature inferred by non-Gaussian diffusion barriers inside the cortical plate during the third trimester. Age-dependent decreases of MK are significant in the precentral and postcentral gyri where no significant changes of cortical FA were found, indicating that DKI offers unique microstructural quantification of internal complexity inside the cerebral cortex. This internal complexity may be indicative of important cellular processes such as changes of neuronal density, dendritic arborization, synaptic formation, and myelination of intracortical axons. With data-driven clustering analysis based on cortical MK and FA measurements, the signature of cortical microstructural maturation was revealed. Furthermore, synchronized maturation of microstructures in the cortex and the associated WM was observed. As effective indicators of cortical microstructural architecture, both cortical MK and FA could serve as sensitive imaging markers for age prediction. With noninvasive properties of cortical MK and FA measurements, longitudinal spatiotemporal mapping of cortical MK and FA could delineate trajectories of cortical differentiation across all cortical regions during typical development. The cortical MK or FA measurement could also be potentially used as clinical biomarker by quantifying aberrant deviation in neurological or psychiatric disorders. Taken together, noninvasively delineating the 4D spatiotemporal pattern of cortical MK and FA in the preterm human brain offers unique insights into cellular processes and associated developmental mechanisms during critical development in the third trimester.

Significant Cortical MK and FA Decreases and Associated Cytoarchitectural Processes.

Significant cortical MK and FA decreases during third trimester have been observed in most cortical regions, as shown in Figs. 1–4. The cortical MK measurement may represent a completely different aspect of the microstructural alterations that take place inside the developing cortical plate but cannot be captured by the cortical FA measurement. The relationship of cortical FA and cytoarchitectural processes has been investigated more extensively in the literature. In an immature cortex, water molecules are more likely to diffuse along the radial glia scaffold (2, 7, 8) instead of perpendicular to it, resulting in high cortical FA measurements (17). Radial glia scaffold in a well-defined columnar structure (2, 7, 8) is considered as the neuronal migration pathway from ventricular and subventricular zones to the cortical surface. Consistent with the literature, high FA values were observed in most of the cortex of the youngest preterm brain around 32 PMW (Fig. 2). Cortical maturational processes including dendritic arborization, synaptic formation, and myelination of intracortical axons lead to the disruption of the well-organized radial glia scaffold, causing FA decreases (17, 28). Significant cortical FA decrease has been observed in early brain development in both human (e.g., refs. 16–24) and animal (e.g., refs. 29–35) studies. Cortical MK is thought to be proportional to the complexity of the microstructure (e.g., refs. 39, 42, and 43), and a higher MK value may indicate more diffusion barriers, that is, more densely packed cells, higher dendritic density, or higher axonal density. Significant correlation between MK and neurite density was observed in a recent mouse brain study (44). Higher MK was found in the occipital cortex in the preterm brain at 33 PMW compared with the frontal cortex, suggesting higher neuronal densities in the primary visual areas. During normal development, neurons are usually overproduced and have to compete among themselves for appropriate synaptic sites (45, 46). Those neurons that are unable to connect with any synaptic targets are fated to die. Such a neuronal process has been observed across the species (47). The reduced diffusion barriers reflected by decreased MK might be associated with the decrease in neuronal density caused by cell death, also known as apoptosis in the cortex (5, 6, 48). Significant MK decreases were found across the entire cortex of the preterm brain from 33 to 40 PMW (Figs. 1 and 3), indicating a continuous decrease in diffusion barriers possibly caused by decreases of neuronal density during this period. Significant decreases of cortical MK and FA thus jointly suggest dendritic arborization, synaptic formation, and cell death. This interpretation of cortical MK and FA changes is illustrated with a schematic depiction in .

Heterogeneous 4D Spatiotemporal Patterns of Cortical MK and FA Measurements.

Time courses of cortical differentiation have been predominantly obtained from histological studies (e.g., refs. 2, 4, 7, and 8), offering insights into cytoarchitectural maturation in brain development. However, knowledge of cortical differentiation based on histological studies was limited to relatively localized and small cortical regions. Mapping the entire cortex with DKI and DTI makes it possible to investigate spatial heterogeneity at each PMW as well as heterogeneity of time courses across cortical regions. The cortical MK (Fig. 1) or FA (Fig. 2) at each individual PMW is heterogeneous across the brain regions. During development from 32 to 41 PMW, the decreases of cortical MK (Fig. 1) and FA (Fig. 2) are also asynchronous across the regions. These heterogeneous cortical MK or FA declines were further indicated by distinguished MK (Fig. 3) or FA (Fig. 4) decreasing clusters obtained from a data-driven clustering analysis. Regional heterogeneity in both cortical MK (Fig. 1) and FA (Fig. 2) indicates asynchronous cortical microstructural maturation across different regions. Slower MK and FA decreases (panel p in the Figs. 1 and 2) in insula and relatively steeper decreases of MK and FA in other cortical regions (other panels in Figs. 1 and 2 and Figs. 1 and 2) suggest a developmental gradient moving from insula to surrounding areas along rostral and caudal direction to the frontal and occipital lobe, respectively. This cortical maturation gradient was consistent to the pattern observed in the histopathological findings (49) and recently in a macaque brain study (50). Heterogeneity of decreases in FA has been reported with DTI of preterm brains (18, 22, 24). During a younger age range in the second trimester, spatiotemporally distinctive cortical FA decrease of human fetal brain has also been recorded (21). Heterogeneous cortical FA decreases across the brain regions from 32 to 41 PMW (Figs. 2 and 4 and ) are consistent to these findings. Relatively low FA values at 32 PMW with almost flat FA time courses from 32 to 42 PMW in the primary sensorimotor cortex (Fig. 2 and cluster 1 in Fig. 4) suggest earlier development of primary sensorimotor cortices than higher-order prefrontal cortex (clusters 2 and 4 in Fig. 4). By quantifying cortical microstructures with a more sophisticated diffusion kurtosis model, we hypothesized that regionally distinctive microstructural changes could also be reflected by MK measurements. The distinctive cortical MK decreases from 33 to 40 PMW were quantitatively characterized and confirmed with Figs. 1 and 3 and . MK clusters with largest MK decline rate were found mainly around the primary visual cortex (cluster 3 in Fig. 3). The cortical plate includes cortical regions unique in cytoarchitecture and function. The brain’s diverse functional capacities are modulated by neural circuits (14), the development of which requires coordination of underlying cellular and molecular processes. Formation of unique cortical regions is likely to be associated with distinctive regional dendritic arborization, synaptic formation, and neuronal apoptosis for differential emergence of functional systems during the third trimester. The heterogeneous 4D spatiotemporal pattern of cortical MK and FA reflects these underlying processes and is consistent with findings on the emergence of functional systems with resting-state fMRI. It has been found that rapid increase of functional connectivity during the third trimester occurs in the primary sensorimotor, auditory, and visual networks (51–53), whereas development of higher-order functional networks, such as executive control network, mostly occurs later in life (54). Furthermore, the heterogeneity of the cortical microstructural development across various regions unifies other lines of experimental findings from neuroanatomical to physiological measurements, including synaptic density with histology (4), regional cerebral glucose utilization with PET (55), and regional cerebral blood flow with arterial spin labeling MRI (e.g., ref. 56).

Cortical MK and FA Measures Are Complementary to Each Other and Sensitive to Age in PMW.

Spatiotemporal dynamics of MK and FA are different across brain regions. Cortical MK and FA measures offer complementary microstructural quantifications. For example, a large decrease in MK yet almost no change in FA in the first phase was found in the primary sensorimotor cortex (Fig. 5), indicating MK could be more sensitive to the microstructural changes in tissues with a relatively isotropic Gaussian diffusion. Therefore, MK and FA measurements should be integrated to get a more comprehensive picture of the underlying cortical microstructural changes in the specific early developmental stage in this study.

Either cortical MK or FA predicts the postmenstrual age of preterm and term-born neonates accurately (Fig. 6), indicating that MK or FA can serve as sensitive imaging markers to characterize normal brain maturation. The R value of the cortical FA from prediction analysis is higher than the R value of the cortical MK. This might be related to the fact of DKI being acquired from a smaller group of preterm neonates. As effective indicators of cortical microstructural architecture, cortical MK and FA may also have the potential to serve as quantitative biomarkers besides surface area (57) or functional connectivity (58) to predict aberrant brain maturation in neurodevelopmental disorders.

The FA measurement of the cerebral cortex after 41 PMW is as low as the noise level and no longer sensitive to cortical microstructural dynamics during later infant and childhood development. In contrast to the limited time window during which significant cortical FA changes can be observed, cortical MK quantifies non-Gaussian diffusion widely occurring in cerebral cortex and is sensitive to alterations of cortical internal complexity across the life span (e.g., ref. 40). Of note, the cortical MK measurement may not decrease monotonically in typical development after 41 PMW. Associated with increase or decrease of cortical microstructural complexity in a certain age period, the MK measurement could increase or decrease, respectively. Besides characterizing typical development, cortical MK measurements could be potentially used for quantifying aberrant deviations in pediatric neurological or psychiatric disorders such as autistic spectrum disorder characterized with excessive number of neurons and synapses.

Synchronized Cortical and White Matter Maturation.

By tracing the WM tracts from the cortical MK clusters and measuring rates of FA changes of these WM tracts, a synchronized cortical and WM maturation was observed as demonstrated in Fig. 7. Specifically, the cluster with the more rapid cortical MK decline is associated with a faster FA increase of corresponding WM tract. Despite that statistical significance on synchronized maturation could not be established with limited cluster numbers, observed synchronized maturation (Fig. 7) reconciles with prior work about coherence between cortical thickness reduction and corresponding WM microstructural enhancement in brain development (e.g., refs. 28, 59, and 60). Although the underlying mechanisms are not completely clear, this synchronized cortical and WM maturation is likely related to the formation and growth of brain circuits. Most WM axons connecting different brain regions are projected from the neurons in the cerebral cortex, and mature coherently with their associated neurons in development (61). Previous studies on cortical thickness (e.g., refs. 62 and 63) also suggested that cortical regions that grow together are perhaps mediated by the development of their underlying axonal connections.

Methodological Considerations and Limitations.

Despite that cortical MK or FA measurements could effectively characterize dynamics of cortical microstructure in typical and atypical brain development, comprehensive validation against neuropathological measurements such as neuronal density and dendritic arborization is warranted. Precise interpretation of underlying cytoarchitectural processes (e.g., model suggested by ) can be established only after the validation. Similarly, the clinical potential of detecting cortical microstructural abnormality by surveying the entire cerebral cortex through 5–10 min of dMRI could be fully exerted only after validation. Preterm birth has been widely recognized as one of the risk factors that might affect normal brain development. Nevertheless, MRI examinations of preterm infants have been predominantly used to understand cortical microstructure development during the third trimester (e.g., refs. 17 and 21–24). Exposure to the extrauterine environment could be relatively subtle compared with dramatic developmental factors during the third trimester (64, 65). Although none of the included neonates were clinically referred and all had been rigorously examined to avoid apparent abnormality, preterm effects could still affect the temporal courses in the present study. It is noteworthy that all neonates (47 out of 76 neonates) who had undergone the 2-y follow-up neurodevelopmental tests demonstrated no apparent neurological disorder or neurodevelopmental delay. MD was also used to measure cellularity and can be relatively more conveniently obtained with single-shell dMRI. The MK measurement, on the other hand, demonstrates unique properties in quantifying neuronal density and cellularity, the latter of which is explored also in other research areas such as differentiating brain tumor grades (e.g., refs. 66 and 67). Despite a higher dMRI resolution of 1.5 × 1.5 × 1.6 mm3 being used, partial volume effects (68) cannot be ignored for measuring cortical MK and FA. With the analysis pipeline (22, 24) demonstrated in , the cortical skeleton was obtained to measure the cortical microstructure at the center or core of the cortical plate. In this way, partial volume effects could be maximally alleviated and measurement accuracy could be enhanced. Last, all data presented in this study are limited by its cross-sectional study design. Longitudinal data in the future are needed for delineating trajectories of cortical microstructural maturation.

Conclusion

This study represents a comprehensive characterization of cortical microstructural changes in preterm brains using both ROI-based and data-driven analyses with combined DKI and DTI datasets. Different from DTI-derived FA measurements, DKI-derived MK measurements provide extra information of cortical architecture by quantifying non-Gaussian diffusion and are especially useful in characterizing the microstructure of cortical tissue with isotropic diffusion. Application of DKI to measure the restricted water diffusion properties of the cerebral cortex provides fresh insight into the dynamic cortical microstructural signature of critical developmental stages in not only the third trimester but also other postnatal developmental periods. Furthermore, cortical MK measurements could be used as potential clinical biomarkers for atypical brain development. The exact neuroanatomical underpinning of cortical MK is not completely known; however, we suspect that MK is more sensitive to neuronal density, which cannot be detected with DTI in the developing brain of the studied age period.

Materials and Methods

Participants.

The study was approved by the Institutional Review Board at the University of Texas Southwestern Medical Center. Neonates were recruited from the Parkland Hospital and scanned at Children’s Medical Center at Dallas. In total, 107 preterm and term neonates were consented, of which 89 neonates were scanned and 76 were included in the study (47 males/29 females; postmenstrual ages, 31.5–41.7 PMW). Written consent was obtained from the parent of every participating neonate. The demographic information, exclusion criteria, and other details of the included subjects are described in .

MRI Data Acquisition.

All neonates were well fed before scanning and remained asleep during scan. No sedation was used. Multishell diffusion-weighted images (DWIs) with b = 1,000 and 1,600 s/mm2 were acquired from a 3T Philips Achieva system (Philips Healthcare) with 1.5 × 1.5 × 1.6 mm3 imaging voxel size. Please see for acquisition details.

Diffusion Kurtosis and Tensor Fitting.

All b = 1,000 s/mm2 DWIs were registered to the nondiffusion-weighted image (i.e., b0 image). The tensor fitting was conducted with DWI of b = 1,000 s/mm2 after motion and distortion correction to obtain the DTI maps including FA, MD, RD, and AD maps. After DWIs of b = 1,600 s/mm2 were corrected for motion and distortion, kurtosis was fitted using in-house software written in MATLAB (2015a; Mathworks) to obtain the DKI maps including MK, RK, and AK maps.

Measurement of DKI and DTI Metrics on the Parcellated Cortical Skeleton.

To alleviate partial volume effects, the cortical DKI- and DTI-derived parameters were measured on the cortical skeleton, that is, the center of the cortical plate, as shown in . Details of measuring DKI and DTI metrics on the parcellated cortical skeletons can be found in .

Clustering Analyses.

To characterize the age-related cortical MK or FA change patterns, we employed the NMF method, an unsupervised, free from regional hypothesis, and data-driven multivariate method recently used for analysis of structural neuroimaging data (69). The approach is part of the NMF family of methods (70) that aim to factorize a tall nonnegative data matrix (constructed by arraying each data sample per column, with a data sample referring to a vectorized MK or FA at cortical skeleton voxels) into two nonnegative matrices:

Details of the clustering analysis can be found in .

SVR Analyses for Age Prediction.

To determine whether FA and MK in different brain regions could serve as a biomarker for brain age prediction, we performed pattern analysis using the SVR algorithm implemented in LIBSVM (A Library for Support Vector Machines 3.21; https://www.csie.ntu.edu.tw/∼cjlin/libsvm/). Postmenstrual age was used as the training measure, and FA or MK values from different brain regions were used as features in the SVR models with linear kernel function. Leave-one-out cross-validation was adopted to evaluate the performance of the SVR models. Pearson correlation coefficient between the actual and predicted ages was computed to assess the prediction accuracy.

Analyses of Cortical MK or FA Time Courses of All Cortical Gyri.

The analyses below were conducted using R statistical software, version 3.0.2 (www.r-project.org). MK time courses demonstrated monophasic linear changes. A monophasic linear model was confirmed with ANOVA for MK time courses. The random-effects linear model was used to fit the cortical MK of each gyrus and age for all subjects:

where t is the postmenstrual age of the jth subject with j from 1 to 26; α and β are the intercepts and slopes for MK of the th cortical gyrus, respectively, with i from 1 to 46 (covering all cortical gyri in the left and right hemisphere of the brain except only six regions around hippocampus and medial temporal area due to poor registration quality for these six regions); ε is the error term. Familywise error rate was controlled with the FDR correction. Values of β1,, obtained in Eq. , were estimated MK change rates.

FA time courses demonstrated nonlinear changes of FA over time. We then tested the random-effects model based on first- and second-order polynomials and a biphasic piecewise linear model. The biphasic piecewise linear model with an inflection point at 36 PMW was selected based on ANOVA analysis. The inflection point was obtained by fitting FA of all cortical gyri in left and right hemisphere to a piecewise linear function and averaging the inflection points of the best fit. The biphasic linear model below was used to fit the cortical FA of each gyrus and age during 31–36 and 36–42 PMW:

where t is the postmenstrual age of the jth subject with j from 1 to 76; c and f are the intercepts and β2,, and β3, are the slopes for FA of the th cortical gyrus, with i from 1 to 46; and are the error term for cortical FA measurement from 31 to 36 PMW and from 36 to 42 PMW, respectively. The analysis was controlled for multiple comparisons using a FDR correction. In addition, an ANOVA F test was performed to prove that FA and MK time courses at each cortical gyrus were not redundant.

Test of the Correlation Between Cortical MK and FA of All Cortical Gyri.

The null hypothesis is that the correlation coefficient between MK and FA from a given cortical gyrus is not significantly different from zero. The tests were conducted with partial correlation analysis between MK and FA measurements of all cortical gyri controlling age effects. Rejection of the null hypothesis indicates a significant correlation between MK and FA at a given cortical gyrus. Familywise error rate was controlled with the FDR correction.

Test of Group Differences Between WM Tracts Traced from the Cortex.

Cortical clusters were dilated into WM and used as the tracing seeds to initialize DTI-based tractography. FA of the WM tract associated to each cluster was then measured. The details of tractography and WM FA measurement are described in . Group differences of the FA change rates between the WM tracts were analyzed using unpaired t test. Familywise error rate was controlled with the Bonferroni correction.