Home | About | Journals | Submit | Contact Us | Français |

**|**HHS Author Manuscripts**|**PMC3571860

Formats

Article sections

Authors

Related links

IEEE Trans Med Imaging. Author manuscript; available in PMC 2013 February 13.

Published in final edited form as:

Published online 2011 September 15. doi: 10.1109/TMI.2011.2168233

PMCID: PMC3571860

NIHMSID: NIHMS405276

Yalin Wang,^{1,}^{2} Xianfeng Gu,^{3} Tony F. Chan,^{1} Paul M. Thompson,^{2} and Shing-Tung Yau^{4}

Yalin Wang: ude.alcu.htam@gnawly; Xianfeng Gu: ude.bsynus.sc@ug; Tony F. Chan: ude.alcu.htam@nahc; Paul M. Thompson: ude.alcu.inol@nospmoht; Shing-Tung Yau: ude.dravrah.htam@uay

The publisher's final edited version of this article is available at IEEE Trans Med Imaging

See other articles in PMC that cite the published article.

In medical imaging, parameterized 3D surface models are of great interest for anatomical modeling and visualization, statistical comparisons of anatomy, and surface-based registration and signal processing. By solving the Yamabe equation with the Ricci flow method, we can conformally parameterize a brain surface via a mapping to a multi-hole disk. The resulting parameterizations do not have any singularities and are intrinsic and stable. To illustrate the technique, we computed parameterizations of cortical surfaces in MRI scans of the brain. We also show the parameterization results are consistent with constraints imposed on the mappings of selected landmark curves, and the resulting surfaces can be matched to each other using constrained harmonic maps. Unlike previous planar conformal parameterization methods, our algorithm does not introduce any singularity points.

Surface-based modeling is valuable in brain imaging to help analyze anatomical shape, to statistically combine or compare 3D anatomical models across subjects, and to map and compare functional imaging parameters localized on anatomical surfaces. Parameterization of these surface models involves computing a smooth (differentiable) one-to-one mapping of regular 2D coordinate grids onto the 3D surfaces, so that numerical quantities can be computed easily from the resulting models [1]. The mesh-based work contrasts with implicit methods, which typically define a surface as the level set of a higher-dimensional function [2]. Relative to level set methods, surface meshes can allow regular 2D grids to be imposed on complex structures, transforming a difficult 3D problem into a 2D planar problem, with simpler data structures, discretization schemes, and rapid data access and navigation. Here we present a new method to parameterize brain surfaces based on algebraic functions. We find a planar conformal parameterization without any singularities by solving the Yamabe equation with the Ricci flow method. This method can compute conformal invariants of brain surfaces which can be used to compare and classify brain surface structures. Compared with previous brain conformal parametrization work [3, 4], the parameterization provided by our algorithm does not have any zero points so there is less area distortion. By solving a harmonic map in the parameter domain, our algorithm provides smooth correspondence fields for matching of different brain surfaces while explicitly matching labeled sets of landmark curves.

Brain surface parameterization has been studied intensively. Schwartz et al. [5], and Timsari and Leahy [6] computed quasi-isometric flat maps of the cerebral cortex. Drury et al. [7] presented a multiresolution method for flattening the cerebral cortex. Hurdal and Stephenson [8] report a discrete mapping approach that uses circle packings to produce “flattened” images of cortical surfaces on the sphere, the Euclidean plane, and the hyperbolic plane. The maps obtained are quasi-conformal approximations of classical conformal maps. Haker et al. [9] implement a finite element approximation for parameterizing brain surfaces via conformal mappings. They select a point on the cortex to map to the north pole of the Riemann sphere and conformally map the rest of the cortical surface to the complex plane by stereographic projection of the Riemann sphere to the complex plane. Gu et al. [10] propose a method to find a unique conformal mapping between any two genus zero manifolds by minimizing the harmonic energy of the map. They demonstrate this method by conformally mapping a cortical surface to a sphere. Ju et al. [11] present a least squares conformal mapping method for cortical surface flattening. Joshi et al. [12] propose a scheme to parameterize the surface of the cerebral cortex by minimizing an energy functional in the *p ^{th}* norm. Wang et al. [3, 4] have used holomorphic 1-forms to parameterize anatomical surfaces with complex (possibly branching) topology. Recently, Ju et al. [13] reported the results of a quantitative comparison of FreeSurfer [14], CirclePack, and least squares conformal mapping (LSCM) with respect to geometric distortion and computational speed.

Suppose *M* is a surface embedded in ^{3}, then it has the natural induced Euclidean metric, denoted by **g**. Suppose is another Riemannian metric on *M*, we say it is conformal to **g**, if the two metrics differ by a scalar function *u* : *M* → , namely = *e*^{2u}
**g**.

For the purpose of brain mapping, we typically want to flatten the cortical surface onto the plane, with a specific set of anatomical landmarks are mapped to specific locations, such as circles in the flattened space;, furthermore, the mapping is required to be conformal. We formulate this problem as finding a conformal metric that induces the prescribed curvature, such that all interior points have zero Gaussian curvature, and the boundary points have constant geodesic curvature. This can be formulated rigorously in terms of the Yamabe equation [15],

$$\{\begin{array}{lll}\stackrel{\sim}{K}& =& 0\\ \mathrm{\Delta}u-K+{e}^{2u}\stackrel{\sim}{K}& =& 0\\ {k}_{\stackrel{\sim}{g}}\mid \partial M& =& \mathit{const}\end{array}$$

(1)

The Yamabe equation can be solved using the Ricci flow method,

$$\frac{\mathit{du}(t)}{\mathit{dt}}=\stackrel{\sim}{K}-K(t).$$

(2)

It has been proven that surface Ricci flow with normalized total area will converge to the desired metric, and the convergence is exponentially fast.

In practice, all surfaces are approximated by discrete piecewise polygonal surfaces so we developed a discrete Ricci flow method that applies to triangulated meshes. We associate a circle of radius *γ _{i}* with vertex

$${l}_{\mathit{ij}}=\sqrt{{\gamma}_{i}^{2}+{\gamma}_{j}^{2}+2\phantom{\rule{0.2em}{0ex}}\mathrm{cos}\phantom{\rule{0.2em}{0ex}}{\varphi}_{\mathit{ij}}\phantom{\rule{0.2em}{0ex}}{\gamma}_{i}\phantom{\rule{0.2em}{0ex}}{\gamma}_{j}}$$

Therefore, the circle radii and the intersection angles determine a metric on the mesh, which is called the *circle packing metric*. By using the circle packing metric, the conformal deformation can be approximated by changing the vertex radii while preserving the intersection angles. The discrete Ricci flow can be defined as

$$\frac{{\mathit{d\gamma}}_{i}}{\mathit{dt}}=-({\stackrel{\sim}{K}}_{i}-{K}_{i}){\gamma}_{i,}$$

(3)

where is the target curvature at the vertex *υ _{i}*.

Let *u _{i}* = ln

$$f(\mathbf{u})={\mathit{\int}}_{{\mathbf{u}}_{0}}^{\mathbf{u}}\sum _{i=1}^{n}{K}_{i}{\mathit{du}}_{i}.$$

(4)

where **u**_{0} = (0, 0, …, 0). The integration path is chosen arbitrarily, namely, the energy is path independent. Because scaling does not affect the curvature, we confine the **u** to lie in the sub-linear space
$\sum {u}_{i}=0$. This energy is strictly convex in this space, therefore it has a unique global minimum. This global minimum is exactly the desired metric for the prescribed curvature. The discrete Ricci flow 3 is the negative gradient flow of the discrete Ricci energy.

In current work, we propose to use the Newton method which is quadratically convergent, and much faster than previous Ricci flow methods. The Newton method requires computation of the Hessian matrix of the energy,

$${H}_{\mathit{ij}}=\frac{\partial {K}_{i}}{\partial {u}_{j}}.$$

The algorithm is equivalent solving Equation 1 that describes a conformal deformation. We use the Ricci flow method [16] to solve this equation.

**Algorithm 1**
*Conformal Mapping to a Multi-hole Punctured Disk*

Input: mesh M, step length ε, energy difference threshold δK;

*Output*:
$\overrightarrow{h}:M\to D$. *Here D R*^{2}, *and D is a multi-hole disk*.

*Computing initial radii γ*${l}_{\mathit{ij}}={\gamma}_{i}^{2}+{\gamma}_{j}^{2}-2{\gamma}_{i}\phantom{\rule{0.2em}{0ex}}{\gamma}_{j}\mathit{cos}{\varphi}_{\mathit{ij}}$._{i}for each vertex, and angle ϕ_{ij}for each edge e_{ij}, such that*Compute boundary loops, denoted as*Γ_{0}, Γ_{1}, …, Γ_{n}.*The*Γ_{0}*is the exterior boundary*.*Set target Gaussian curvature of each interior vertex to be zero,*= 0._{i}*For any vertex on υ*Γ_{k}_{0},*set its target Gaussian curvature to*${\stackrel{\sim}{K}}_{k}=\frac{2\pi}{\left|{\mathrm{\Gamma}}_{i}\right|}$,*where*|Γ_{0}|*denotes the number of vertices in*Γ_{0}.*For any vertex on υ*Γ_{k}≠ 0,_{i}, i*set its target Gaussian curvature to*${\stackrel{\sim}{K}}_{k}=-\frac{2\pi}{\left|{\mathrm{\Gamma}}_{i}\right|}$,*where*|Γ_{i}|*denotes the number of vertices in*Γ_{i}.*Update the vertex radii with the Ricci flow*,$${\gamma}_{i}(t+1)={\gamma}_{i}(t)+\epsilon \times ({\stackrel{\sim}{K}}_{i}(t)-{K}_{i}(t))\times {\gamma}_{i}(t).$$*Update the target Gaussian curvature for boundary vertices, suppose υ*Γ_{k}_{i}, suppose e_{k−1,k},*e*_{k,k+1}Γ= Σ_{i}, then let S_{i}_{epq Γi}*l*, then_{pq}$${\stackrel{\sim}{K}}_{k}=\frac{{l}_{k-1,k}+{l}_{k,k+1}}{2{S}_{i}}\times {C}_{i},$$*where C*${C}_{i}=\{\begin{array}{cc}2\pi ,& i=0\\ -2\pi ,& i\ne 0\end{array}$._{i}is*Repeat step 6 and 7 until the maximal Gaussian curvature error, max*|_{i}*K*−_{i}|,_{i}*is less than δK*.

After the computation of conformal parameterizations for open boundary genus zero surfaces with a multiple-hole punctured disk, we can compute the direct correspondence of two surfaces by solving a constrained harmonic mapping problem [4]. Given two surfaces *S*_{1} and *S*_{2}, their punctured disk parameterizations are *τ*_{1} : *S*_{1} → *R*^{2} and *τ*_{2} : *S*_{2} → *R*^{2}, we want to compute a map, *ϕ* : *S*_{1} → *S*_{2}. Instead of directly computing of *ϕ*, we can easily find a harmonic map between the parameter domains. We look for a harmonic map, *τ* : *R*^{2} → *R*^{2}, such that *τ* *τ*_{1}(*S*_{1}) = *τ*_{2}(*S*_{2}), *τ* *τ*_{1} (*S*_{1}) = *τ*_{2}(*S*_{2}), Δ*τ* = 0. Then the map *ϕ* can be obtained by
$\varphi ={\tau}_{1}\circ \tau \circ {\tau}_{2}^{-1}$. Since *τ* is a harmonic map while *τ*_{1} and *τ*_{2} are conformal map, the resulting *ϕ* is a harmonic map.

We tested our algorithm with cortical surfaces extracted from 3D MRI scans of the brain and we also tested its ability to accommodate constraints with different landmark sets. Specifically, each of the rows in Figure 1 shows a cortical left hemisphere cortex labeled by four, seven and twelve landmarks, respectively. After cutting along these landmark curves, the cortical surface becomes an open boundary high genus surface. Our algorithm conformally maps the surface to a 3-hole disk (first row), 6-hole disk (second row) and 11-hole disk (third row). The perimeter of the corpus callosum is mapped to the exterior circular disk boundary and other landmark curves are mapped to the disk’s inner circle boundaries.

Illustrates conformal maps of the same brain with various landmark sets. Each of the three rows shows the brain with four, seven and twelve landmarks, respectively. The landmark curves are labeled by thick blue lines. The last column show their parameterization **...**

Figure 2 illustrates how our algorithm is used to match two left hemisphere cortical surfaces. As shown in Figure 2(a), (b), (d) and (e), we selected four major landmark curves on two different cortices, to illustrate the approach (thick lines show the precentral and postcentral sulci, the superior temporal sulcus, and the corpus callosum boundary at the mid-sagittal plane). By cutting the surface along these landmark curves, we obtain two genus-3 open-boundary surfaces. Figure 2(c) and (f) show their conformal map to a 3-hole disk. Because of the shape difference between two cortices, the centers and the radii of inner circles are different. By computing a constrained harmonic map from (f) to (c), we have a new parameterization (Figure 2(g)) of the cortex in the second row ((d) and (e)). The inner circle centers and radii of the new parameterization are identical to the parameterization in (c). With the new 3-hole disk as the canonical space, we can easily compute a direct surface correspondence between two surfaces (a) and (d). Because the inner circles and exterior circle are identical for the two parameterizations, the landmark curves lying in the surface are exactly matched to each other. Figure 2 (h)-(m) illustrate the direct surface correspondence by morphing between these two cortical surfaces. Figure 2(h) and (m) and surfaces (a) and (d) respectively, viewed from a different viewpoint. (j), (k) and (l) are the intermediate shapes when linearly interpolating the surface correspondence vector field between (h) and (m). We can see that although surface shape changes substantially from (h) to (m), the relative locations of the three landmarks remain the same. This verifies that our algorithm provides a method to perform surface matching, while explicitly matching sulcal curves or, potentially, any other landmarks lying in the surface. Although some previous approaches [3, 4] had the same motivation as ours, because their work introduced singularities at the so-called zero points, so their surface matching results have some errors and inevitable distortions in the areas around the zero points. Our approach provides an improved method for global surface matching with exact landmark matching capability.

In this paper, we presented a brain surface conformal parameterization method based on algebraic functions. With the Ricci flow method, we solved the Yamabe equation to obtain a conformal deformation that conformally maps open boundary surfaces to multi-hole disk. We tested our algorithm on the hippocampus and surface models of the cerebral cortex, including major cortical sulci as anatomical landmark constraints. Used as a canonical space, the multi-hole disk conformal parameterization provides a brain surface matching approach that can exactly match landmark curves lying on the surfaces. Compared with other work that conformally maps brain surfaces to parallelograms, our algorithm offers some advantages because it does not introduce any singular points. Our future work will include empirical application of the Ricci flow concept to medical applications in computational anatomy, including the detection of population differences and the tracking of brain change over time.

This work was funded by the National Institutes of Health through the NIH Roadmap for Medical Research, Grant U54 RR021813 entitled Center for Computational Biology (CCB). The work was performed while the third author was on leave at the National Science Foundation as Assistant Director of the Directorate for Mathematics & Physical Sciences.

1. Thompson PM, Giedd JN, Woods RP, MacDonald D, Evans AC, Toga AW. Growth patterns in the developing human brain detected using continuum-mechanical tensor mapping. Nature. 2000 Mar;404(no. 6774):190–193. [PubMed]

2. Osher S, Sethian JA. Fronts propagating with curvature dependent speed: Algorithms based on Hamilton-Jacobi formulations. J Comput Phys. 1988;79:12–49.

3. Wang Y, Gu X, Hayashi KM, Chan TF, Thompson PM, Yau ST. MICCAI. II. Palm Springs; CA, USA: 2005. Brain surface parameterization using Riemann surface structure; pp. 26–29.

4. Wang Y, Gu X, Hayashi KM, Chan TF, Thompson PM, Yau ST. Surface parameterization using Riemann surface structure. 10th IEEE ICCV. 2005:1061–1066.

5. Schwartz EL, Shaw A, Wolfson E. A numerical solution to the generalized mapmaker’s problem: Flattening nonconvex polyhedral surfaces. IEEE TPAMI. 1989 Sep;11(no. 9):1005–1008.

6. Timsari B, Leahy RM. An optimization method for creating semi-isometric flat maps of the cerebral cortex. Proceedings of SPIE, Medical Imaging; San Diego, CA. Feb 2000.

7. Drury HA, Van Essen DC, Anderson CH, Lee CW, Coogan TA, Lewis JW. Computerized mappings of the cerebral cortex: A multiresolution flattening method and a surface-based coordinate system. J Cognitive Neurosciences. 1996;8:1–28. [PubMed]

8. Hurdal MK, Stephenson K. Cortical cartography using the discrete conformal approach of circle packings. NeuroImage. 2004;23:S119–S128. [PubMed]

9. Angenent S, Haker S, Tannenbaum A, Kikinis R. Conformal geometry and brain flattening. Med Image Comput Comput -Assist Intervention. 1999 Sep;:271–278.

10. Gu X, Wang Y, Chan TF, Thompson PM, Yau S-T. Genus zero surface conformal mapping and its application to brain surface mapping. IEEE TMI. 2004 Aug;23(no. 8):949–958. [PubMed]

11. Ju L, Stern J, Rehm K, Schaper K, Hurdal MK, Rottenberg D. IEEE ISBI. Arlington, VA, USA: 2004. Cortical surface flattening using least squares conformal mapping with minimal metric distortion; pp. 77–80.

12. Joshi AA, Leahy RM, Thompson PM, Shattuck DW. IEEE ISBI. Arlington, VA, USA: 2004. Cortical surface parameterization by p-harmonic energy minimization; pp. 428–431. [PMC free article] [PubMed]

13. Ju L, Hurdal MK, Stern J, Rehm K, Schaper K, Rottenberg D. Quantitative evaluation of three surface flattening methods. NeuroImage. 28(no. 4):869–880. 2005. [PubMed]

14. Fischl B, Sereno MI, Dale AM. Cortical surface-based analysis II: Inflation, flatteningm and a surface-based coordinate system. NeuroImage. 1999;9:179–194. [PubMed]

15. Schoen R, Yau S-T. Lectures on Differential Geometry. International Press; 1994.

16. Hamilton RS. The Ricci flow on surfaces. Contemp Math. 71:237–262. 1988.

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's national science library in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |