About the Authors:
Guorong Wu
Affiliation: Department of Radiology and BRIC, University of North Carolina at Chapel Hill, Chapel Hill, NC, 27599, United States of America
Xuewei Peng
Affiliations Department of Radiology and BRIC, University of North Carolina at Chapel Hill, Chapel Hill, NC, 27599, United States of America, Department of Electrical & Computer Engineering, Texas A&M University, College Station, TX, 77843, United States of America
Shihui Ying
Affiliation: Department of Mathematics, School of Science, Shanghai University, Shanghai, 200444, China
Qian Wang
Affiliation: Med-X Research Institute, Shanghai Jiao Tong University, Shanghai, 200240, China
Pew-Thian Yap
Affiliation: Department of Radiology and BRIC, University of North Carolina at Chapel Hill, Chapel Hill, NC, 27599, United States of America
Dan Shen
Affiliation: Department of Mathematics and Statistics, University of South Florida, Tampa, FL, 33620, United States of America
Dinggang Shen
* E-mail: [email protected]
Affiliations Department of Radiology and BRIC, University of North Carolina at Chapel Hill, Chapel Hill, NC, 27599, United States of America, Department of Brain and Cognitive Engineering, Korea University, Seoul, Republic of Korea
Introduction
Since the advent of Magnetic Resonance Imaging (MRI), many imaging-based studies have been initiated to study structural variations within a population, between populations, and across different times of the same population [1–5]. In these studies, image registration is key for removing structural variations that confound the analysis of group differences or changes, e.g., those associated with brain disorders [6–8]. In particular, subtle changes may be elusive in cases where image registration is not performed with sufficient accuracy.
Many pairwise registration methods have been developed for registering a group of images to a template image [9–14]. However, template selection is not a trivial task [15, 16] and, if done improperly, will bias the subsequent statistical analysis. To deal with this issue, groupwise registration methods have been recently proposed to align all images jointly onto a common space, without the need of explicitly specifying the template. The objective function in groupwise registration aims to minimize either the overall intensity difference [17] or the entropy of joint intensity distributions across all images [18, 19]. To optimize the large-scale objective function in groupwise registration, efficient gradient-based Gauss-Newton optimization is proposed in [20]. A hierarchical groupwise registration mechanism is also used in [21] by selecting the key points in the image and only letting the key points drive the entire groupwise registration. The advantage of groupwise registration over pairwise registration has been widely demonstrated in the literature [15, 17, 19, 20, 22–24].
Many existing groupwise registration methods require an explicitly defined target of registration. For example, Joshi et al. [25] proposed an efficient groupwise registration algorithm by alternating between two steps: (1) register all images, separately, to a tentatively estimated group-mean image, and (2) update the group-mean image by averaging all registered images. This groupwise registration method works well with the datasets involving only small structural variations. However, for the datasets with large and complex variations, this method could result in a fuzzy group-mean image and thus decrease the registration accuracy. To solve this problem, Wu et al. (2011) proposed a sharp-mean based groupwise registration method by using a patch-based weighted-averaging method. However, a common limitation of these methods is that registration needs to be performed with respect to the group-mean image, regardless of whether the appearance of an individual image under registration is significantly different from the group-mean image or not. These methods are hence limited in dealing with the datasets with large and complex structural variations.
On the other hand, the image distribution information has been shown useful in guiding groupwise registration. For example, Wang et al. (2010) has demonstrated that the improved registration accuracy can be obtained by clustering the images in a population into several subgroups and then registering each of the subgroups with groupwise registration. Since images in a subgroup are similar in appearance, accurate registration can be obtained with relative ease. After within-subgroup registration, a representative image from each subgroup can be used for between-subgroup registration. However, one limitation of this method is that the within-subgroup distribution information is not specifically used for groupwise registration.
We developed an approach, called ABSORB (Atlas Building by Self-Organized Registration and Bundling) [17], to leverage image distribution information to guide groupwise registration. In ABSORB, each image is registered to a small number (e.g., m = 3∼5) of neighboring images, which are similar in appearance. Each image is then spatially transformed using the average of the m resulted deformations. This procedure is applied to all images iteratively. In this way, ABSORB avoids the need to register two images with large structural differences. However, the value of m is empirically determined in ABSORB and then fixed throughout the whole registration process. In addition, ABSORB is not aware of the image distribution of the entire dataset, but only a few neighbors.
To combine the advantages of above-mentioned hierarchical registration framework and ABSORB, we further proposed an approach called HUGS (Hierarchical Unbiased Graph Shrinkage) for groupwise registration [23]. In HUGS, we characterize the image distribution of the entire dataset by using a graph with its nodes representing the images and its edges representing the degree of similarity between the images. Only similar images are connected on the graph. The registration of all images to a hidden common space is formulated as a dynamic graph shrinkage problem, where all nodes move closer to each other along the edges. Since both global and local image distribution information is encoded in the graph, the topology of the entire image manifold can be preserved more accurately throughout registration.
Despite being more accurate, HUGS is still limited in dealing with heterogeneous datasets since a single simple graph is often insufficient in modeling complex image distributions. HUGS uses a simple threshold-based method to construct the graph by connecting image pairs with distances smaller than the given threshold. Given a heterogeneous dataset, the threshold has to be significantly relaxed in order to ensure that all images are connected on the graph, which often results in an over-connected and less-efficient graph. This will also cause many unnecessary registrations between dissimilar image pairs, which not only significantly increases the computation time but also results in considerable registration errors.
In this paper, we present a novel groupwise registration method that inherits the advantages of all our previously developed methods but, at the same time, addresses their limitations. Specifically, we introduce an improved version of HUGS, called enhanced HUGS or eHUGS, to deal with the registration of heterogeneous datasets by constructing a hierarchical graph to model the complex image distribution. The main idea is that we first represent the manifold of the population of brain images by a graph that characterizes the affinity or similarity between any two images. If two images are structurally similar, they are connected by an edge on the graph. This graph is used to help avoid the need to register two images with large difference. Since on the graph the images are represented as nodes and the distances between images are represented as edges, the groupwise registration problem then becomes a graph shrinkage problem, where the goal is to shrink the edges of the graph by progressively registering images adjacent on the graph. To construct a hierarchical graph, the images are first clustered into several subgroups via affinity propagation [26], which can group images into subgroups based on a set of automatically determined representative images. A representative image is identified to represent each subgroup, and the distribution in each subgroup is modeled using a graph. Then, the distribution of all these “low-level” graphs is captured by using “high-level” graphs that model the distribution of the representative images from all subgroups. After building the hierarchal graph, dynamic graph shrinkage is finally employed to register all images to the common space.
In the experiments, we first tuned the parameters and validated the registration performance of eHUGS by using 20 longitudinal series of images [27], each with 4 time points. After that, we applied eHUGS to 285 images of infants with ages ranging from 2 weeks to 2 years. We evaluated the quality of infant atlases generated by averaging those registered images in the common space. The results indicate that eHUGS yields greater registration accuracy than other groupwise registration methods, reflected by its built infant atlases with much clearer anatomical details.
Ethics Statement
This study does not involve human subjects or samples. The individuals in this manuscript have given written informed consent at the time of enrollment for imaging and completed questionnaires approved by each participating sites. The authors have obtained approval from the ADNI Data Sharing and Publications Committee to use the data. The authors confirm that the data was analyzed anonymously.
Materials and Methods
Given a large dataset of heterogeneous brain MR images, our goal is to simultaneously register all images to a hidden common space. The group-mean image does not need to be explicitly specified, thus avoiding concerns of how the low quality of group-mean images might undermine registration accuracy [22]. At the same time, the topological distribution of the entire image dataset will also be preserved throughout registration.
To achieve the above goal, we propose a hierarchical graph to model the distribution of all images in the dataset, and then attain groupwise registration by shrinking the graph dynamically until all images are warped to the hidden common space. For clarity and conciseness, we use I = {Ii|i = 1, …, N} to denote a group of N images that are distributed on the manifold. The term di,j denotes the distance between image Ii and image Ij. Since our focus is on deformable image registration, all images in the dataset are first affine-registered using FLIRT [28].
Overview of eHUGS
An overview of eHUGS is provided in Fig 1. eHUGS is an iterative method to deform all images together gradually as the graph shrinks. First, image similarity is computed based on a current set of deformed (or initially the affine transformed) images. Then, a hierarchical graph representing the distribution of the deformed images is built based on manifold M, which is determined using image clustering. In contrast to HUGS, eHUGS partitions the heterogeneous image group into several subgroups. Here, each subgroup contains only the homogeneous data with images of similar appearance. A low-level graph is then constructed to capture the distribution of images in each subgroup. These low-level graphs are further connected using the high-level graphs. Each image is allowed to deform only towards its connected counterparts in this hierarchical graph, progressively bringing the images closer to each other.
[Figure omitted. See PDF.]
Fig 1. An illustration of Enhanced Hierarchical Unbiased Graph Shrinkage (eHUGS) method.
In this example, the images are clustered into 5 subgroups based on their similarity. Then a hierarchical graph is constructed with each low-level graph (indicated by solid lines) describing the image distribution in each subgroup, and each high-level graph (indicated by dashed lines) encodes the interaction between the subgroups. Groupwise registration is formulated as the dynamic shrinkage of this hierarchical graph.
https://doi.org/10.1371/journal.pone.0146870.g001
Groupwise registration is carried out by shrinking the graph gradually to deform all images to the hidden common space. Graph shrinkage is accomplished by dynamically and synchronously warping each image towards its connected neighbors in the graph. To cater to the hierarchical structure of the graph, we apply graph shrinkage sequentially for different levels. It is worth noting that more levels may be used to better represent the image distribution if needed, but will increase the computation time significantly. For our datasets, we found that using two levels is sufficient for reasonable representation of data, as well as for reasonable registration accuracy. More specifically, in the low level, each subgroup shrinks as all member images agglomerate. In the high level, the representative images of the subgroups are used to drive each subgroup closer to the common space. Graph shrinkage is executed alternatively throughout different levels until all images are registered to the common space. In what follows, details for both constructing the hierarchical graph and executing graph shrinking are provided.
Construction of Hierarchical Graph
Graph construction is key to the success of eHUGS. As shown by Ying et al. [23], the graph needs to satisfy three criteria: (1) all images must be connected by the graph to guarantee that all images can be warped to a common space; (2) the number of edges in the graph needs to be minimized to reduce the computation load and also make graph shrinkage efficient; (3) the graph needs to capture the image distribution sufficiently well to ensure good registration accuracy.
We employ the hierarchical graph to capture the distribution of heterogeneous image data more accurately. Our graph reflects the relatively small variation within each subgroup as well as the large variation between subgroups. This is inspired by the work of Wang et al. [29], where the entire dataset of images are clustered into several subgroups based on the similarity between images. We follow a similar approach and derive the hierarchical graph based on the outcome of image clustering. First, the whole dataset of images are separated into several subgroups, each comprising similar images. Second, a representative image is selected from each subgroup to produce a high-level graph.
As illustrated in Fig 2, a set of N images I = {Ii|i = 1, …, N} reside in a manifold M. The distance di,j for each pair of images Ii and Ij is defined as the sum of squared differences (SSD) di,j = ‖Ii − Ij‖2. Note that there are many other choices for defining image similarity/distance, such as mutual information (MI) [30] and normalized mutual information (NMI) [31]. Our choice of SSD is simply for the sake of computational efficiency. Images with greater similarity tend to be spatially close to each other on the manifold. Since it is easier to register images with similar anatomical structures, the graph is constructed by arranging similar images into the common subgroups. The construction of the graph is accomplished by data-driven clustering as detailed next.
[Figure omitted. See PDF.]
Fig 2. Illustration of images residing on a high dimensional manifold and connected via the geodesic paths.
The geodesic distances between image Ii and other images are shown by dash arrows.
https://doi.org/10.1371/journal.pone.0146870.g002
Data-Driven Image Clustering.
A data-driven approach is used to cluster images into different subgroups based on their distances on the manifold. Image pairs with small distances tend to be clustered into the same subgroup, while images with large distances tend to be clustered into different subgroups. Methods such as k-means clustering [32] and affinity propagation [26] can be used. We opt for affinity propagation because, unlike k-means clustering, it does not require the explicit specification of the number of clusters.
An affinity matrix of similarities between all pairs of images is required by affinity propagation. We first apply affine registration to align all images to a common space, and then populate the N × N affinity matrix S whose elements si,j are defined as the negative SSD, i.e., si,j = −di,j = −‖Ii − Ij‖2. We also set the preference value p in affinity propagation as below, which controls the likelihood of an image as a subgroup exemplar (representative image) and also eventually influences the number of exemplars:(1)
In general, small preference value results in the small number of representative images. Since we do not have any priori knowledge on the observed images, we use the same value of p for all images, so that all images have equal possibility of becoming the representative image in each subgroup.
We further construct a N × N symmetric connection matrix E to encode the graph. If image Ii and image Ij is connected, we set the connection matrix E’s element ei,j = 1; otherwise, ei,j = 0. All diagonal entries in E are set to 0 so that no image is self-connected.
Specifically, we apply the following steps to build the graph:
1. Cluster the N images I = {Ii|i = 1, …, N} into Ω subgroups Gα (α ∈ {1, …, Ω} and ) by using affinity propagation.
2. Determine the global center image Io of the dataset,(2)to approximate the hidden common space.
3. Within each subgroup Gα, a representative image is determined as the one nearest to the global center image Io. Thus, with Ω subgroups, we have Ω representative images: . Particularly, when the subgroup includes the global center image Io, its representative image which is the nearest to the global center image Io is the global center image Io itself. Thus, we have . The advantage of choosing the image nearest to the global center for delegating the whole subgroup is that this (representative) image can make all other images in the whole subgroup deform easily towards the global center.
4. For each subgroup Gα, intra-group connections are made between image and all other member images. For any pair of images Ii and Ij ∈ Gα with j ≠ i, we have ei,j = 1 if i = iα or j = iα, and ei,j = 0 otherwise.
5. For the high-level connections between subgroups, the global center image Io is connected to all other representative images that are determined in Step (3). For any pair of images and j ≠ i, the conclusion can be drawn that ei,j = 1 if i = o or j = o, and ei,j = 0 otherwise. An example graph is shown in Fig 1.
We can assume that all images within each subgroup are structurally similar. Thus, the image closest to the global center can be chosen as the representative image to guide deformations of all images in the subgroup to the global center image.
Graph-Shrinking-Guided Image Normalization
We assure that all images are connected, and therefore the graph has a total of N − 1 edges. Based on this graph, all images in the manifold can be warped in accordance to their connected images. The registration process can then be viewed as a dynamic graph shrinkage process, as described in our previous HUGS (Hierarchical Unbiased Graph Shrinkage) for groupwise registration [23]. A brief description of HUGS is given here.
In general, we can assume the deformation of each individual image as a dynamic procedure of time variable t. Thus, Ii(t) can be used to represent the deformed image Ii at time t. Next, we use a graph defined in the brain image manifold as follows. Let be the graph nodes and E = {eij: i,j = 1, …, N} be the edges between two nodes in the graph. eij = 1 means existence of a link between Ii(t) and Ij(t). Otherwise, there is no direct link between Ii(t) and Ij(t) in the graph. Also, we define a N × N weighted adjacency matrix where each element exp (vij(t)) describes the geodesic pathway between two images, such that vij(t) > 0 if eij = 1 and vij(t) = ∞ otherwise. Note that vij(t) denotes the velocity vector of geodesic pathway, which indicates the distance between Ii(t) and Ij(t). Here, deformable image registration [12] is performed to estimate each velocity vector vij(t) and further calculate the deformation pathway between Ii(t) and Ij(t) by exp (vij(t)), where ‘exp’ is the exponential map [33]. Intuitively, the goal of our graph-based groupwise registration is to minimize the velocity vectors along all graph edges, as defined below:(3)
The principle behind F(t) is demonstrated in Fig 3. First, all images are assumed to be sitting in a high-dimensional manifold. Then, the topology of their distribution is described by a graph, where the graph edges denote the local connectivity between graph nodes. Specifically, the velocity vector vij(t) is associated with each graph edge, where the integration along vij(t) forms the geodesic distance from Ii(t) to Ij(t). Thus, the minimization of F(t) can be regarded as a dynamic graph shrinking procedure, which deforms each image from Ii(t) to Ii(t + Δt) with the decreased overall geodesic distance, while keeping the topology of the entire graph. As time t increases, all Ii(t)s are supposed to meet at the population center, with properly determined velocity vector vij(t) and time increment Δt.
[Figure omitted. See PDF.]
Fig 3. Registration via graph shrinkage.
At time t, all the images in the graph (denoted by blue points) are connected to each other by the geodesic path (denoted by blue dash lines). At time t+Δt, all the images are deformed according to their connected images through geodesic paths (denoted by purple dash lines). The deformed images (denoted by red points) become closer to each other and a new graph describing their relationship is generated (denoted by red dash lines).
https://doi.org/10.1371/journal.pone.0146870.g003
As we formulate the problem of groupwise registration as the dynamic shrinkage of graph, it is critical to determine the deformation of each image Ii(t) at time t, which can minimize energy function F(t). Suppose each image has been deformed from Ii(t0) to Ii(tk), where {tk} is the discretization of time t (k = 0, … K, t0 = 0, tK → ∞). Regarding the local connectivity of each node Ii(tk) in the graph, it is reasonable to move Ii(tk) along the average direction according to its connected nodes. Since the velocity vector sits on the tangent space of Ii(tk) on the manifold, it can be efficiently calculated by linear averaging as , where is the number of connections for Ii(tk). Given the direction of each velocity vector at time point tk, it can be proven that the entire energy function F(t) is strictly and monotonously decreased if each node Ii(tk) moves along the direction of velocity vector from tk to tk + Δtk, where the time increment Δtk is bounded by . Eventually, the geodesic pathway φi from each image Ii to the population center can be obtained by concatenating the deformation segments from t0 to tK, i.e., , where ‘∘’ denotes the deformation composition [34]. During graph shrinkage, all images are gradually deformed to the hidden common space with preservation of graph topology. For numerical implementation, we discretize t using tk (k = 0,1,2,…, with t0 = 0 and tk → ∞ as k → ∞) based on the interval Δt. Then we use as the new deformed image at time tk, and also as the composited deformation field of image at time tk.
Summary of Our Method
Our method is briefly summarized as follows:
1. Perform affine registration to obtain a set of N affine-aligned images I = {Ii|i = 1, …, N}.
2. Compute all distances di,j = ‖Ii − Ij‖2, ∀i,j ∈ {1, …, N}.
3. Cluster images I into Ω subgroups Gα () by using affinity propagation.
4. Construct a null connection matrix E with elements ei,j = 0, ∀i,j ∈ {1, …, N}.
5. Obtain the global center image Io by Eq 2.
6. Determine in each subgroup Gα the image nearest to the global center image Io as the representative image .
7. Establish intra-subgroup connections in the same subgroup by connecting all images Ii {Ii ∈ Gα and i ≠ iα} with the representative image . Then, establish inter-subgroup connections by connecting all representative images and the global center image Io.
8. Perform registration via graph shrinkage iteratively until convergence.
Results and Discussion
eHUGS offers an efficient and accurate solution to tackle heterogeneous image datasets (We found eHUGS did not achieve significant improvement over HUGS on small dataset such as IXI (30 subjects) and NIREP (16 subjects) since single graph is sufficient to describe the image distribution for small number of subjects). We evaluate the proposed method on both adult and infant datasets. Specifically, we tune the parameters of our eHUGS registration method using the 4-time-point longitudinal brain MR scans of 20 elderly normal subjects, obtained from the ADNI dataset [27]. Then, using these tuned parameters, we apply eHUGS to the infant dataset, which includes 285 MR brain images of infant subjects with 2 weeks to 2 years of age. The proposed method is compared with the standard group-mean method [25], ABSORB [17], and HUGS [35].
Experiments on Elderly Brains
In this experiment, 20 normal subjects are chosen from the ADNI database. In addition to the baseline scan, each subject was also scanned after 6 months, 1 year, and 2 years. Sample baseline scans are shown in Fig 4. This dataset is preprocessed using the following steps. First, anterior commissure (AC) posterior commissure (PC) correction is applied to all images, which are further resampled to size 256×256×256 with voxel resolution 1mm×1mm×1mm. After applying the N3 algorithm [36] for intensity inhomogeneity correction, Brain Surface Extractor (BSE) [37] and Brain Extraction Tool (BET) [38] are used for skull stripping. Next, each image is segmented into three types of tissues: white matter (WM), gray matter (GM), and cerebrospinal fluid (CSF) by using FAST [39] in FSL software package (http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/). We further visually inspect and manually correct the segmentation results. Hence, we can regard these segmentation results as the ground truth in evaluating the registration accuracy. Finally, affine registration using FLIRT [28] is applied.
[Figure omitted. See PDF.]
Fig 4. Image samples of elderly subjects obtained from the ADNI database.
In this figure, large variations across subjects can be observed.
https://doi.org/10.1371/journal.pone.0146870.g004
In eHUGS, the distance matrix is calculated using the pre-processed brain images and then affinity propagation is applied to cluster all the images. In this way, we obtain 20 subgroups based on the observed image dataset. The clustering result is shown in Fig 5 (A), where all images have been projected to a 3-dimensional space. Each dimension of the projected space denotes one of the 3 largest eigenvectors obtained by performing Principal Component Analysis (PCA) on the image distribution. Both intra- and inter-subgroup connections are also shown in Fig 5 (B).
[Figure omitted. See PDF.]
Fig 5. Clustering and graph construction.
(a) Images are clustered into 20 subgroups, represented with different markers and colors, by using affinity propagation. (b) Intra-subgroup (black color) and inter-subgroup (orange color) connections in the graph.
https://doi.org/10.1371/journal.pone.0146870.g005
Evaluations of eHUGS, HUGS, ABSORB, and the standard group-mean method are performed on a Dell workstation with 2 Xeon E74850, 2.0Ghz, 10 core CPU, and 256GB/1066Mhz RAM. The computation times of all methods are summarized in Table 1. There are in total 79 edges in the graph constructed in eHUGS, compared with 1265 edges in the graph given by HUGS. This causes HUGS to be nearly 20 times slower than eHUGS. The standard group-mean method gives 79 edges, and ABSORB gives more than 79 edges. Generally, in terms of computation time, eHUGS ≅ group-mean < ABSORB < 1-week computation time< HUGS (~1 month).
[Figure omitted. See PDF.]
Table 1. Computation times of the standard group-mean method, ABSORB, HUGS, and eHUGS on the elderly brain image dataset, obtained from ADNI.
https://doi.org/10.1371/journal.pone.0146870.t001
To quantitatively measure registration accuracy, the Dice ratio is used to measure the degree of overlap of a specific type of tissue between different images after registration. The Dice ratio is defined as(4)
Here, | ∙ | is the voxel count within the region-of-interest. A and B denote the regions spanned by a specific type of tissue, respectively, in the two images. Better registration should yield higher Dice ratios [40]. Since no explicit template image is selected in groupwise registration, the ground-truth segmentation in the common space is obtained by majority voting from all aligned images. Then, the Dice ratio associated with each tissue type (WM, GM, and CSF) is computed with respect to the estimated ground-truth segmentation.
Before (non-rigid) groupwise registration, the average Dice ratios of WM, GM, CSF, and all tissues in overall (after affine registration) are 70.01±2.99%, 55.20±3.00%, 60.64±2.81%, and 61.95±2.93%, respectively. Table 2 shows the average Dice ratios of each tissue type given by the standard group-mean method, ABSORB, HUGS, and eHUGS, respectively. The lower quartiles, medians, and upper quartiles of Dice ratios on WM, GM and CSF are also shown in Fig 6(A)–6(C), respectively, for different registration approaches. It can be observed that our eHUGS achieves the highest Dice ratio for each tissue type. Since the conventional group-mean registration method independently registers each individual image to the population center, such outlier images undermine the group mean and also mislead the registration result of whole group. On the contrary, our graph-based groupwise registration method adaptively deforms each image based on the neighbouring images in the high-dimensional image manifold, with the goal to improve the overall registration accuracy while minimize the influence from the outlier images.
[Figure omitted. See PDF.]
Fig 6. The box plots of Dice ratios of WM, GM and CSF for the elderly brain dataset (ADNI).
(a) Dice ratio produced by the standard group-mean method. (b) Dice ratio produced by ABSORB. (c) Dice ratio produced by HUGS. And (d) Dice ratio produced by eHUGS.
https://doi.org/10.1371/journal.pone.0146870.g006
[Figure omitted. See PDF.]
Table 2. Dice ratios of WM, GM, CSF and all tissues, obtained by the standard group-mean method, ABSORB, HUGS, and eHUGS, respectively, for the elderly brain dataset (ADNI).
https://doi.org/10.1371/journal.pone.0146870.t002
The average Dice ratios of WM, GM and CSF obtained by the proposed method (eHUGS), HUGS, ABSORB, and the standard group-mean method in each iteration are shown in Fig 7(A)–7(D), respectively. For each method, the Dice ratio becomes stable after a certain number of iterations. From Fig 7(A)–7(D), we can observe that the average Dice ratios for WM, GM, CSF and all tissues in overall given by eHUGS are higher than any other method. We apply the paired t-test between our eHUGS and all other three counterpart methods, and our eHUGS achieves significant improvement over all three methods in terms of Dice ratio with p < 0.01. More strictly, we further apply Mann-Whitney U-test, after the false discovery rate (FDR) correction, and also find that the improvements by eHUGS over all other three groupwise registration methods are still statistically significant (p < 0.01) in all tissue types, as indicated by the red ‘*’ in Table 2.
[Figure omitted. See PDF.]
Fig 7. Changes of Dice ratios with progress of groupwise registration for the elderly brain dataset (ADNI).
https://doi.org/10.1371/journal.pone.0146870.g007
It is worth noting that our graph-based groupwise registration is very flexible to work with other similarity measurements. Here, we also use mutual information (MI) to compute the pairwise similarity and construct the hierarchical graph. The Dice ratios in WM, GM, and CSF after groupwise registration using MI are 88.54%, 77.24%, and 82.49%, respectively, which are comparable to the results using SSD (last row in Table 2).
Experiments on 285 Infant Brains (Aged from Birth to 2-year-old)
In this experiment, we evaluate the groupwise registration result on 95 longitudinal infant image series. Each infant subject has been scanned at birth, 1 year old, and 2 years old. All these images were acquired using a Siemens head-only 3T scanner for 95 subjects. Due to WM myelination, we use T2-weighted images for neonates and T1-weighted images for 1–2 year(s) old infants. For neonatal subjects (scanned at birth), we use T2-weighted images with the following parameters: TR = 7380ms, TE = 119ms, flip angle = 150°, acquisition matrix = 256×128, and voxel resolution = 1.25×1.25×1.95mm3. For 1-year-old and 2-year-old subjects, we use T1-weighted images, with the following imaging parameters: TR = 1900ms, TE = 4.38ms, flip angle = 7°, acquisition matrix = 256×192, and voxel resolution = 1×1×1mm3. T2-weighted images are affine aligned onto the respective 1-year-old T1-weighted images of the same subject and further resampled to 1×1×1mm3. All images are cropped to the size 256×256×198. Then skull stripping, inhomogeneity correction are performed to all T1- and T2-weighted MR images. Since we have each infant subjects scanned in 0-year-old, 1-year-old, and 2-year old, we employ the state-of-the-art infant segmentation method [41] to obtain reasonable segmentation results for infant images by letting the segmentation of later time points guide the challenging segmentation in the 0 year old. We further manually inspect and correct the segmentation result, in order to guarantee the segmentation quality in evaluating the registration accuracy. After that, we perform FLIRT in FSL package to linearly register all infant subjects to the common space. As shown in Fig 8, structural variations across different infants are large, and also the appearances across different times of same infant are also large.
[Figure omitted. See PDF.]
Fig 8. Sample images of two subjects from the infant database.
We can observe large structural variation across subjects, and also large appearance variation across different images of same subject.
https://doi.org/10.1371/journal.pone.0146870.g008
In eHUGS, again the distance matrix is calculated for the pre-processed brain images, and then affinity propagation is applied to cluster the images, resulting in 2 subgroups. As can be seen from Fig 8, the images of 1-year-olds and 2-year-olds are relatively similar, especially when compared with their differences with those of neonates. This is consistent with the fact that the total brain volume increases 101% in the first year of life, in contrast to only 15% in the second year (Knickmeyer, et al., 2008). The clustering result is shown in Fig 9 (A), where all images are projected to a 3-dimensional space. Again, each dimension denotes one of the 3 largest eigenvectors by performing PCA on the image distributions. The intra- and inter-subgroup connections are also shown in Fig 9 (B).
[Figure omitted. See PDF.]
Fig 9. Clustering and graph construction.
(a) Images are clustered into 2 subgroups, represented with different markers and colors, by using affinity propagation. (b) Intra-subgroup (black color) and inter-subgroup (orange color) connections in the graph.
https://doi.org/10.1371/journal.pone.0146870.g009
There are in total 284 edges in the graph constructed by eHUGS, compared with 24,494 edges in the graph given by HUGS. This causes HUGS to be nearly 84 times slower than eHUGS. In fact, for this dataset, based on our estimation, HUGS needs more than a year to complete the registration. In this case, we only compare eHUGS with the standard group-mean method and ABSORB. For all methods, the computation times are summarized in Table 3.
[Figure omitted. See PDF.]
Table 3. Computation times of the standard group-mean method, ABSORB, HUGS, and eHUGS for the infant brain image dataset.
https://doi.org/10.1371/journal.pone.0146870.t003
Before (non-rigid) groupwise registration, the average Dice ratios of WM, GM, CSF, and all tissues in overall are 60.93±9.15%, 72.29±11.11%, 39.19±7.57%, and 57.47±9.28%, respectively. Table 4 shows the average Dice ratio of each tissue type. The lower quartiles, medians, and upper quartiles of Dice ratios on WM, GM and CSF are shown in Fig 10(A)–10(C), respectively, for different registration approaches. It can be observed that eHUGS again achieves the highest Dice ratio in all comparisons.
[Figure omitted. See PDF.]
Fig 10. The box plots of the Dice ratios of WM, GM and CSF for the infant dataset.
(a) By the standard group-mean method. (b) By ABSORB. (c) By eHUGS.
https://doi.org/10.1371/journal.pone.0146870.g010
[Figure omitted. See PDF.]
Table 4. Dice ratios of WM, GM, CSF and all tissues, obtained by the standard group-mean method, ABSORB, and eHUGS, respectively, for the infant dataset.
https://doi.org/10.1371/journal.pone.0146870.t004
The average Dice ratios of WM, GM and CSF obtained by the proposed method (eHUGS), the standard group-mean method, and ABSORB in each iteration are shown in Fig 11(A)–11(D). From Fig 11(A)–11(D), we can draw similar conclusions as for the elderly brain dataset. That is, the average Dice ratios on WM, GM, CSF and the overall by the proposed method are higher than any other method. Under paired t-test, we find that our eHUGS achieves significant improvements (p < 0.01) over group-mean, and ABSORB in terms of Dice ratio on all tissue types. By considering all 285 aligned images together, we also perform the Mann-Whitney U test between the Dice ratios of eHUGS and all other methods. After FDR correction, the multiple comparison tests show that eHUGS achieves significant improvement in CSF with p < 0.01, as indicated by the red ‘*’ in Table 4. Similarly, we also replace the SSD measurement with MI in constructing the graph. The final groupwise registration results by using two image similarity measurements are comparable, where the tissue overlap ratios in WM, GM, and CSF are (77.44%, 78.88%, 52.79%) by using SSD and (77.96%, 78.64%, 52.96%) by using MI.
[Figure omitted. See PDF.]
Fig 11. Changes of Dice ratios with progress of groupwise registration for the infant brain dataset.
https://doi.org/10.1371/journal.pone.0146870.g011
Since the infant atlas is very important in many brain development studies, we specifically evaluate the registration accuracy by inspecting the quality of the infant atlases, constructed from all aligned infant images. Here, we use a simple averaging to get the mean image, and use majority voting to generate the tissue probability maps for WM, GM, and CSF, respectively. From top to bottom, in Fig 12, we show the constructed atlas (mean image, WM, GM, and CSF probability maps) by the standard group-mean method, ABSOB, and eHUGS. Through visual inspection, the atlases constructed by eHUGS show greater anatomical details than any other method.
[Figure omitted. See PDF.]
Fig 12. Illustration of the atlases of infant brains.
A-C are obtained by the standard group-mean method, ABSORB, and eHUGS, respectively. In each panel, from left to right are the mean image, and tissue probability maps for CSF, GM and WM, respectively.
https://doi.org/10.1371/journal.pone.0146870.g012
Conclusion
In this paper, we have presented an effective and efficient method for unbiased groupwise registration of large dataset of images. Specifically, we propose an enhanced version of HUGS, namely eHUGS, to model the complex image distribution by using the hierarchical graph and thus overcome the limitation of using a single graph in HUGS. Experiments using both elderly and infant brain images indicate that eHUGS gives the improved registration performance. In all experiments, eHUGS yields the best Dice ratio when compared with the standard group-mean method, ABSORB, and HUGS.
In our proposed eHUGS, our main idea is to use hierarchical graph for representing image distribution, and then register images progressively to each other for avoiding direct registration between different images. For our current datasets, using simple image similarity metric, we found that 2 levels are enough to represent the image distribution. But, for other datasets, more levels may be needed for representing more complex image distribution. On the other hand, when affine aligning all images to the comment space, although we select a population center as a template, pairwise registration is still used. Potentially, we can use groupwise affine registration for alignment of the whole dataset, which could improve the registration results of all affine-aligned images. All these will be our future work.
Besides, our other future work also includes (1) evaluating the registration performance of eHUGS on large datasets; (2) applying eHUGS to other applications, such as group comparison for discovering AD imaging biomarkers.
Author Contributions
Conceived and designed the experiments: GW XP SY QW DS DGS. Performed the experiments: GW XP SY QW DS DGS. Analyzed the data: GW XP SY QW DS DGS. Contributed reagents/materials/analysis tools: GW XP SY QW DS DGS. Wrote the paper: GW XP SY QW PY DS DGS.
Citation: Wu G, Peng X, Ying S, Wang Q, Yap P-T, Shen D, et al. (2016) eHUGS: Enhanced Hierarchical Unbiased Graph Shrinkage for Efficient Groupwise Registration. PLoS ONE 11(1): e0146870. https://doi.org/10.1371/journal.pone.0146870
1. Toga AW, Thompson PM. The role of image registration in brain mapping. Image and Vision Computing. 2001;19(1–2):3–24. pmid:19890483
2. Toga AW, Thompson PM. TEMPORAL DYNAMICS OF BRAIN ANATOMY. Annual Review of Biomedical Engineering. 2003;5:119–45. pmid:14527311
3. Zitova B. Image registration methods: a survey. Image and Vision Computing. 2003;21(11):977–1000. pmid:citeulike:494149.
4. Auzias G, Colliot O, Glaunes JA, Perrot M, Mangin JF, Trouve A, et al. Diffeomorphic brain registration under exhaustive sulcal constraints. IEEE Trans Med Imaging. 2011;30(6):1214–27. pmid:21278014.
5. Jia H, Yap P-T, Wu G, Wang Q, Shen D. Intermediate templates guided groupwise registration of diffusion tensor images. NeuroImage. 2011;54(2):928–39. pmid:20851197
6. Hill DLG, Batchelor PG, Holden M, Hawkes DJ. Medical image registration. Physics in Medicine and Biology. 2001;46:R1–R45. pmid:11277237
7. Durrleman S, Fillard P, Pennec X, Trouve A, Ayache N. Registration, atlas estimation and variability analysis of white matter fiber bundles modeled as currents. Neuroimage. 2011;55(3):1073–90. pmid:21126594.
8. Fonov V, Evans AC, Botteron K, Almli CR, McKinstry RC, Collins DL. Unbiased average age-appropriate atlases for pediatric studies. NeuroImage. 2010;54(1):313–27. pmid:20656036
9. Shen D, Davatzikos C. HAMMER: Hierarchical attribute matching mechanism for elastic registration. Medical Imaging, IEEE Transactions on. 2002;21(11):1421–39.
10. Beg F, Miller M, Trouve A, Younes L. Computing large deformation metric mappings via geodesic flows of diffeomorphisms. International Journal of Computer Vision. 2005;61(2):139–57. pmid:Beg05.
11. Avants B, Epstein CL, Grossman M, Gee JC. Symmetric diffeomorphic image registration with cross-correlation: Evaluating automated labeling of elderly and neurodegenerative brain. Medical Image Analysis. 2008;12(1):26–41. pmid:17659998
12. Vercauteren T, Pennec X, Perchant A, Ayache N. Symmetric log-domain diffeomorphic registration: a demons-based approach. Medical Image Computing and Computer-Assisted Intervention–MICCAI 20082008. p. 754–61.
13. Vercauteren T, Pennec X, Perchant A, Ayache N. Diffeomorphic demons: efficient non-parametric image registration. NeuroImage. 2009;45(1, Supplement 1):S61–S72. pmid:19041946
14. Cootes TF, Twining CJ, Petrovic VS, Babalola KO, Taylor CJ. Computing accurate correspondences across groups of images. IEEE Trans Pattern Anal Mach Intell. 2010;32(11):1994–2005. pmid:20847389.
15. Joshi S, Davis B, Jomier M, Gerig G. Unbiased diffeomorphic atlas construction for computational anatomy. Neuroimage. 2004;23 Suppl 1:S151–60. pmid:15501084.
16. Park H, Bland PH, Hero AO, Meyer CR. Least biased target selection in probabilistic atlas construction. Medical Image Computing and Computer-Assisted Intervention–MICCAI 20052005. p. 419–26.
17. Jia H, Wu G, Wang Q, Shen D. ABSORB: Atlas building by self-organized registration and bundling. NeuroImage. 2010;51(3):1057–70. pmid:20226255
18. Learned-Miller EG. Data driven image models through continuous joint alignment. Pattern Analysis and Machine Intelligence, IEEE Transactions on. 2006;28(2):236–50.
19. Spiclin Z, Likar B, Pernus F. Groupwise Registration of Multimodal Images by an Efficient Joint Entropy Minimization Scheme. IEEE Transactions on Image Processing. 2012;21(5):2546–58. pmid:22294031
20. Wachinger C, Navab N. Simultaneous Registration of Multiple Images: Similarity Metrics and Efficient Optimization. IEEE Trans on Pattern Anal and Mach Intell. 2013;35(5):1221–33.
21. Wu G, Wang Q, Jia H, Shen D. Feature-based Groupwise Registration by Hierarchical Anatomical Correspondence Detection. Human Brain Mapping. 2012;33(2):253–71. pmid:21391266
22. Wu G, Jia H, Wang Q, Shen D. SharpMean: Groupwise Registration Guilded by Sharp Mean Image and Tree-based Registration. NeuroImage. 2011;56(4):1968–81. pmid:21440646
23. Ying S, Wu G, Wang Q, Shen D. Hierarchical unbiased graph shrinkage (HUGS): A novel groupwise registration for large data set. Neuroimage. 2014;84:626–38. pmid:24055505.
24. Orchard J, Mann R. Registering a multisensor ensemble of images. IEEE Transactions on Image Processing. 2010;19(5):1236–47. pmid:20040419
25. Joshi S, Davis B, Jomier M, Gerig G. Unbiased diffeomorphic atlas construction for computational anatomy. pmid:NeuroImage. 2004;23(Supplement 1):S151–S60.
26. Frey BJ, Dueck D. Clustering by passing messages between data points. Science. 2007;315:972–6. pmid:17218491
27. ADNI. http://www.loni.ucla.edu/ADNI/.
28. Jenkinson M, Smith S. A global optimisation method for robust affine registration of brain images. Medical Image Analysis. 2001;5(2):143–56. pmid:11516708
29. Wang Q, Chen L, Yap PT, Wu G, Shen D. Groupwise registration based on hierarchical image clustering and atlas synthesis. Human Brain Mapping. 2010;31(8):1128–40. pmid:20063349.
30. Pluim JPW, Maintz JBA, Viergever MA. Mutual-information-based registration of medical images: a survey. IEEE Transactions on Medical Imaging. 2003;22(8):986–1004. pmid:12906253
31. Studholme C, Hill DLG, Hawkes DJ. An overlap invariant entropy measure of 3D medical image alignment. Pattern Recognition. 1999;32(1):71–86. pmid:WOS:000078298000006.
32. MacQueen J, editor Some methods for classification and analysis of multivariate observations. Proceedings of the fifth Berkeley symposium on mathematical statistics and probability; 1967: California, USA.
33. Helgason S. Differential Geometry, Lie Groups and Symmetric Space: Academic Press; 2001.
34. Ashburner J. A fast diffeomorphic image registration algorithm. NeuroImage. 2007;38(1):95–113. pmid:17761438
35. Ying S, Wu G, Wang Q, Shen D. Hierarchical unbiased graph shrinkage (HUGS): a novel groupwise registration for large data set. NeuroImage. 2014;84(1):626–38.
36. Sled JG, Zijdenbos AP, Evans AC. A nonparametric method for automatic correction of intensity nonuniformity in MRI data. Medical Imaging, IEEE Transactions on. 1998;17(1):87–97.
37. Shattuck DW, Sandor-Leahy SR, Schaper KA, Rottenberg DA, Leahy RM. Magnetic resonance image tissue classification using a partial volume model. NeuroImage. 2001;13(5):856–76. pmid:11304082
38. Smith SM. Fast robust automated brain extraction. Human Brain Mapping. 2002;17(3):143–55. pmid:12391568
39. Zhang Y, Brady M, Smith S. Segmentation of brain MR images through a hidden Markov random field model and the expectation maximization algorithm. IEEE Transactions on Medical Imaging. 2001;20(1):45–57. pmid:11293691
40. Rohlfing T. Image Similarity and Tissue Overlaps as Surrogates for Image Registration Accuracy: Widely Used but Unreliable. IEEE Transactions on Medical Imaging. 2012;31(2):153–60. pmid:21827972
41. Shi F, Shen D, Yap P-T, Fan Y, Cheng J, An H, et al. CENTS: Cortical Enhanced Neonatal Tissue Segmentation. Human Brain Mapping. 2011;32(2):382–96.
You have requested "on-the-fly" machine translation of selected content from our databases. This functionality is provided solely for your convenience and is in no way intended to replace human translation. Show full disclaimer
Neither ProQuest nor its licensors make any representations or warranties with respect to the translations. The translations are automatically generated "AS IS" and "AS AVAILABLE" and are not retained in our systems. PROQUEST AND ITS LICENSORS SPECIFICALLY DISCLAIM ANY AND ALL EXPRESS OR IMPLIED WARRANTIES, INCLUDING WITHOUT LIMITATION, ANY WARRANTIES FOR AVAILABILITY, ACCURACY, TIMELINESS, COMPLETENESS, NON-INFRINGMENT, MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. Your use of the translations is subject to all use restrictions contained in your Electronic Products License Agreement and by using the translation functionality you agree to forgo any and all claims against ProQuest or its licensors for your use of the translation functionality and any output derived there from. Hide full disclaimer
© 2016 Wu et al. This is an open access article distributed under the terms of the Creative Commons Attribution License: http://creativecommons.org/licenses/by/4.0/ (the “License”), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.
Abstract
Effective and efficient spatial normalization of a large population of brain images is critical for many clinical and research studies, but it is technically very challenging. A commonly used approach is to choose a certain image as the template and then align all other images in the population to this template by applying pairwise registration. To avoid the potential bias induced by the inappropriate template selection, groupwise registration methods have been proposed to simultaneously register all images to a latent common space. However, current groupwise registration methods do not make full use of image distribution information for more accurate registration. In this paper, we present a novel groupwise registration method that harnesses the image distribution information by capturing the image distribution manifold using a hierarchical graph with its nodes representing the individual images. More specifically, a low-level graph describes the image distribution in each subgroup, and a high-level graph encodes the relationship between representative images of subgroups. Given the graph representation, we can register all images to the common space by dynamically shrinking the graph on the image manifold. The topology of the entire image distribution is always maintained during graph shrinkage. Evaluations on two datasets, one for 80 elderly individuals and one for 285 infants, indicate that our method can yield promising results.
You have requested "on-the-fly" machine translation of selected content from our databases. This functionality is provided solely for your convenience and is in no way intended to replace human translation. Show full disclaimer
Neither ProQuest nor its licensors make any representations or warranties with respect to the translations. The translations are automatically generated "AS IS" and "AS AVAILABLE" and are not retained in our systems. PROQUEST AND ITS LICENSORS SPECIFICALLY DISCLAIM ANY AND ALL EXPRESS OR IMPLIED WARRANTIES, INCLUDING WITHOUT LIMITATION, ANY WARRANTIES FOR AVAILABILITY, ACCURACY, TIMELINESS, COMPLETENESS, NON-INFRINGMENT, MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. Your use of the translations is subject to all use restrictions contained in your Electronic Products License Agreement and by using the translation functionality you agree to forgo any and all claims against ProQuest or its licensors for your use of the translation functionality and any output derived there from. Hide full disclaimer