Group‐guided individual functional parcellation of the hippocampus and application to normal aging

Abstract Aging is closely associated with cognitive decline affecting attention, memory and executive functions. The hippocampus is the core brain area for human memory, learning, and cognition processing. To delineate the individual functional patterns of hippocampus is pivotal to reveal the neural basis of aging. In this study, we developed a group‐guided individual parcellation approach based on semisupervised affinity propagation clustering using the resting‐state functional magnetic resonance imaging to identify individual functional subregions of hippocampus and to identify the functional patterns of each subregion during aging. A three‐way group parcellation was yielded and was taken as prior information to guide individual parcellation of hippocampus into head, body, and tail in each subject. The superiority of individual parcellation of hippocampus is validated by higher intraregional functional similarities by compared to group‐level parcellation results. The individual variations of hippocampus were associated with coactivation patterns of three typical functions of hippocampus. Moreover, the individual functional connectivities of hippocampus subregions with predefined target regions could better predict age than group‐level functional connectivities. Our study provides a novel framework for individual brain functional parcellations, which may facilitate the future individual researches for brain cognitions and brain disorders and directing accurate neuromodulation.

AD in elderly individuals (O'Brien et al., 2010;Salami, Pudas, & Nyberg, 2014). At the subregion-level, recent studies found lower functional connectivity between posterior hippocampus and medial prefrontal cortex in venerable age compared to younger age population (Damoiseaux, Viviano, Yuan, & Raz, 2016;Wang et al., 2010). All these studies demonstrated that to delineate the functional patterns of hippocampus is important to reveal the neural basis of cognition and memory decline during aging (Andrewshanna et al., 2007;Koch et al., 2010).
Hippocampus is a functionally heterogamous area with different functional subregions having different functions (Das et al., 2011;Zarei et al., 2013;Zhong et al., 2019). The anterior hippocampus is mainly involved in emotion processing while the posterior hippocampus primarily participates in spatial memory (Blum, Habeck, Steffener, Razlighi, & Stern, 2014). Thus, to uncover the aging effects at the subregion-level may better identify the neural mechanism of aging.
Although subregions of the hippocampus have been defined with cytoarchitecture or connectivity-based parcellation approaches, however, different parcellation studies proposed different parcellation schemes for hippocampus (Amunts et al., 2005;Cheng & Fan, 2014;Cheng, Zhu, Zheng, Liu, & He, 2020;Robinson et al., 2015;Zhong et al., 2019). Moreover, all the previous parcellation schemes for hippocampus are based on group-level mapping and thus are hardly to reflect the individual differences.
Individual differences in cognition and behaviors have been widely reported in neuroscience researches (Brown, 2017;Tavor et al., 2016). The brain structural and functional variability across individuals has also been identified in recent literatures (Mueller et al., 2013). The individual brain functional parcellation and functional networks mapping were developed and employed to reveal the functional underpinnings for brain cognitions and disorders (Han et al., 2019;Wang, Buckner, et al., 2015). However, all the existing individual brain mapping approaches are based on the predefined prior brain atlas to guide parcellation, and thus they are constrained by the parcellation schemes of the used atlases (Han et al., 2019;Li et al., 2019;Wang et al., 2018).
Although consistent anatomical and functional topological architecture of some brain areas have been documented, the topological variations between structural, anatomical, and functional parcellations have been demonstrated (Caspers et al., 2008;Ruschel et al., 2014;Wang et al., 2012;Wang et al., 2017). Thus, atlases-free individual parcellation and functional mapping are essential and may better characterize the individual functional organization to delineate the individual differences compared to group-level mapping.
In this study, we proposed a new individual brain parcellation approach using affinity propagation clustering which is able to automatically define the number of clusters (Frey Brendan & Detbert, 2007). Then, we validated the superiority of our proposed method by measuring the functional similarities compared to groupbased parcellation results. In addition, the relationship between individual variations and coactivation patterns of hippocampus mostly related cognitive functions were evaluated. Finally, we applied this approach to define individual functional subregions of hippocampus in aging population and to test whether individual parcellation could better predict age than group-level parcellation.

| Subjects
A public adult lifespan resting-state functional magnetic resonance imaging (rs-fMRI) data (Southwest University Adult lifespan Dataset, SALD) was accessed through f1000 project (http://fcon_1000. projects.nitrc.org/indi/retro/sald.html). By excluding the bad fMRI data with large head motion (see the following fMRI preprocessing), a total of 262 healthy subjects (101 males/161 females, age range of 19-75 years, mean and SD = 42.03 ± 16.73 years) with high quality rs-fMRI were finally used in this study.

| Resting-state fMRI data acquisition
The rs-fMRI data were scanned using a Siemens 3T Tim Trio MRI scanner with echo planar imaging sequence. The subjects were instructed to lie down, close their eyes, and rest without thinking anything and not fall asleep. The acquisition parameters were the followings: TR = 2,000 ms, TE = 30 ms, voxel size = 3.4 Â 3.4 Â 3 mm 3 with 1 mm gap, 32 axial slices, and 242 volumes. The detailed information about subjects and rs-fMRI scanning parameters can be found in a previous study (Wei et al., 2018).

| Resting-state fMRI data pre-processing
The preprocessing of resting-state fMRI data including the following steps: discarding the first 10 volumes to facilitate magnetization equilibrium; realigning all the remained volumes to the first volume to correct head motion; normalizing to the standard EPI template in MNI space; regressing out Friston 24-parameter model of head motion, white matter, cerebrospinal fluid and global mean signals; filtering with a temporal band-pass of 0.01-0.1 Hz. To exclude the head motion effects, the data were discarded if the head-movement exceeded 2 mm or 2 in any direction. Moreover, scrubbing was further used to eliminate the bad images (before two time points and after one time point) exceeding the pre-set criteria (frame displacement [FD], FD < 0.5) for excessive motion.
If the deleted number of volumes exceeding half of the time points, that is, 116 volumes for each subject, this subject was discarded for the following analysis. Under these criteria, 232 subjects were excluded and the remained 262 subjects were used for analyses.

| Hippocampus definition
To define the individual functional subregions of hippocampus, the human bilateral hippocampus seed masks were defined using Harvard-Oxford cortical atlas with 25% probability (Desikan et al., 2006). After obtaining the hippocampus masks, the bilateral masks were downsampled into 3 mm cubic voxel to calculate functional connectivity for parcellation.

| Group-guided individual parcellation of hippocampus
To achieve individual parcellation of hippocampus, we first performed the group-level parcellation of hippocampus across all the subjects.
For group-level parcellation, we first calculated the whole brain functional connectivity map for each voxel within hippocampus. Next, the similarity for the functional connectivity maps of every pair of voxels within the hippocampus was defined using eta 2 (Wang, Yang, et al., 2015;Wang et al., 2016), and a similarity matrix S for each subject was obtained. Then, all the similarity matrixes were averaged to obtain one similar matrix which was clustered to obtain the parcellation results of hippocampus. To automatically identify the optimal number of clustering, the affinity propagation method (Frey & Dueck, 2007) was employed to segment the similarity matrix.
The AP algorithm requires two input parameters: input similarity matrix S, in which S(i,k) is the negative value of the squared Euclidean distance between points i and point k in the similarity matrix (for simi- and preference "P," recommended as the median of S(i,k) when there is no prior, determines the number of clustering (Frey & Dueck, 2007).
The algorithm takes all data points as potential exemplars and calculates the following two messages including the "availability" a(i,k) and the "responsibility" r(i,k)) to characterize the appropriateness of exemplars selection for the data points iteratively to obtain the optimal clustering results (Zhang, Li, et al., 2011;Zhang, Tuo, et al., 2011). The number of clusters was finally determined when the iterative process converges. The details of AP clustering are as follows (Xia et al., 2008): Step 1: initialize the availabilities: a i, k ð Þ¼0: Step 2: update the responsibilities: Step 3: update the availability: a k, k ð Þ X i 0 s:t::i 0 ≠ k max 0,r i 0 , k À Á È É Step 4: iterations convergence terminate: For parameter of preference "P," the recommended value usually does not obtain a satisfactory result. In our study, we used a semisupervised algorithm for step optimization to search the optimal "P." The optimal "P" value was determined by identifying the largest silhouette value across different "P" values (Kaufman & Rousseeuw, 1990). At a specific clustering solution, silhouette value, sil, is defined as follows: where a represents intracluster similarity, that is, the average distance between data point i and the other points in a specific cluster, and b is the interclusters similarity, that is, average distance between data point i in one cluster and all the data points in the other cluster. The average of the silhouette values is defined: where N is the total number of data points. "AS" reflects the quality of clustering results varying from 0 to 1, and a larger "AS" represents better clustering quality. To initialize the "P" value, S med and S min as median and min of S(i,k) were calculated, respectively. The "P" changed from S med recommended by the original paper to the values of the end of search (S med + step Â number of steps, where step = S min / T, T is a random value to determine the length of step, T = 100 used in this study). If the "AS" values go smooth from one step to all the other steps behind it, the corresponding "P" value is considered to be optimal. The smooth of the changing "AS" values was characterized by the gradient values calculated as the latter value minus the previous value, and no blunt change of the gradient value was considered to be smooth. The number of clusters corresponding to the optimal "P" was taken as the optimal group-level parcellation results.
Next, the group-level parcellation of hippocampus was used to guide the individual parcellation using Litekemans method which has fast speed and high accuracy (http://www.cad.zju.edu.cn/home/ dengcai/Data/code/litekmeans.m). The mean whole brain functional connectivity map of each subregion of hippocampus yielded by group-level parcellation was first calculated and taken as initial clustering centers across all the voxels within each subregion and across all the subjects. Finally, the whole brain functional connectivity maps of all the voxel in hippocampus in each subject were clustered to achieve individual parcellation of hippocampus.

| Validation of the individual parcellation approach
To test and validate the reliability of the developed group-guided individual functional parcellation approach, we applied this method to parcellate supplementary motor area (SMA), a brain region that has been widely used to test the applicability and reliability of resting-state fMRI-based functional parcellation because of its established functional architecture (Johansen-Berg et al., 2004;Kim et al., 2010;Zhang et al., 2015). SMA mask was defined using automated anatomical labeling template (Tzourio-Mazoyer et al., 2002). Then, individual parcellation of SMA was performed using our developed groupguided individual parcellation approach. Both two-and three-way parcellations of SMA were performed in both group and individual levels referring to previous studies.

| Functional similarity validation
To validate the superiority of individual parcellation results compared to group-level parcellation, intraregional functional connectivities and time series similarities were calculated in this study. The intraregional functional connectivity similarity is the average Pearson correlation coefficients between any pair of voxels in each subregion. The intraregional time series similarity is characterized using Kendall's coefficient concordance (KCC) of all the voxels' time series in a specific subregion (Kendall, 1990).
where R i is the sum rank of the ith voxel in a specific subregion of hip-

| Age prediction
To further validate the superiority of individual functional parcellation, the individual and group-level functional connectivities were separately taken as features to predict individual ages. To define hippocampus involved functional networks, the voxel-wise whole brain functional connectivity analysis of bilateral hippocampus was first performed to define the target brain areas. One-sample t test was used to identify the significantly connected brain areas with hippocampus.
The significant level was corrected using family discovery rate method with p < .001, and minimum cluster size > = 100.
Next, the functional connections of each hippocampus subregions yielded by individual and group-level parcellation with these target brain areas were calculated and taken as features for prediction. Relevance vector regression method which showed better performances than other prediction methods was applied to predict the age (Cui & Gong, 2018;Cui et al., 2016). A 10-fold cross-validation strategy was used to estimate the generalization ability of the prediction mode (Tang et al., 2018). Pearson correlation coefficient between the real age and predicted age was calculated to depict the prediction performance. The mean absolute error (MAE) values between real and predicted age were also calculated to evaluate prediction.

| Age-related functional connections with hippocampus
To explore age associated functional connections with hippocampus, correlation analyses were performed between seed-to-target functional connections and age. The seed to targets functional connec- tions measured with Pearson's correlation coefficients between time series were calculated for individual and group-level hippocampus subregions. The Fisher r-to-z transformation was used to change the r value to z value. Finally, correlation analysis between each individual and group-level functional connectivity and age was performed, and the significant level was set at p < .05 with Bonferroni correction.

| Comparisons with iteration adjusted individual parcellation approach
We also compared our method with the previous iteration adjusted individual parcellation approach developed by Wang, Buckner, et al. (2015). 3 | RESULTS

| Individual functional parcellation of the hippocampus
To determine the number of clusters, the maximum average silhouette value and smooth variation of gradient of the average silhouette were taken as the criteria to terminate searching the optimal preference parameter (Figure 1a,b). Based on the criteria, the optimal three-way group-level parcellation of bilateral hippocampus was found. The three-way group-level parcellation identified head, body, and tail of hippocampus ( Figure 2a). Using the group-level parcellation results of hippocampus as prior information, the individual parcellation results for hippocampus were obtained, and five randomly selected individual parcellation results are shown in Figure 2b. The obvious individual variations of hippocampus parcellation results could be observed.
To quantitatively describe the individual variations, we calculated the overlap degree at the voxel level for each subregion across all the subjects ( Figure 2c). The tail part in bilateral hippocampus show high overlap degree, while the overlap degree in body and head part of hippocampus is relatively lower than tail part, especially at the intersected part.
To validate our developed group-guided individual parcellation approach, we adopted the same procedures to parcellate SMA in which the functional topography has been well established. The group-guided individual functional parcellation identified similar group-level functional topography of this areas identified by previous studies (Crippa et al., 2011;Johansen-Berg et al., 2004;Kim et al., 2009;Zhang et al., 2015). The individual variations of functional subregions in SMA could be obviously observed in five randomly selected subjects ( Figure S1).

| Age prediction
Based on the whole brain functional connectivity mapping of hippocampus, eight target brain areas showing significantly functional connections with hippocampus were defined ( Figure S4 and Table S1).
After obtain the eight target brain areas, the functional connectivity between each individual and group-level hippocampus subregion and each target was calculated.
To explore whether individual functional connectivities can better  Finally, individual functional connections of individual hippocampus subregions defined using iteration adjusted individual parcellation approach were taken as features to predict age. Correlation analysis showed that individual functional connections can well predict age but with lower correlation coefficient compared to our method ( Figure S5).

| Associations between hippocampus connectivities and age
F I G U R E 3 Functional similarities comparison between the group and individual parcellation results. The intraregional functional connectivity (a) and time series (b) similarities between any pair of voxels within each hippocampus subregion were calculated. Almost all the functional connectivity and time series similarities of hippocampus subregions at the individual level are higher than that at group level except the head part of right hippocampus in which the intraregional functional connectivity and time series similarities at group level are higher than individuals. (*represents p < .05 after Bonferroni corrected)

| DISCUSSION
This study proposed a novel individual brain parcellation approach based on the semisupervised affinity propagation clustering using resting state fMRI. The proposed method was subsequently applied to parcellate hippocampus to study the aging effects on hippocampus.
The hippocampus was successfully subdivided into three parcels in each individual guided by a three-way group parcellation results.

| Affinity propagation algorithm for brain parcellation
Affinity propagation algorithm could automatically determine the number of clusters by setting the preference "P" value which is recommended as the minimum or median value of the input similarity matrix to obtain small or moderate number of clusters in the original affinity propagation method (Frey Brendan & Detbert, 2007). However, given that the recommended input parameter "P" usually does not work to achieve the ideal clustering results, especially for brain parcellation, thus, in this study, we set a range of "P" values and the optimal "P" value was selected with high silhouette values and small or smooth gradient variations across these silhouette values. Our findings demonstrated that the method to select the optimal number of F I G U R E 4 Correlations between individual variations and activities patterns of cognitive functions within hippocampus. Voxel-wise correlation analyses between individual variations of hippocampus parcellation results defined using our proposed approach and coactivation patterns, that is, activation likelihood estimation (ALE) values, of cognitive memory, episodic recall, and emotion of angry of hippocampus were performed. Significant correlations between individual variations and ALE values of cognitive memory, episodic recall, and emotion of angry were found F I G U R E 5 Age prediction using individual and group-level functional connections. Both individual and grouplevel functional connectivities were separately taken as features to predict age. The Pearson correlation coefficient was used to characterize the prediction performance. The individual functional connections showed better performance, that is, higher correlation coefficient, than group-level functional connections-based prediction clusters could obtain reliable parcellation results by comparing our findings with previous hippocampus parcellation results.

| Individual brain parcellation
The human brain is a complex system with huge individual differences in brain structure and functions (Mueller et al., 2013). To map the individual functional atlas is prevalent to reveal the individual differences in behaviors and cognitions (Goulas et al., 2012;Nebel et al., 2014). With predefined 18 cortical networks as prior information, Wang, Buckner, et al. (2015) developed an individual cortical parcellation approach and demonstrated that the functional networks mapped by this approach were highly reproducible within subjects and effectively captured the variability across subjects.
Subsequently, Wang et al. (2018)  To avoid these problems, we developed a fully datadriven approach without prior information as guide to map the individual functional atlas. Our method was applied to hippocampus and obtained reliable group and individual parcellation results.
Although the hippocampus was parcellated with different techniques, almost all of these parcellation results are at group level not considering the individual variations. In our study, we developed a group-level guided individual parcellation of hippocampus, which could not only capture the population information but also characterize the individual variations.

| Hippocampus and aging
The human brain undergoes a series of complex structural and functional decline during aging (Damoiseaux, 2017). Existing researches demonstrate that the overall gray matter volume of hippocampus in the elderly is significantly smaller than in healthy adults (Apostolova et al., 2012;Fleischman et al., 2014;Golomb et al., 1993;Miller & O'Callaghan, 2005), and the anterior hippocampus shrinks faster than the posterior hippocampus during aging (Chen et al., 2010).
Aging is also closely related to functional integrity of hippocampus, but the detailed relationship is still controversial (Blum et al., 2014). The functional connectivities between hippocampus and the default mode network including posterior cingulate cortex, medial prefrontal cortex, and lateral parietal cortex are negatively correlated with age. The functional connectivity between the left and right posterior hippocampus decreases with age (Damoiseaux et al., 2016).
However, in another study, the coupling between the bilateral hippocampus was found to increase with age (Salami et al., 2014). In addition, the precuneus and bilateral medial temporal lobes were found to have stronger connections to the anterior hippocampus in young people, but have stronger connections to the posterior hippocampus in the elder (Blum et al., 2014). In spite of inconsistency of these findings, all the evidence suggests an important role of hippocampus in aging. In our study, we found the functional connectivities between hippocampus and medial prefrontal cortex, superior frontal gyrus, retrosplenial cortex, and cerebellum were negatively correlated with ages. Moreover, we found that the functional connections of hippocampus with its involved network could well predict age. Our findings provide further supporting information that to characterize the functional patterns of hippocampus is crucial to reveal the neural basis of cognitive decline during aging and disease states.

| Limitations and prospects
The current research has several limitations. First, how to choose the preference "P" value is very hard. Although we set a range of "P" values to find the optimal solution in our study, whether a better method to determine the optimal "P" should be further investigated.
Second, we did not applied this approach to clinical patient data, the stability of our proposed method needs to be further validated since the functional connectivity pattern may changes greatly during disease conditions.

DATA AVAILABILITY STATEMENT
All data and code used for data analysis are available upon request.