A multi-resolution scheme for distortion-minimizing mapping between human subcortical structures based on geodesic construction on Riemannian manifolds

Abstract

In this paper, we deal with a subcortical surface registration problem. Subcortical structures including hippocampi and caudates have a small number of salient features such as heads and tails unlike cortical surfaces. Therefore, it is hard, if not impossible, to perform subcortical surface registration with only such features. It is also non-trivial for neuroanatomical experts to select landmarks consistently for subcortical surfaces of different subjects. We therefore present a landmark-free approach for subcortical surface registration by measuring the amount of mesh distortion between subcortical surfaces assuming that the surfaces are represented by meshes. The input meshes can be constructed using any surface modeling tool available in the public domain since our registration method is independent of a surface modeling process. Given the source and target surfaces together with their representing meshes, the vertex positions of the source mesh are iteratively displaced while preserving the underlying surface shape in order to minimize the distortion to the target mesh. By representing each surface mesh as a point on a high-dimensional Riemannian manifold, we define a distance metric on the manifold that measures the amount of distortion from a given source mesh to the target mesh, based on the notion of isometry while penalizing triangle flipping. Under this metric, we reduce the distortion minimization problem to the problem of constructing a geodesic curve from the moving source point to the fixed target point on the manifold while satisfying the shape-preserving constraint. We adopt a multi-resolution framework to solve the problem for distortion-minimizing mapping between the source and target meshes. We validate our registration scheme through several experiments: distance metric comparison, visual validation using real data, robustness test to mesh variations, feature alignment using anatomic landmarks, consistency with previous clinical findings, and comparison with a surface-based registration method, LDDMM-surface.

Introduction

Shape analysis of human subcortical structures has received much attention to detect functional abnormalities of human brains in early stages of brain diseases. Reduction of hippocampus volumes was observed for subjects with Alzheimer's disease (Wang et al., 2006, Apostolova et al., 2006, Devanand et al., 2007), and structural change of the basal ganglia was reported for subjects with schizophrenia (Mamah et al., 2007, Ellison-Wright et al., 2008). Inter-subject correspondences between anatomical features play a key role in shape analysis on human subcortical structures as in cortical surface analysis (Grenander and Miller, 1998, Fischl et al., 1999, Joshi et al., 2004, Gholipour et al., April 2007, Miller and Qiu, 2009). Individual human subcortices can be represented in the form of point sets, surfaces, or volumes in their local coordinate systems (Rettmann et al., 2002, Mangin et al., 2004, Wang et al., 2005, Van Essen, 2005).

Intuitively, volume deformation is accompanied by surface deformation. However, the neuroanatomical significance of subcortical surface deformity has not been fully explored unlike that of cortical surface deformity. In other words, there have been few attempts to explore it since the spatial resolution of MR images has not been high enough to detect the regional surface deformity of subcortical structures. Tools for their surface deformity analysis have not been mature enough, either. However, subcortical structures such as basal nuclei and thalami have consistent connection with cortices. Thus, deformations of certain cortical regions could influence on localized regions of subcortical structures, which may result in subcortical surface deformation in a disease-specific way. In subcortical structure studies, volume-based approaches have dominated surface-based approaches. It is not until recently that surface-based approaches have appeared. For example, Csernansky et al. (2004) reported that local deformity of thalamic surfaces showed more significant effect than volume deformity in diagnosis of schizophrenia. Local surface deformity analysis on hippocampi demonstrated that surface deformity provided neuroanatomical significance (Thompson et al., 2004, Wang et al., 2006). In particular, Becker et al., 2006, Frisoni et al., 2008 reported that surface-based morphometry is more reliable for hippocampus shape analysis than tensor-based morphometry relying on volumetric registration. Our choice of representation is surface meshes to further support surface-based approaches along this line. Subcortical structures including hippocampi and caudates have a small number of salient features such as heads and tails unlike cortical surfaces. Therefore, it is a challenging problem to establish inter-subject surface correspondences of human subcortical structures with only such features (Shi et al., 2009).

In this paper, we adapt the multi-resolution framework (Kilian et al., 2007) on Riemannian manifolds to our problem setting in order to present a multi-resolution scheme for distortion-minimizing diffeomorphic mapping between human subcortical surfaces without using predefined landmarks. The source and target surface meshes are first preprocessed so that they have the same number of vertices and the same connectivity. Our scheme then iteratively displaces the vertex positions of the source mesh at each level in a hierarchical fashion in order to minimize the distortion toward the target mesh. By representing source and target surface meshes as their respective points on a high-dimensional Riemannian manifold, we reduce the distortion minimization problem to the problem of finding the point on the manifold that is closest to the target point. The source mesh is allowed to move on the source surface, and the target mesh is fixed on the target surface. To estimate mesh distortion, we develop a distortion metric to measure the deviation from an isometric deformation while also penalizing triangle flipping. Using this metric, the amount of mesh distortion from the (current) source mesh to the target mesh is estimated as the total distortion along the deformation path between the source and the target points on the manifold. We minimize the distortion between the two meshes by constructing a geodesic curve from the source point to the target point on the manifold. To do this, we present a multi-resolution framework on the high-dimensional Riemannian manifold. This scheme can robustly construct the inter-surface mapping that naturally aligns the anatomical features of the subcortical surfaces. The efficacy of the proposed mapping method is demonstrated using several experiments, including the anatomic landmark-based validation and a group analysis on clinical data.

Compared to the work in Kilian et al. (2007), the contributions of our work in the methodological aspect are three-fold: First, we deal with a surface registration problem for subcortical structures, where the source mesh is allowed to move on a source surface. In the work of Kilian et al. (2007), however, both the source and target meshes are fixed during optimization. Specifically, our objective is to find the isomorphic mesh on the source surface, which minimizes the amount of distortion toward the (target) mesh on the target surface. On the other hand, their objective is to interpolate the fixed source and target meshes in order to generate a sequence of in between meshes. Second, we propose a new metric to measure the amount of distortion between the two meshes, which aims at preserving the isometry between the (final) source mesh to the target mesh while avoiding triangle flipping. Our metric consists of two terms: an isometry term and a Laplacian preservation term. Originated from their metric, the former term is to preserve the inter-mesh isometry as much as possible. Unlike their corresponding term, it prevents excessively high penalty to long edges. The latter is to avoid triangle flipping, which is newly introduced to our metric. Third, we present a new multi-resolution solver for our problem, which is quite different from their solver because of the constraint that the source mesh is allowed to move on the source surface.

From the practical aspect, the contributions of our work are two-fold: First, our method is independent of a surface modeling process. Subcortical surface meshes can be constructed from three-dimensional T1 images by using any surface modeling tool available in the public domain. With these meshes as input data, our method performs surface registration without any assumption on the meshes except that the meshes represent genus-zero surfaces. Second, our method requires no predefined landmarks, such as feature points or curves on surfaces. Since it is non-trivial for neuroanatomical experts to select landmarks consistently for subcortical surfaces of different subjects, a landmark-free approach is suitable for large-scale group analysis on subcortical structures.

In order to prepare the source and target meshes, we collected a group of 49 subjects composed of 26 healthy controls and 23 patients with Alzheimer's disease (AD). The subjects underwent high-resolution T1-weighted volume magnetic resonance imaging (MRI) at the Samsung Medical Center, Seoul, South Korea. Specifically, three-dimensional T1-weighted spoiled-gradient (SPGR) echo images were acquired using a 1.5-T MRI scanner (GE Sigma, Milwaukee, WI, USA) with the following imaging parameters: coronal slice thickness = 1.5 mm; echo time = 7 ms; repetition time = 30 ms; the number of excitations = 1; flip angle = 45^∘; field of view = 22 × 22 cm; and matrix size = 256 × 256. The healthy controls had no history of neurological or psychiatric abnormalities. Also, the cognitive functioning of the healthy subjects was confirmed to be within normal limits as assessed by a series of neuropsychological tests including Mini-Mental State Examination (MMSE) (Folstein et al., November 1975). The AD patients met NINCDS/ADRDA (McKhann et al., 1984) criteria for probable AD, in which diagnostic procedures included a clinical interview, a neurological examination, and a series of neuropsychological tests. Laboratory tests including complete blood count, blood chemistry, vitamin B12/foliate, syphilis serology, and thyroid function tests were also performed to verify the absence of secondary causes of the dementia in any of the patients. We obtained informed consents from all the subjects, and the study was approved by the Institutional Review Board of the Samsung Medical Center. Table 1 presents the demographic characteristics of the participants.

The T1 images of each subject were processed to extract a subcortical surface. The surface extraction method consists of three parts: image parcellation, manual correction, and surface mesh generation. First, the anatomical parcellations of every subcortical structure were generated from the MRI images of each individual subject using version 5.0 of the FreeSurfer software package (Athinoula A. Martinos Center at the Massachusetts General Hospital, Harvard Medical School; http://www.surfer.nmr.mgh.harvard.edu/). The boundaries of resulting regions may be inaccurate and may have topological errors such as holes (Qiu and Miller, 2008, Khan et al., 2008). The second step is to refine those regions by performing manual correction for the parcellated images: A neuroanatomist manually edited the automatically parcellated regions slice-by-slice using version 3.6.1 of 3D Slicer (Pieper et al. (2006)). The boundary of each parcellated region in a binary image was smoothed using a low pass filter. The kernel size was set to be small (3 ⁎ 3 ⁎ 3 mm³) to smooth out noisy features along the boundary while preserving the original boundary shape. Finally, we converted the smoothed image to a binary image and generated subcortical surface meshes from the boundaries of the parcellated region using the FreeSurfer software package. Note that although we employed this sequence of surface extraction steps for our experiments, our registration method is independent of a surface modeling process.

Many methods have been proposed for establishing inter-subject correspondences of human brain models. From an algorithmic viewpoint, methods for brain registration can be categorized into two groups: volume registration and surface registration (Gholipour et al., April (2007)). With advancement of cortical surface extraction techniques, surface-based methods are becoming increasingly popular. For our purpose, we focus on reviewing surface-based methods. For volume-based methods, we refer readers to a survey paper (Gholipour et al. (2007)). Surface-based methods exploit sulcal landmarks that are defined on a cortical surface to guide the registration process (Towle et al., 2003, Fischl et al., 2004). For example, manually identified surface landmarks such as a set of labeled cortical sulci (Van Essen et al., 1998a, Zhong and Qiu, 2010), sulcal lines (Gu et al., 2003, Durrleman et al., 2008) and a set of sulcal points (Glaunès et al., 2004) were used for cortical surface registration. Surface features such as mean curvature, convexity, or sulcal depth have also been utilized (Van Essen et al., 1998b, Fischl et al., 1999, Robbins, 2004, Tosun and Prince, 2005, Lyttelton et al., 2007, Yeo et al., 2010). However, subcortical structures have a small number of salient features such as the hippocampal head and tail unlike cortical surfaces. Therefore, the cortex registration techniques based on local surface features cannot be directly applied to subcortical surface registration.

The geometry of subcortical surfaces can be represented using spherical harmonic (SPHARM) basis functions (Brechbühler et al., 1995). Gerig et al. (2001) proposed a registration method by minimizing the mean squared distance between corresponding points of two SPHARM-based surfaces, which was later extended by Huang et al. (2007) to present a fast surface alignment algorithm. Shi et al. (2009) recently presented a subcortical surface registration technique by exploiting more general harmonic basis functions on arbitrary manifolds defined by eigenvectors of the Laplace–Beltrami(LB) operator (Reuter et al., 2006, Levy, 2006, Qiu et al., 2008). The large deformation diffeomorphic metric mapping (LDDMM) has been a popular surface registration method for both subcortical surfaces (Vaillant and Glaunès, 2005, Qiu and Miller, 2008, Qiu et al., 2009, Miller and Qiu, 2009) and regions of interest (ROI) on a cortical surface (Vaillant et al., 2007, Qiu et al., 2007). The LDDMM method measures the quality of a mapping based on surface normals (Vaillant and Glaunès, 2005), and its usefulness was validated for neuroanatomical analysis (Vaillant et al., 2007, Qiu et al., 2007, Zhong et al., 2010). In particular, for subcortical surface analysis, the LDDMM method was successfully employed as a preprocess step (Qiu and Miller, 2008, Qiu et al., 2009, Miller and Qiu, 2009). However, none of these methods (Vaillant and Glaunès, 2005, Huang et al., July 2007, Shi et al., 2009) explicitly addressed the distance preservation issue during registration.

The importance of distance preserving mapping is well-known for various applications (Timsari and Leahy, 2000, Mémoli and Sapiro, 2004, Bronstein et al., 2005, Lipman and Funkhouser, 2009, Balasubramanian et al., 2010). However, many distance preserving methods have focused on mapping of surfaces to their flattened ones on the plane (Timsari and Leahy, 2000, Sander et al., 2001, Zigelman et al., 2002, Balasubramanian et al., 2010). For example, a cortical surface was flattened by exploiting an area and angle preserving mapping (Timsari and Leahy, 2000) or an edge length preserving mapping (Balasubramanian et al., 2010). The inter-surface mapping method in Schreiner et al. (2004) presented a distance preserving mapping between two arbitrary surfaces of the same topology. However, this method requires a set of predetermined point correspondences to establish an initial coarse mapping. Eckstein et al. (2007) presented a deformable matching framework by exploiting a pseudo-Hausdorff distance to explicitly define matching errors. Although this method employed an isometry-based metric to measure distance distortion during deformation, their optimization scheme does not guarantee to provide an optimal solution in the sense that the deformation may not necessarily follow geodesic paths under the proposed metric.