Parcellation of the primary cerebral cortices based on local connectivity profiles

Connectivity-based parcellation using diffusion MRI has been extensively used to parcellate subcortical areas and the association cortex. Connectivity profiles are vital for connectivity-based parcellation. Two categories of connectivity profiles are generally utilized, including global connectivity profiles, in which the connectivity information is from the seed to the whole brain, and long connectivity profiles, in which the connectivity information is from the seed to other brain regions after excluding the seed. However, whether global or long connectivity profiles should be applied in parcellating the primary cortex utilizing connectivity-based parcellation is unclear. Many sources of evidence have indicated that the primary cerebral cortices are composed of structurally and functionally distinct subregions. Because the primary cerebral cortices are rich in local anatomic hierarchical connections and possess high degree of local functional connectivity profiles, we proposed that local connectivity profiles, that is the connectivity information within a seed region of interest, might be used for parcellating the primary cerebral cortices. In this study, the global, long, and local connectivity profiles were separately used to parcellate the bilateral M1, A1, S1, and V1. We found that results using the three profiles were all quite consistent with reported cytoarchitectonic evidence. More importantly, the results using local connectivity profiles showed less inter-subject variability than the results using the other two, a finding which suggests that local connectivity profiles are superior to global and long connectivity profiles for parcellating the primary cerebral cortices. This also implies that, depending on the characteristics of specific areas of the cerebral cortex, different connectivity profiles may need to be adopted to parcellate different areas.


Introduction
Exploring the structural and functional organization of the brain is one of the most prominent research endeavors in neuroscience. Mapping its subdivisions is a very important aspect of exploring the human brain (Cloutman and Lambon Ralph, 2012). The primary cerebral cortices, including the primary motor (M1), primary somatosensory (S1), primary auditory (A1), and primary visual (V1) cortices, are specialized processing centers for dealing with related primary information (O'Leary et al., 2007). Exploring the fine structure of the primary cerebral cortices has long been an important goal for researchers in the neural sciences. A large amount of work has been done using post-mortem histological tools (Geyer et al., 2000;Schönwiesner et al., 2002) or functional methods (Binkofski et al., 2002;Blankenburg et al., 2003;Formisano et al., 2003). For example, the M1 was identified as having two subregions using the quantitative cytoarchitectonic method (Geyer et al., 1996). The S1 was divided cytoarchitectonically into four areas, 3a, 3b, 1, and 2 (Schönwiesner et al., 2002), and the A1 was segregated into two areas (Galaburda and Sanides, 1980) or three areas (Morosan et al., 2001;Upadhyay et al., 2007). These studies provide evidence that the primary cerebral cortices are indeed composed of structurally and functionally distinct subregions.
Structural connectivity-based parcellation using diffusion MRI can non-invasively subdivide regions in vivo and has been extensively applied to parcellate the subcortical areas and association cortices (see Cloutman and Lambon Ralph, 2012, for a review). The results obtained from structural connectivitybased parcellation tend to resemble the results from functional and cytoarchitectonic studies (Knosche and Tittgemeyer, 2011;Cloutman and Lambon Ralph, 2012). This provides important evidence that structural connectivity-based parcellation is indeed capable of representing the fine substructure of cortical regions. In structural connectivity-based parcellation studies, the formation of connectivity profiles is crucial for the eventual parcellation. Two categories of connectivity profiles, that is, global connectivity profiles and long connectivity profiles, are generally utilized. Global connectivity profiles encompass connectivity information from the seed to the whole brain (Klein et al., 2007;Mars et al., 2011) and long connectivity profiles encompass connectivity information from the seed to the other brain regions after excluding the seed region (Anwander et al., 2007;Beckmann et al., 2009;Wang et al., 2012;Fan et al., 2014). Since the subcortical areas and association cortex receive widespread projections from the distributed brain system, global and long connectivity profiles can adequately represent their characteristics (Mesulam, 1990(Mesulam, , 2008. However, the primary cerebral cortices are heavily myelinated (Glasser and Van Essen, 2011;Glasser et al., 2014), rich in local hierarchical connections (Felleman and Van Essen, 1991), and possess a high density of local functional connectivity profiles (Sepulcre et al., 2010). Therefore, using connectivity-based parcellation to form global or long connectivity profiles and then applying these to parcellate the primary cerebral cortices may not be ideal given the nature of the primary cortices. In contrast, parcellating the primary cerebral cortices using local connectivity profiles, that is, the connectivity information within the primary cerebral cortices, could be quite interesting.
To investigate this, we chose the M1, S1, A1, and V1 as regions of interest (ROIs) and compared them using the three kinds of connectivity profiles (global, long, and local). We first constructed a connectivity matrix for each of the three profiles and then clustered the profiles which subdivided the seed ROIs on two independent datasets. We also parcellated the four primary cerebral cortices based on resting-state functional magnetic resonance imaging (fMRI). The results for each seed ROI were compared with the cytoarchitectonic results from FZ Jülich's SPM Anatomy toolbox (Eickhoff et al., 2005) and with the results from the resting-state fMRI. In addition we compared the results between the three kinds of connectivity profiles. Using this technique we were able to evaluate the usefulness of the different types of profiles for parcellating the primary cerebral cortices and, in particular, to test our hypothesis that local connectivity might be particularly useful in these regions.

Subjects and MRI Data Acquisition
Two independent datasets were used in this study. Dataset 1 included 20 healthy, right-handed subjects (10 males and 10 females, mean age: 18.5 years, range: 17-20 years). All participants provided written informed consent. None of the participants had ever suffered from any psychiatric or neurological disease, and none had any contraindications for MRI scanning. Diffusion-weighted images, T1weighted images, and resting-state functional magnetic resonance imaging (fMRI) were acquired on a 3.0 T GE MR scanner (General Electric, Milwaukee, WI, USA). The diffusion-weighted images were acquired using spin-echo echo-planar imaging (TR = 8.5 s, T = minimum, 75 axial slices, resolution = 2 mm × 2 mm × 2 mm, FOV = 256 mm × 256 mm) in non-collinear 64 directions (b = 1000 s/mm 2 ), along with one non-diffusion weighted volume (b = 0 s/mm 2 ). The T1-weighted images were acquired using a Sag 3D BRAVO sequence (TR = 1.9 s, TE = 3 ms; FOV = 256 mm × 256 mm; in-plane resolution = 1 mm × 1 mm; slice thickness = 1 mm; 192 slices). The resting-state fMRI images were obtained using a gradient-echo single-shot echoplanar imaging sequence (GE-EPI) with the following imaging parameters: TR/TE = 2000/30 ms; FOV = 240 mm × 240 mm; matrix = 64 × 64; FA = 90 • , slice thickness = 3.4 mm; 0.6 mm gap; 39 transversal slices; 255 volumes. All aspects of the study were approved by the Ethics Committee of Tianjin Medical University.

Preprocessing of DTI and fMRI Images
The diffusion and structural MR images were preprocessed using FMRIB's Diffusion Toolbox (FSL 4.1 2 ) and included the following steps: (1) correcting for eddy currents and head motion; (2) co-registering the skull-stripped T1-weighted image to the b = 0 images in native DTI space and then transforming to Montreal Neurological Institute (MNI) space; and (3) transforming the seed masks from MNI space to the native DTI space with nearest-neighbor interpolation using the inverted transformation parameters obtained above.
The resting-state fMRI images were preprocessed using the following steps: (1) discarding the first ten volumes; (2) correcting for slice timing and head motion; (3) intensity scaling of the fMRI images; (4) spatially smoothing and temporally band-pass filtering; and (5) removing nuisance signals, including the signal averaged from the whole-brain mask, the signal averaged from the white matter mask, and the signal averaged from the ventricular mask, six motion parameters, and their first derivatives.

Definition of ROIs
There are many different definitions for the location of the primary cerebral cortices, and the number of subregions in the primary cerebral cortices has varied in previous studies. Because some researchers have suggested that architectonic and probability maps must currently be considered as state of the art for cortical parcellation (Wasserthal et al., 2014), in this study we used the ROIs obtained using FZ Jülich's SPM toolbox 3 (Eickhoff et al., 2005) as the seed masks and calculated the overlap between our parcellation results and the parcellation results this toolbox provided.
The detailed definitions of the M1, S1, A1, and V1 that we used were as follows. The M1 corresponded to BA 4 (Geyer et al., 1996) and two subregions were identified. We chose BA 3 as a seed region for S1, based on a recent study (Keysers et al., 2010), which described this region as being homogenous with the other sensory cortices. Also using other's results (Schönwiesner et al., 2002) as guidance, we subdivided BA 3 into two subregions. A1 was chosen as Te1 (Morosan et al., 2001) and three subregions were identified. We defined BA 17 (Amunts et al., 2000) as V1. Based on the research by Boucard et al. (2009) and Yu et al. (2014), we identified the number of subregions of V1 as 2. Bilateral seed ROIs for these primary cerebral cortices were obtained using FZ Jülich's toolbox. These seed ROIs were transformed into MNI space by an affine registration using FSL (Jenkinson and Smith, 2001).

Parcellation of the Primary Cerebral Cortices based on Diffusion MRI
Voxelwise estimates of the fiber orientation distribution were calculated using Bedpostx within FSL. The probability distributions for two fiber directions at each voxel were calculated 2 http://www.fmrib.ox.ac.uk/fsl 3 http://www.fz-juelich.de/SharedDocs/Downloads/INM/INM-1/DE/Toolbox/ Toolbox 2 0.html;jsessionid=0B82DC28FF06C133EF5ACAD6D841DE54? nn=563092 using a multiple fiber extension  of a previously published diffusion modeling approach (Triarhou, 2007). Thus, the fiber tracts between each voxel in the seed region and every voxel of the whole brain were estimated. We sampled 5000 streamline fibers per voxel to estimate the connectivity probability. The connectivity profile of each voxel in the seed ROI was then used to construct the connectivity matrix (Johansen-Berg et al., 2004). If the number of voxels in the ROI mask is defined as M and the number of voxels in the brain mask is N, the method for calculating the global connectivity profiles included all the connectivity information from the ROI mask to the whole brain, resulting in an M × N connectivity matrix. Using the same designations, the method for calculating the long connectivity profiles included all the rest of the information from the whole brain after excluding the ROI information, resulting in an M × (N-M) connectivity matrix. Similarly, the method using local connectivity profiles included only the connectivity information within the ROI mask, resulting in an M × M connectivity matrix. In this way, three types of cross-correlation (CC) matrices were formed for each subject. Subsequently, spectral clustering was introduced to automatically cluster the CC matrices in order to obtain distinct subregions (Ng et al., 2002). Finally, all the subregions were transformed back into MNI space and a maximum probability map (MPM; Eickhoff et al., 2005) was obtained. This entire process was applied to all the selected ROIs.

Parcellation of the Primary Cerebral Cortices using Resting-State fMRI
The M1, A1, S1, and V1 were parcellated using resting-state fMRI to test the viability of local connectivity profiles. In order to be able to compare these results with those based on connectivity-based parcellation using the local connectivity profiles, we again parcellated the M1, S1, and V1 into two subregions and the A1 into three subregions. The parcellation based on the resting-state fMRI data included the following steps: (1) calculating the Pearson correlation between voxels in the seed region; (2) constructing the correlation matrix; (3) clustering the correlation matrix; and (4) calculating the MPM of the parcellation results for each parcellated region.

Calculation of Consistency and Stability
The consistency and stability of the results obtained from various subjects are often used to justify the acceptance of a parcellation (Kahnt et al., 2012;Li et al., 2013). In this study, we determined consistency with the cytoarchitectonic results by measuring the overlap ratio between our MPM results using connectivity-based parcellation and the results from FZ Jülich's SPM Anatomy Toolbox. We also calculated the consistency between the results based on the local connectivity profiles and the results based on the resting-state fMRI data by measuring their overlap. We determined the stability by transforming the result for each individual subject into MNI space and then calculating the stability as the overlap ratio for each subject's result with the MPM.

Connectivity-Based Parcellation of the Primary Cerebral Cortices
Two subregions, the M1 anterior and M1 posterior, were identified within the M1 for all three different connectivity profiles using connectivity-based parcellation (Figure 1, M1). The parcellation results for both the independent datasets were consistent with the cytoarchitectonic subdivision of M1 (Geyer et al., 1996) from FZ Jülich's SPM Anatomy Toolbox (Eickhoff et al., 2005; Figure 1, M1). Quantitatively, the overlap ratio between our identified regions and the SPM anatomy toolbox were 79.61, 79.25, and 79.41% in the left hemisphere and 77.10, 78.74, and 78.90% in the right hemisphere using the local, global, and long connectivity profiles in dataset 1, respectively (Figure 2). The A1 was divided into three subregions along the medio-lateral axis of Heschl's gyrus (Figure 1, A1). The correspondence between our parcellation results and the cytoarchitectonic results (Morosan et al., 2001) from FZ Jülich's SPM Anatomy Toolbox (Eickhoff et al., 2005) for A1 was good, as it had been for the M1. The consistency between the two results was 82.33, 82.13, and 81.96% in the left hemisphere and 82.63, 82.98, and 80.68% in the right hemisphere using the local, global and long connectivity profiles, respectively. The S1 was parcellated into two subregions, which were arranged medially to laterally (Figure 1, S1). Our results showed 54.11, 53.95, and 58.79% overlap in the left hemisphere and 53.36, 51.39, and 51.86% overlap in the right hemisphere using the local, global, and long connectivity profiles, respectively, with the cytoarchitectonic results (Schönwiesner et al., 2002; Figure 2). We parcellated the V1 into two areas, the anterior and posterior areas. The parcellation of the V1 was performed to ensure the completeness of our work on the parcellation of the primary cerebral cortices and was also used later in the validation of the stability of the local connectivity profiles. Our parcellation results for the M1, A1, S1, and V1 in dataset 2 showed a similar pattern to those from dataset 1 (Supplementary Figure S1).

Resting-State fMRI based Parcellation of the Primary Cerebral Cortices
The results based on resting-state fMRI images (Figure 3) were highly consistent with the results obtained by connectivity-based parcellation using local connectivity profiles. The consistency was 96.74% for the left M1 and 98.96% for the right M1; 74.00% for the left S1 and 77.67% for the right S1; 94.82% for the left A1 and 90.00% for the right A1; and 95.45% for the left V1 and 92.94% for the right V1.

Stability of the Three Parcellation Results
We further compared the parcellation results using the local, global, and long connectivity profiles separately. The parcellation results of the group MPM and four randomly selected subjects from dataset 1 are shown (Figure 4 for the left M1, Figure 5 for the left S1, Figure 6 for the left A1 and Figure 7 for the left V1, Supplementary Figure S3 for the right M1, Supplementary Figure  S4 for the right S1, Supplementary Figure S5 for the right A1, and Supplementary Figure S6 for the right V1). For the M1, V1, and S1, the borders of the subregions that were identified using local  connectivity profiles exhibited significantly higher stability than the results obtained by using global and long connectivity profiles in both hemispheres for dataset 1 (Wilcoxon signed rank test, p < 0.001; Figure 8). For the A1, the parcellation results obtained using local connectivity profiles were more stable than  FIGURE 5 | The primary somatosensory cortex (S1) in the left hemisphere was parcellated into two subregions using local, global, and long connectivity profiles. The MPM from dataset 1 of the S1 (first column) and four individual example results using local, global, and long connectivity profiles are shown.
the results obtained using long connectivity profiles, but no difference was found between the local and global connectivity profiles (p = 0.1942 for the left hemisphere, p = 0.7738 for the right hemisphere). In dataset 2, the borders of the subregions of all four primary cerebral cortices obtained using local connectivity profiles exhibited significant higher stability than the results obtained FIGURE 6 | The primary auditory cortex (A1) in the left hemisphere was parcellated into three subregions using local, global, and long connectivity profiles. The MPM from dataset 1 of the A1 (first column) and four individual example results using local, global, and long connectivity profiles are shown.
FIGURE 7 | The primary visual cortex (V1) in the left hemisphere was parcellated into two subregions using local, global, and long connectivity profiles. The MPM from dataset 1 of the V1 (first column) and four individual example results using local, global and long connectivity profiles are shown.
FIGURE 8 | Comparison of the stability between the results using local, global, and long connectivity profiles from dataset 1; groups with * means that the differences were significant, p < 0.001; error bars indicate mean ± SD.
using global and long connectivity profiles in both hemispheres (Wilcoxon signed rank test, p < 0.001; Supplementary Figure S2), with an exception that there was no significant difference between the results using the local and global connectivity profiles for A1 in the right hemisphere (p = 0.0018).

Discussion
In this study, we proposed that local connectivity profiles might be advantageous for parcellating the primary cerebral cortices (M1, S1, A1, and V1). Global, long, and local connectivity profiles were used to test this hypothesis. We found that the results of the three profiles were all quite consistent with the cytoarchitectonic results. More importantly, the results obtained using the local connectivity profiles showed less inter-subject variability than the results obtained using the other two. This finding suggests that local connectivity profiles are suitable for parcellating the primary cerebral cortices and superior to global and long connectivity profiles. This may also imply that different connectivity profiles should be adopted depending on the structural characteristics of different areas of the cerebral cortices.
Subdivision of the primary cerebral cortices is essential for improving our understanding of the human brain. Since the important milestone work, Brodmann's map (Brodmann, 1909), was constructed, intensive research based on cytoarchitecture, myeloarchitecture, and functional brain imaging has been performed for the primary cerebral cortices. The M1 is primarily responsible for descending motor commands for voluntary movement. Many studies in both monkeys and humans have suggested that the M1 contains subdivisions. For example, Geyer et al. (1996) subdivided the M1 into the M1 anterior and the M1 posterior based on quantitative cytoarchitecture, quantitative distributions of transmitter-binding sites, and positron emission tomography. Rathelot and Strick (2009) found two subdivisions of the M1 in rhesus monkey, based on different distributions of cortico-motoneuronal cells. Nebel et al. (2014) used patterns of correlations in fMRI data to localize five functional subdivisions of the M1. However, previous research into the M1 was primarily based on post-mortem analysis or functional brain imaging. In this study, we subdivided the M1 into the M1 anterior and the M1 posterior (Figure 1, M1) using a non-invasive, connectivity-based parcellation method with local, global, and long connectivity profiles. The borders of the two subregions obtained using local connectivity profiles were almost the same as the results obtained using global and long connectivity profiles and corresponded well with previous results (Geyer et al., 1996). The consistency between our results and cytoarchitectonic results indicates that connectivity-based parcellation is appropriate for the M1.
The A1 is primarily responsible for processing auditory information in humans, and research has indicated that the A1 is tonotopically organized (Talavage et al., 2004), but the precise locations of these tonotopic areas remains unclear. Much research (Kaas and Hackett, 2000;Morosan et al., 2001) has been devoted to investigating the subdivisions of the A1 because of the good correspondence between these subdivisions and the tonotopic organization within it (Kaas and Hackett, 2000;Formisano et al., 2003). Currently, the human A1 is difficult to identify in vivo and the information about it has primarily been defined from post-mortem brains. For example, the A1 has been variously described as BA 41 (Brodmann, 1909), area TC (von Economo andHorn, 1930), or Te1 (Morosan et al., 2001). Additionally, the number of subdivisions of the A1 has varied between different studies (Iannetti et al., 2003;Sewards, 2011;Wasserthal et al., 2014). In this study, we chose the cytoarchitectonic map of Te1 in Morosan et al. (2001) as a seed ROI for the A1 and identified it as having three subregions (Figure 1, A1). The results obtained using local connectivity profiles were in good agreement with the results obtained using global and long connectivity profiles and were similar to the Te1.1, Te1.0, and Te1.2 in Morosan et al. (2001). Our parcellation results were also structurally in agreement with the functionally mirror image tonotopic organization of the A1, which consists of a high-low-high frequency gradient from the rostrolateral to the caudomedial direction (Formisano et al., 2003;Upadhyay et al., 2007;Langers and van Dijk, 2012). Thus, our parcellation of the A1 into three subregions based on the structural connectivity-based method is reasonable.
The S1 is primarily involved in processing tactile and nociceptive stimuli and is complex and hierarchically organized (Ploner et al., 2000). Identifying the subregions of the S1 should be helpful for studying brain representations of the body. Traditionally, the S1 has been described as an area consisting of 3a, 3b, 1, and 2 (Brodmann, 1909). In this study, we followed the recent definition by Keysers et al. (2010) that defined the S1 as BA 3, because of its homogeneity with the other sensory cortices. We found that the two subregions of the S1 (Figure 1, S1) identified using connectivity-based parcellation were different from Geyer's results (Schönwiesner et al., 2002), a difference that may be attributable to our using different imaging modalities. Additionally, the BA 3a region is so small (Schönwiesner et al., 2002) that current low-resolution diffusion MRI brain images cannot parcellate it well. In any case, our results obtained using local connectivity profiles were still in accord with the results obtained using global and long connectivity profiles, and the medial-lateral organization was in good agreement with the somatotopic maps of the fingers that have been identified for 3b (Sanchez-Panchuelo et al., 2012. This finding suggests that connectivity-based parcellation is applicable to the S1. Area V1 is primarily used to respond to the visual world and is thought to be relatively functionally monomodal. However, increasing evidence indicates that the V1 is a more complex and comprehensive area than previously thought (Petro et al., 2014). For example, Sugita (1999) and Lee and Nguyen (2001) found that the V1 could respond to stimuli that were not directly presented to the retina. Zhu et al. (2014) found that the V1 was hyperactive in individuals with post-traumatic stress disorder even during eyes-closed, resting-state fMRI scanning. Further, some recent studies suggested that the V1 may not be a single anatomic region. Yu et al. (2014) found that the anterior and posterior subregions of the V1 in patients with primary openangle glaucoma differed in their changes in cortical thickness compared to the corresponding subregions of normal subjects. The FZ Jülich's SPM Anatomy toolbox does not provide information about the subregions of the V1. In the present study, we used local, global, and long connectivity profiles to identify subdivisions of the V1. We found that the results obtained using the three types of profiles were highly consistent and were in line with the anterior-posterior pattern of the V1 provided by Boucard et al. (2009) and Yu et al. (2014). The results obtained in this current study may provide a new view of the V1. Considering all the above evidence, we conclude that the parcellation results of the primary cerebral cortices obtained using local connectivity profiles correspond very well with the results obtained using global and long connectivity profiles as well as with the probabilistic cytoarchitectonically defined cortical areas.
In this study, the M1, A1, S1, and V1 were also parcellated using resting-state fMRI to test the viability of local connectivity profiles. The results (Figure 3) obtained from resting-state fMRIbased parcellation and the results obtained using local connectivity profiles-based parcellation corresponded very well. These findings further validate the viability of the local connectivity profiles-based parcellation.
The present study found that the borders of the subregions of the primary cerebral cortices obtained using local connectivity profiles showed less inter-subject variability compared to the borders of those obtained using global and long connectivity profiles (Figure 8). These findings have many potential explanations. First, convergent evidence indicates that the primary sensory and motor areas are heavily myelinated (Sigalovsky et al., 2006;Sereno et al., 2013;Glasser et al., 2014;Wasserthal et al., 2014). Some studies have suggested that horizontal fibers within the cortex consist primarily of local axonal ramifications of pyramidal neurons (Levitt et al., 1993). These horizontal fibers are important intracortical fibers between myelinated cells (Hellwig, 2002). Thus, it is reasonable for us to assume that the heavily myelinated primary cerebral cortices are rich in local axonal ramifications. From the above, we can infer that the primary cerebral cortices may be rich in local structure and function, which could lead to our finding local connectivity profiles, which are statistical descriptions of local axonal ramifications. On the other hand, the information processing steps of the primary cerebral cortices are fixed, which is beneficial for dealing with continuingly recurring standard tasks (Hellwig, 2002). This characteristic was partially confirmed by results using BOLD fMRI that showed that the primary cerebral cortices exhibited high local functional connections based on an intrinsic activity correlation (Sepulcre et al., 2010). Therefore, unlike the subcortical areas and the association cortex, the primary cerebral cortices may primarily exchange information that is constrained to local regions. Finally, the strength of the connectivity profiles obtained using tractography, as well as the certainty about the orientation measurements, decreases with distance between the source and target areas (Raj and Chen, 2011). Thus, local connectivity profiles, especially for the primary cerebral cortices, which are heavily myelinated and rich in local functional connectivity, may be less influenced by these deficiencies than global and long connectivity profiles. In summary, we therefore suggest that local connectivity profiles are more suitable for parcellating the primary cerebral cortices.

Conclusion
We found that local connectivity profiles seemed to do a better job of parcellating the primary cortices. Prior to our study, many subcortical areas and the association cortex have been subdivided by connectivity-based parcellation using global or long connectivity profiles. Our findings and previous studies imply that, in order to obtain accurate and useful parcellations, frameworks should be adopted that correspond to the characteristics of the various regions of the cerebral cortices. Our work provides a new perspective for completing the parcellation framework.