Both cortical folding and structural connection patterns are more elaborated during the evolution of primate neocortex. For instance, cortical gyral shapes and structural connection patterns in humans are more complex and variable than those in chimpanzees and macaques. However, the intrinsic relationship between gyral folding and structural connection and their coevolution patterns across primates remain unclear. Here, our qualitative and quantitative analyses of in vivo diffusion tensor imaging (DTI) and structural magnetic resonance imaging (MRI) data consistently demonstrate that structural fiber connection patterns closely follow gyral folding patterns in the direction “tangent” to the cortical sphere, and this close relationship is well preserved in the neocortices of macaque, chimpanzee, and human brains, despite the progressively increasing complexity and variability of cortical folding and structural connection patterns. The findings suggest a hypothesis that a common axonal fiber pushing mechanism sculpts the curved patterns of gyri in the tangent direction during primate brain evolution. Our DTI/MRI data analysis provides novel insights into the structural architecture of primate brains, a new viewpoint of the relationship between cortical morphology and connection, and a basis for future elucidation of the functional implication of coevolution of cortical folding and structural connection patterns.
Brain evolution has been an intriguing research topic for centuries (Zilles et al. 1988, 1989; Deacon 1990; Schoenemann 2006; Rakic 2009; Rogers et al. 2010). Over the past decades, there has been increasing interest in applying noninvasive neuroimaging methods (e.g. Rilling and Insel 1999; Rilling and Seligman 2002; Woods et al. 2011) to study brain evolution because 3D morphological and anatomical attributes of mammalian brains can be effectively visualized and measured. In particular, advanced neuroimaging techniques such as diffusion tensor imaging (DTI) (Mori 2006) and functional magnetic resonance imaging (fMRI) (Friston 2009) have been recently used to study the connection and function of primate brains. For instance, a recent comparative DTI study (Rilling et al. 2008) reported a prominent temporal lobe projection of the human arcuate fasciculus whose counterpart is significantly smaller or absent in nonhuman primates (Rilling et al. 2008). This work suggests that DTI is a powerful approach in revealing interesting structural connectivity patterns of brain evolution and that DTI-derived fiber bundles can be potentially linked to brain function (Rilling et al. 2008). In another recent study (Vincent et al. 2007), fMRI was employed to examine the intrinsic functional architecture in anesthetized macaque brains and resting-state human brains. It was found that 4 functional networks including the oculomotor, somatomotor, visual, and default-mode systems are well preserved across macaque and human brains (Vincent et al. 2007). These studies suggest that neuroimaging such as DTI and fMRI is a powerful approach to studying brain evolution.
In general, during the evolution of the neocortex of primate brain, the variability and complexity of gyral shape patterns are gradually pronounced, as illustrated in Figure 1a. Meanwhile, the axonal fiber connection patterns are much more elaborated during the evolution of primate cortex (Rilling et al. 2008; Krubitzer 2009). As an example, Figure 1b illustrates the DTI-derived fibers emanating from the postcentral gyrus (Fig. 1a) in 3 primate brains. It is evident that there are more DTI-derived axonal fibers connecting to the frontal and occipital lobes in the human brain than those in the macaque and chimpanzee brains. Based on the patterns in Figure 1a,b and existing neuroscience knowledge of the relationship between corticogenesis and axonal wiring (Rakic 1988; Welker 1990; Van Essen 1997; Sur and Rubenstein 2005; Rash and Grove 2006; Nie et al. 2011), we hypothesize that gyral folding and structural fiber connection patterns “coevolve” in primate brains; we test this hypothesis via quantitative modeling and analysis of cortical surfaces reconstructed from MRI data and axonal fibers reconstructed from DTI data. Specifically, we developed novel computational methods to quantitatively measure the complexities of cortical gyral shape and fiber connection patterns and performed statistical correlation analyses for 3 primates. Our results derived from quantitative modeling clearly showed the increasing complexities of gyral folding and fiber connection patterns in all of the 3 primates of macaques, chimpanzees, and humans. In particular, the outcome of our analysis showed a well-preserved, strong, positive correlation between the complexities of gyral folding and fiber connection patterns across the 3 primates, suggesting the coevolution of gyral folding and structural connection patterns in primates.
We further explored the possible underlying mechanism that can explain the clear coevolution of gyral folding and structural connections patterns across the 3 primate neocortices. Our previous studies already demonstrated that significantly denser DTI-derived axonal fibers are connected to the gyral regions than to the sulcal regions (Nie et al. 2011), as shown in Figure 2a. Our extensive quantitative DTI data analysis in human, chimpanzee, and macaque brains replicated this observation (Nie et al. 2011), and based on this observation, we hypothesized that axonal fiber pushing (blue arrows in Fig. 2a) is a mechanism that induces or regulates the “radial” formation of gyral regions (dashed yellow curve in Fig. 2a) in the cerebral cortex (Nie et al. 2011). This hypothesis is specifically supported by the biomechanical experiment in the fruit fly and by the cortex growth simulation experiment in Nie et al. (2011). In this paper, we report “another” interesting finding: the gyral folds in the “tangent” direction of the cortical sphere (dashed green curve in Fig. 2b) are connected by dense axonal fibers that orient perpendicularly to the gyral crest curves, as demonstrated in Figure 2b (red arrows in Fig. 2b). We replicated this close relationship between tangent gyral folds and axonal fiber orientation directions by extensive qualitative and quantitative analyses in all of the 3 primate brains in this paper. Thus, we are inspired to extend the hypothesis in Nie et al. (2011) to that of axonal fiber pushing. This can also be a mechanism that induces or regulates the curved folding patterns in the direction tangent to the gyral crest line, as illustrated in Figure 2b. Our results show that this hypothesis is supported by a series of experiments in the 3 primates based on DTI/MRI data. Therefore, the work reported in this paper, together with our previous studies in Nie et al. (2011), suggests a novel axonal fiber pushing mechanism of cortical gyral folding. This hypothesis can well explain not only the folding patterns of gyri along the radial direction (dashed yellow curve in Fig. 2a) (Nie et al. 2011) of the roughly spherical-shaped cortex (e.g. Fig. 2a), but also the tangent folding of gyral crests on the cortical surface (e.g. Fig. 2b).
The major contributions of this work are in the following 2 aspects. First, we developed, evaluated, and applied a novel pipeline of computational approaches to quantify cortical gyral folding and fiber orientation patterns in order to enable and facilitate the joint analysis of their coevolution patterns. As a result, the correlation between fiber orientation patterns and tangent gyral folding patterns can be quantitatively assessed via a joint representation framework. These novel approaches have been shown to be effective in describing the 3 types of primate brains and have revealed a well-preserved, strong, positive correlation between the complexities of gyral folding and fiber connection patterns across the 3 primates. Secondly, we revealed a novel finding that tangent gyral folds are connected by dense DTI-derived fibers orienting perpendicularly to the gyral crest curves and replicated and confirmed across the 3 species of humans, chimpanzees, and macaques we studied. This interesting finding provides novel insights into the underlying mechanism of coevolution patterns of gyral folding and fiber connection in primate brains: a common axonal fiber pushing mechanism that sculpts the tangent folding patterns of cortical gyri during primate brain evolution.
Materials and Methods
Data Acquisition and Preprocessing
Human Brain Imaging
Nineteen healthy volunteers were scanned in a GE 3 T Signa MRI system (GE Healthcare, Milwaukee, WI, USA) using an 8-channel head coil at the Bioimaging Research Center of the University of Georgia under Institutional Review Board approval. DTI data were acquired using a spatial resolution of 2 × 2 × 2 mm3, matrix size = 128 × 128, 60 slices, repetition time (TR) = 15.5 s, echo time (TE) = min-full, b-value = 1000 with 30 DWI gradient directions and 3 B0 volumes acquired. All slices were aligned to the AC-PC line. T1-weighted MRI images were acquired using a fast spoiled gradient-recalled echo protocol: TE = min-full, TR = 7.5 ms, ﬂip angle = 20°, 154 axial slices, slice thickness = 1.2 mm, and field of view (FOV) = 256 × 256 mm. For the anatomical data, preprocessing included brain skull removal, gray matter (GM), and white matter (WM) tissue segmentations (Liu et al. 2007). Then, the GM/WM cortical surfaces were generated (Liu et al. 2008). For the DTI data, preprocessing included brain skull removal, motion correction, and eddy current correction. Fiber tracts were generated on the basis of the DTI data using a method detailed elsewhere (Zhang, Guo, Li et al. 2010).
Chimpanzees Brain Imaging
All chimpanzees were members of a colony in the Yerkes National Primate Research Center (YNPRC) at Emory University in Atlanta, GA, USA. All imaging studies were conducted under the IACUC approval of Emory University. MRI and DTI scans were obtained from adult female chimpanzees. Twenty-four cases of MRI/DTI scans were used for this study. Prior to MRI/DTI scanning, the chimpanzee subjects were immobilized with ketamine injections (2–6 mg/kg, i.m.) and subsequently anesthetized with an intravenous propofol drip (10 mg/kg/h), following standard veterinary procedures used at the YNPRC. The subjects remained sedated for the duration of the scans as well as the time needed for transportation between their home cage and the scanner location. After completing the MRI/DTI scans, the chimpanzees were temporarily housed in a single cage for 6–12 h to allow the effects of anesthesia to wear off before being returned to their home cage and cage mates. The veterinary staff and research staff assessed the general well-being (i.e. activity and food intake) of the chimpanzees twice daily after the scan for possible distress associated with anesthetic accesses.
Both anatomical MRI and DTI scans were acquired on a Siemens 3 T Trio scanner (Siemens Medical System, Malvern, PA, USA), with a standard birdcage coil. Foam cushions and elastic straps were used to minimize head motion. DTI data were collected with 2 different imaging protocols: multi-shot double spin-echo echo planar imaging (MS-EPI) sequence and single-shot double spin-echo EPI (SS-EPI) sequence. For MS-EPI, a dual spin-echo technique combined with bipolar gradients was used to minimize eddy current effects. The parameters used for diffusion data acquisition were as follows: diffusion-weighting gradients applied in 60 directions with a b-value of 1000 s/mm2, TR/TE of 5740/91 ms, FOV of 230 × 230 mm2, matrix size of 128 × 128, resolution of 1.8 × 1.8 × 1.8 mm2, 41 slices with no gap, covering the whole brain. For each diffusion direction, 2 diffusion-weighted images were acquired, each with 1 of the possible left–right phase-encoding directions and 2 averages, allowing for correction of susceptibility-related distortion using an algorithm described in the literature (Andersson et al. 2003). For each average of diffusion-weighted images, 6 images without diffusion weighting (b = 0 s/mm2) were also acquired with matching imaging parameters. The total MRI/DTI scan time was approximately 50 min. For SS-EPI, an SS-EPI sequence with reduced FOV in the phase-encoding direction and a partial Fourier imaging technique were used to scan the chimpanzees. The parameters used for diffusion data acquisition were as follows: diffusion-weighting gradients applied in 60 directions with a b-value of 1000 s/mm2, TR/TE of 5900/84 ms, FOV of 129 × 230 mm2, matrix size of 72 × 128, resolution of 1.8 × 1.8 × 1.8 mm3, 41 slices with no gap, covering the whole brain. For each average of diffusion-weighted images, 1 image without diffusion weighting (b = 0 s/mm2) was also acquired with matching imaging parameters. The total MRI/DTI scan time was 53.4 min. High-resolution T1-weighted MRI images were acquired with a 3D magnetization-prepared rapid gradient echo (MP-RAGE) sequence for all participants. For subjects scanned using the MS-EPI sequence, the scan protocol, optimized at 3 T, used a TR/inversion time/TE of 2400/1100/4.13 ms, a flip angle of 8°, a volume of view of 256 × 256 × 154 mm3, a matrix of 256 × 256 × 192, and a resolution of 1.0 × 1.0 × 0.8 mm3, with 2 averages. The total T1 scan time was approximately 20 min. For subjects scanned using SS-EPI, the scan protocol is similar, despite that the resolution is 0.8 × 0.8 × 0.8 mm3 isotropic and volume of view is 154 × 154 × 154 mm3. Preprocessing steps were similar to those used in the processing of human MRI/DTI data.
Macaques Brain Imaging
All macaques were members of a colony at YNPRC. All MRI and DTI scans were conducted under the IACUC approval of Emory University. MRI/DTI data from 8 subjects were used for this study. Prior to scanning, the subjects were immobilized with ketamine injections (2–6 mg/kg, i.m.) and subsequently anesthetized with an intravenous propofol drip (10 mg/kg/h), following standard veterinary procedures used at the YNPRC. The macaques remained sedated for the duration of the scans as well as the time needed for transportation between their home cage and the scanner location. After completing the MRI scan, the macaques were temporarily housed in a single cage for 6–12 h to allow the effects of anesthesia to wear off before being returned to their home cage and cage mates. The veterinary staff and research staff observed the general well-being (i.e. activity and food intake) of the macaques twice daily after the scan for possible distress associated with anesthetic accesses.
Both anatomical MRI and DTI scans were performed on a Siemens 3T Trio scanner (Siemens Medical System), with a standard knee coil. Foam cushions and elastic straps were used to minimize head motion. A specially designed holding device was used to stabilize macaque's head during scanning, with 2 plastic screws anchoring in the macaque's ear canals tightly. High-resolution T1-weighted MRI images were acquired with a 3D MP-RAGE sequence for all participants. The scan protocol used a TR/inversion time/TE of 2500/950/3.49 ms, a flip angle of 8°, a volume of view of 128 × 128 × 96 mm3, a matrix of 256 × 256 × 192, and a resolution of 0.5 × 0.5 × 0.5 mm3, with 3 averages. The total T1 scan time was approximately 33 min. Diffusion MRI data were collected with a diffusion-weighted, multi-shot (3 segments), spin-echo EPI sequence. A dual spin-echo technique combined with bipolar gradients was used to minimize eddy current effects. The parameters used for diffusion data acquisition were as follows: diffusion-weighting gradients applied in 60 directions with a b-value of 1000 s/mm2, TR/TE of 6970/104 ms, FOV of 141 × 141 mm2, matrix size of 128 × 128, resolution of 1.1 × 1.1 × 1.1 mm3, 41 slices with no gap, covering the whole brain. Similar to that for chimpanzees, diffusion-weighted images were acquired with phase-encoding directions of opposite polarity (left–right), each with 4 averages, to correct for susceptibility-induced distortion. For each average of diffusion-weighted images, 5 images without diffusion weighting (b = 0 s/mm2) were also acquired with matching imaging parameters. The total MRI and DTI scan time was approximately 90 min. Preprocessing steps were similar to those used in the processing of human data.
Joint Representation of Gyral Shape and Structural Fiber Connection
Since the major objective of this work is to assess the coevolution of gyral folding and structural fiber connection patterns, we employed a joint representation of gyral shape and structural fiber connection in this work. Specifically, our extensive recent studies (Li, Guo, Li et al. 2010; Li, Guo, Nie et al. 2010; Nie et al. 2011) have shown that DTI-derived WM fiber endpoints closely follow the gyral folding patterns, as highlighted by the orange arrows in Figure 2f–h. This observation has been replicated in all of the DTI data sets of human, chimpanzee, and macaque brains we analyzed, indicating that the gyral shape pattern and the structural connection pattern have complementary information that can facilitate quantitative analysis of their coevolution brain patterns across primate species. Therefore, we performed whole-brain fiber tracking (via MEDINRIA, http://www-sop.inria.fr/asclepios/software/MedINRIA/) and aligned cortical surfaces reconstructed from T1-weighted MRI image to DTI space with linear transform matrix obtained from the volumetric registration (via FSL FLIRT). After linear alignment, the misalignment between anatomical space and DTI space caused by distortion is relative small, as shown in Supplementary Figure S1. Specifically, the average surface distances caused by the misalignment are 2.21 mm for humans, 1.76 mm for chimpanzees, and 1.35 mm for macaques, respectively. Then, we performed segmentation of 4 major gyri, including precentral gyrus, postcentral gyrus, superior temporal gyrus, and superior frontal gyrus in each hemisphere, from the above-reconstructed cortical surfaces based on experts’ interactive labeling. The open-source software Paraview (http://www.paraview.org/) was used for this interactive gyrus segmentation. Supplementary Figure S2a shows an example of the 4 segmented gyri in 3 primate brains. After each gyrus is segmented from the cortical surface in the DTI image space, the emanating fibers were extracted from the results of whole-brain streamline tractography via a similar method detailed elsewhere (Zhu, Li, Faraco, Deng, Zhang et al. 2011; Zhu, Li, Faraco, Deng, Zhu et al. 2011). As a result, each gyral shape and its structural fiber connections are colocalized and jointly represented on the “same” gyrus in the DTI space, as demonstrated in Supplementary Figure S2b. It should be noted that the joint representation of gyral folding and structural connection patterns treats the brain as a whole via applying in vivo MRI and DTI to study the brain's GM and WM architecture. By exploring the complementary information provided by gyral folding and connectivity patterns, this joint representation methodology offers advantages over other approaches that treat GM and WM separately (Li, Guo, Li et al. 2010; Li, Guo, Nie et al. 2010; Nie et al. 2011) and could provide important insights into the regularity and variability of the brain structure.
Quantification of Fiber Orientation Patterns
To quantitatively measure the relationship between fiber orientation patterns and gyral folding patterns based on the above joint representation, a coordinate system needs to be set up. However, the orientation of a gyrus can vary in different ways; for example, the precentral gyrus and the postcentral gyrus roughly follow a superior–inferior orientation, whereas the superior temporal gyrus and the superior frontal gyrus have a rough anterior–posterior orientation. This makes it very difficult to set up a global coordinate system along the tangent direction for each gyrus. Instead, we defined a 3D local Cartesian coordinate system for each gyrus, as shown in Figure 3a,d. In each local coordinate system, the average normal direction of vertices on the gyral crest curves, which was obtained by thresholding principal curvatures on the cortical surface, is treated as the Z-axis (Fig. 3b). Let p and q denote 2 points manually selected at each end of a gyrus, Y-axis is defined by, and X-axis is obtained by the right-hand rule (Fig. 3d). By default, the Y-axis is identified as the tangent direction of a gyrus in this paper.
After the local coordinate system has been constructed, we projected the DTI-derived fibers that intersect with the gyrus to the Y-Z plane. The fiber orientation pattern is defined by the angle between the Z-axis and the fiber's main direction vector on the end that intersects with the gyrus, as denoted by in Figure 3e. As shown in Figure 3e, if the projection of the main direction vector on Y-axis is positive, is positive. Otherwise, it is negative. It should be noted that this paper mainly focusses on the correlation between fiber orientation and gyral folding in the tangent direction, that is, along the Y-axis. Although we may lose certain fiber orientation information in other directions by projecting 3D fibers to the 2D Y-Z plane, this approach enables us to jointly examine the correlation between fiber orientation pattern and gyral folding pattern defined in the next section in the same local coordinate system.
To facilitate quantitative measurements of the correlation between the fiber orientation pattern and the gyral folding pattern, we mapped the fiber orientation pattern (an example shown in Fig. 3f) to the gyrus that it intersects with, as shown in Figure 3g. Subsequently, the fiber orientation pattern is color-coded on the gyral crest for both visualization and quantification. It should be noted that the selection of gyral endpoints p and q may influence the direction of X- and Y-axes, and thus the fiber orientation patterns may vary slightly due to this influence. However, our extensive experience with this approach showed that this type of influence on our final analysis result is relatively small. For instance, Figure 3i(i and ii) illustrates one example of how different selections of endpoints for the local coordinate system construction influence the computed fiber orientation patterns. It is evident that the fiber orientation maps in both cases are similar. But if we revert (p, q) as shown in Figure 3i(iii), the fiber orientation pattern is flipped due to the flipping of X- and Y-axes. Thus, in order to unify the fiber orientation pattern and to make the results reliable, we selected the (p, q) pair on 2 ends of a gyral ridge in all primate brains by following a superior–posterior to inferior–anterior rule, as shown in Figure 3h.
Quantification of Gyral Folding Patterns
The gyral folding pattern is quantified by the shape of its crest curve on the gyral ridge. We used the fast marching method in Li, Guo, Li et al. (2010) and Li, Guo, Nie et al. (2010) to find the shortest geometric path weighted by , which is inversely proportional to the maximum principal curvature of the surface between 2 manually selected gyral endpoints p and q on the surface (Fig. 3j,k) as follows:is a section of line on the surface of gyrus such that and are connected to each other. Denoting as the projection of on the X-Y plane of the gyral local coordinate system defined in the previous section, the shape pattern of each vertex on curve is defined by the average determinant of the orientation matrix in each vertex's neighborhood on (Fig. 3n): the orientation matrix of vertex of . Denote as the coordinate of vertex : 3o, we showed the folding pattern values of the simulated curve of sine functions with linearly increasing amplitudes. We chose the average determinant of the orientation matrix as the folding pattern descriptor of a gyrus for the following reasons. First, it includes the orientation information of a gyral crest curve, for example, negative value for clockwise or positive value otherwise, as illustrated by the blue arrows in Figure 3o. Secondly, it includes local curvature information. The larger the absolute value, the larger the curvature of the curve. In other words, the higher the average absolute value of this descriptor is, the more curved the gyral folding pattern will be. Thirdly, it is robust to noise. As we are considering a local neighborhood when computing this value, it is less influenced by local noises. With these properties, we use the determinant of the orientation matrix on the gyral crest curve as a cortical folding pattern descriptor and compare its value with the fiber orientation patterns defined in Figure 3e.
The Correlation Between Fiber Orientation Patterns and Gyral Folding Patterns
To measure the correlation between fiber orientation patterns and gyral folding patterns, we used the gyral crest curve as a common reference along which we obtained the fiber orientation patterns (see Quantification of Fiber Orientation Patterns) and computed its Pearson's correlation with the gyral folding patterns along the same gyral crest curve (see Quantification of Gyral Folding Patterns). The computational pipeline is illustrated in Figure 4a–d. Before comparison, the patterns in 2 ends of a gyrus were discarded since the defined folding pattern cannot be reliably computed at the end of a curve.
It should be noted that as the defined local coordinate system is relatively constant for each gyrus, for certain gyri (e.g. the superior temporal gyrus in this paper), the normal directions on their top ridges may rotate from one end to the other. This could result in a global increasing or decreasing of fiber orientation patterns as shown in Figure 4e,f, which may cause inaccuracy when comparing them with the defined gyral folding patterns. As our gyral folding pattern descriptor contains only local shape information, the global trend of fiber orientation pattern has to be removed before computing their correlations with the gyral folding pattern. Our previous studies demonstrated that this type of trend can be effectively removed by a linear regression model, as shown in Figure 4g,h. It is evident that after this de-trending, the fiber orientation pattern and the folding pattern are closely correlated along the gyrus, as highlighted by the purple arrows in Figure 4h,i. In this paper, our experiments showed that only the superior temporal gyrus needs this de-trending preprocessing.
Based on the neuroimaging data sets and computational approaches mentioned earlier, we investigated the coevolution of gyral folding and structural fiber connection patterns in the 3 primate species in 3 manners, allowing us to visually and quantitatively assess the relationship between gyral folding in the tangent direction and fiber orientations patterns. The specifics of these investigations and the results are described subsequently.
Coevolution of Cortical Folding and Structural Fiber Connection Patterns
We compared the complexities of shape patterns of the precentral and postcentral gyri (Supplementary Fig. S2a) and their fiber connection patterns in 3 primate species. The average folding pattern description of the gyral crest defined in equation (3) was used as the gyral shape complexity. The fiber connection pattern complexity was quantified by the diversity of fiber endpoint destinations, that is, the ratio of endpoints of fiber connections to the frontal and occipital lobes (purple and blue fiber bundles in Fig. 1b) to the total number of fiber connections. Our rationale is that the precentral and postcentral gyri in the macaque brains have much less fiber connections to the frontal and occipital lobes, whereas the precentral and postcentral gyri in chimpanzee and human brains have more and more elaborated fiber connections to the frontal and occipital lobes (Fig. 1b). Therefore, we postulated that the percentage of fibers connecting the primary motor and sensory areas to the frontal and occipital lobes is a good indicator of structural connectivity complexity for these areas and plotted the complexities of both gyral folding patterns and fiber connection patterns in a 2D diagram for each of the 4 gyri on both hemispheres in Figure 5. For instance, Figure 5a shows the complexities of the left precentral gyrus folding patterns (vertical axis) and fiber connection patterns (horizontal axis) for macaque (triangles), chimpanzee (rectangles), and human (diamonds) brains. It is evident that there is an approximately linear increasing trend for both complexities of gyral folding patterns and fiber connection patterns, likely reflecting the coevolution of cortical gyral shape and fiber connection patterns. In Figure 5a, the trend is regressed by a black line. This linear trend can also be observed in the right precentral gyri in Figure 5b and in the postcentral gyri on both hemispheres, as shown in Figure 5c,d.
Based on the results shown in Figure 5, we conjecture that the complexity of the gyral folding pattern coevolutes with that of the fiber connection pattern in a close positive relationship, suggesting a deep-rooted evolutionarily preserved mechanism underlying the structural fiber connection pattern and the cortical gyral folding pattern. In further analysis described subsequently, we further investigated the possible underlying mechanisms that underlie the clear coevolution of cortical folding and structural connections patterns across the 3 primate neocortices and tested the hypothesis that the gyral folding pattern is induced or regulated by axonal fiber pushing.
Visual Examination of the Relationships Between Gyral Folding and Fiber Orientation Patterns
We examined the relationship between gyral folding and fiber orientation patterns both qualitatively and quantitatively. First, we visually assessed the fiber orientation pattern and the gyral shape pattern for 4 segmented gyri in the right hemisphere in Figure 6. Specifically, Figure 6a–d illustrates the precentral gyrus, postcentral gyrus, superior frontal gyrus, and superior temporal gyrus, respectively. In each part figure, the first column demonstrates the joint visualization of gyral shape and the emanating DTI-derived fibers, the second column shows the color-coded fiber orientation in the local coordinate system, and the last column illustrates the color-coded gyral shape pattern quantified by the descriptor (see Quantification of Gyral Folding Patterns). Each row shows the result for one example of macaque, chimpanzee, and human brains, respectively. By visual examination of these 4 gyri, it is clear that there is a close relationship between the fiber orientation and cortical folding along the gyral crest curve. For instance, in Figure 6a, the blue arrows highlight negative fiber orientations in the second column, which colocalize well with the negative folding patterns in gyral shapes, as highlighted by the blue arrows in the third column as well. Meanwhile, the positive fiber orientations highlighted by the orange arrows in the second column also colocalize well with the positive folding patterns in the third column, which are pointed to by the orange arrows. For comparison, the corresponding blue and orange arrows are also provided in the first column of joint visualizations in Figure 6. These similar positive correlations between fiber orientation patterns and gyral folding patterns have been replicated in other gyri, as shown in Figure 6b–d. Our extensive observations in all of the 3 primate cases we studied overwhelmingly demonstrate that there is an evolutionarily preserved positive correlation between the cortical folding and fiber connection orientation patterns, which support the concept of coevolution of cortical folding and fiber connection patterns.
Similarly, Supplementary Figure S3 shows the visualization of fiber connection orientation and cortical gryal folding patterns for 4 gyri in the left hemisphere. Again, the positive correlation between fiber orientation and gyral folding patterns is evident, as highlighted by the blue and orange arrows in Supplementary Figure S3a–d. By joint visualizations shown in Figure 6 and Supplementary Figure S3, it is clear that for all 8 gyri on both hemispheres of 3 primates, we observed that there exists a positive relationship between fiber orientations and gyral folding patterns. Therefore, we conjecture that axonal fiber pushing forces induce or regulate the tangent folding patterns of gyral crest curves, as illustrated in Figure 2b. This conjecture can not only explain the close positive relationship between fiber orientations and folding patterns as shown in Figure 6 and Supplementary Figure S3, but also well interpret the increasing complexities of cortical gyral folding patterns across primates. For instance, the precentral gyrus in Figure 6a exhibits increasing complexity of shapes across primates, for example, 3 clear positive folds and 3 clear negative folds in humans, 1 clear negative fold and 2 positive folds in chimpanzees, and no clear folds in macaques. These increasing complexities of cortical gyral folding are exactly mirrored by the increasing complexities of fiber connection patterns, for example, the denser and complex fiber connections in humans and much simpler and sparser fiber connections in chimpanzees and macaques, as shown by the first column in Figure 6a. From a biomechanical perspective, it is reasonable to hypothesize that the increasingly complex fibers “push” the gyri into increasingly complex shape patterns, such as more and more positive and negative folds during the evolution of primate brains.
Similarly, axonal fiber pushing of tangent gyral folding can well explain other gyri on both hemispheres of 3 primates, as shown in Figure 6 and Supplementary Figure S3, for example, the postcentral gyrus in Figure 6b and the precentral gyrus on the left hemisphere in Supplementary Figure S3a. It should be noted that although the gyral folding patterns and fiber orientation patterns are increasingly complex across the species, the same observation that fiber orientations are positively correlated with tangent gyral folding patterns is well preserved across the 3 primates suggests that axonal fiber pushing might be an evolutionarily preserved mechanism of brain architecture formation and organization. In addition, our approaches have been tested on the same group of chimpanzee subjects scanned with both scan sequences (Supplementary Fig. S4). The results are similar in different scans, thus confirming that our approaches are not sensitive to scan parameters. Finally, the influence of DTI parameters is tested by using different diffusion-weighting gradient directions for tractography for the same subject (Supplementary Fig. S5). Specifically, for the DTI data with 60 directions scanned, we separated the image into 2 groups by the order of direction. One group contains 30 directions with odd index from the original file, whereas the other group contains the rest of the 30 directions. As shown in Supplementary Figure S5, the results are very similar in different groups and are very similar to the results in Supplementary Figure S4 as well, thus further confirming that our approaches are not sensitive to scan protocols and parameters.
Quantification of the Relationship Between Gyral Folding and Fiber Orientation Patterns
In addition to the visual examinations aided by Figure 6 and Supplementary Figure S3, we performed quantitative measurements of the correlation between gyral folding and fiber connection patterns. Specifically, Pearson's correlations between the gyral folding patterns described using the methods in “Quantification of Gyral Folding Patterns” and the fiber orientation patterns described in “Quantification of Fiber Orientation Patterns” were measured via the approaches provided in “The Correlation Between Fiber Orientation Patterns and Gyral Folding Patterns,” and the results are shown in Figure 7. Figure 7a–d shows the results for the right precentral gyrus, the right postcentral gyrus, the left precentral gyrus, and the left postcentral gyrus, respectively. For each part figure, the results for chimpanzee and human brains were provided in the first and second rows, respectively. The average correlations for the right precentral gyrus, the right postcentral gyrus, the left precentral gyrus, and the left postcentral gyrus in chimpanzee brains are 0.69, 0.68, 0.75, and 0.73, respectively. The average correlations for the right precentral gyrus, the right postcentral gyrus, the left precentral gyrus, and the left postcentral gyrus in human brains are 0.72, 0.66, 0.69, and 0.64, respectively. Statistical measurements of mean values and standard deviations are displayed on the top of each set of results in Figure 7. In total, the average Pearson correlation for 4 gyri is 0.67 for 19 human brains and 0.71 for 24 chimpanzee brains we tested. These relatively high positive correlations between gyral folding patterns and fiber orientation patterns in both chimpanzee and human brains further suggest that axonal fiber pushing might be an evolutionarily preserved mechanism of brain architecture formation and organization.
It should be noted that due to the relatively straight shapes of gyri in macaques, for example, as shown in the third columns in Figure 6 and Supplementary Figure S3, their gyral folding pattern quantifications are more prone to noises and not reliable. Hence, we did not include Pearson's correlations of macaques’ gyral folding and fiber orientation patterns in our result reports. The positive relationship between gyral folding patterns and fiber orientation patterns in macaque brains can still be clearly appreciated in the first rows of Figure 6 and Supplementary Figure S3, although the closeness of this positive relationship is less pronounced than that in the chimpanzee and human brains.
Discussion and Conclusion
In this work, we developed, evaluated, and applied a novel computational pipeline of algorithms for the quantification of tangent gyral folding patterns and fiber orientation patterns based on a joint representation approach. These quantifications are then mapped onto the same gyral crest curves to enable and facilitate the quantitative measurement of correlations between fiber orientation patterns and tangent gyral folding patterns. Application of these approaches to experimental data obtained from human, chimpanzee, and macaque subjects revealed a well-preserved strong, positive correlation of the complexities of tangent gyral folding and fiber connection patterns across the 3 primate species. Our qualitative and quantitative analyses also found that tangent gyral fold shapes are connected by dense DTI-derived fibers orienting perpendicularly to the gyral crest curves, which has been replicated and confirmed in 4 major gyri on both hemispheres across the 3 species. Both these findings motivated us to hypothesize an underlying mechanism of axonal fiber pushing that drives the coevolution patterns of gyral folding and fiber connection in primate brains and sculpts the tangent folding patterns of cortical gyri during primate brain evolution. This axonal fiber pushing conjecture well explains not only the folding patterns of gyri along the radial direction of the spherically shaped cortex, but also the folding of gyral crest along the tangent direction on the cortical surface.
Nevertheless, it should be pointed out that our current study only examined the association between axonal orientation patterns and gyral folding patterns, and the potential causal relationship has not been investigated yet. Future studies of the potential causal relationships entail a combination of different experimental and computational approaches. We conjecture that at least the following 2 lines of research efforts are worth pursuit in the future. First, DTI and MRI data sets of human genetic disorders with abnormal axonal wiring patterns can be analyzed and modeled to examine whether or not disrupted axonal connections would be associated with abnormal cortical folding patterns. If it turns out that abnormal axonal wiring caused disrupted cortical folding pattern, it can be considered as a supporting evidence of the causal relationship. For example, the agenesis of the corpus callosum (AgCC) (Tyszka et al. 2011) is a congenital disorder in which corpus callosum is absent. Then, we can compare the differences between the cortical folding pattern of AgCC patients and that of normal controls to see whether or not abnormal axonal wiring would result in disrupted gyral folding pattern. Secondly, biomechanical experiments and simulations can be performed to further elucidate the possible pushing forces within the wiring axons and their mechanical effects on the growing cortical layers. For instance, by using the atomic force microscopy tip as a scalpel to cut axons in cultured neurons, we can release the stress and thus will be able to directly discern between positive and negative tensions by observing the relaxation of the axon. The experimental settings could be similar to those in Xiong et al. (2009). In addition, extensive mechanical modeling and simulations, for example, those in Supplementary Figure S6, can be performed in the future to gain further mechanical evidence and feasibility of the axonal pushing mechanism of cortical folding.
To date, our axonal fiber pushing theory of cortical folding has been based on the qualitative and quantitative analyses of in vivo DTI data sets of macaque, chimpanzee, and human brains. In the future, our hypothesis can be further tested using independent microscale microscopic bioimaging studies. For instance, the axonal fiber wiring patterns in gyral and sulcal regions could be verified by microscopic histology studies using animal models in the future, which can provide independent and reliable confirmation of our macroscale neuroimaging studies. It can even be considered to apply the similar electron microscopic techniques in Bock et al. (2011) to examine the axonal wiring patterns at the single axon resolution and to test whether the axonal terminations are concentrated on cortical gyri. In addition, time-lapsed confocal microscopic imaging studies could be performed on fast-growing animal models or cultured neurons to probe whether or not there is real axonal pushing forces induced by the wiring axonal cones.
As for data analysis, in the present study, we analyzed only 4 major gyri including the precentral gyrus, post central gyrus, superior frontal gyrus, and superior temporal gyrus on both hemispheres to demonstrate the major findings and hypotheses in this paper. In the future, other major gyri such as middle temporal gyrus, middle frontal gyrus, and cingulated gyrus will be investigated for additional confirmation of our findings and theory. Also, other cortical surface shape descriptors and fiber bundle pattern descriptors can be incorporated to measure the complexities of cortical folding and fiber connection patterns. For instance, the curvedness of the gyral crest curve can be quantitatively modeled by a combination of accumulative integral of angles of tangents, bend, centroid ratio, curve bend angle, and concave/convex via the methods in Zhang, Guo, Li et al. (2010). Meanwhile, the fiber connection pattern can be quantitatively described by the novel trace-map models (Zhu, Li, Faraco, Deng, Zhang et al. 2011; Zhu, Li, Faraco, Deng, Zhu et al. 2011). Once both gyral folding patterns and fiber connection patterns are quantitatively modeled, computational simulation models (Nie et al. 2011) based on these quantitative measurements can be adopted to verify the axonal pushing mechanism of cortical folding.
To summarize, our DTI and MRI studies have revealed a positive correlation between axonal fiber orientation and tangent gyral folding pattern across 3 primate species (macaque, chimpanzee, and human), which suggest a hypothesis of axon fiber pushing mechanism of cortical folding. To further test this hypothesis, many future studies should be performed to replicate the existing findings in other primate species via different imaging modalities, to reproduce the existing findings in multiple scales, and to investigate the causal relationship between axonal pushing and cortical folding. Once these biological and biomechanical supports are obtained, we can then possibly establish the axonal pushing mechanism as an underlying theory in the future to uniformly explain the coevolution of gyral folding and axonal fiber connection patterns across primate or mammalian brains. Finally, the potential interactions between the proposed axonal pushing mechanism and other cortical folding mechanisms such as the brain skull constraint (Nie et al. 2010) and behaviors of neurons (Rakic 1988, 2009) and glial cells (Hevner and Haydar 2012) should be studied in the future.
T.L. was supported by the NIH Career Award EB 006878, NIH R01 HL087923-03S2, and the University of Georgia start-up research funding. K.L. was supported by the NWPU Foundation for Fundamental Research. L.L. and X.H. were supported by NIH PO1 AG026423. K.L. and T.Z. were supported by the China Government Scholarship.
The mechanical simulations in Supplementary Figure S6 were assisted by Xin Chen, NWPU. The authors would like to thank the anonymous reviewers for their constructive comments that have helped to significantly improve this paper. Conflict of Interest: None declared.