Group Title: Department of Computer and Information Science and Engineering Technical Reports
Title: An Intrinsic geometric framework for simultaneous non-rigid registration and segmentation of surfaces
Full Citation
Permanent Link:
 Material Information
Title: An Intrinsic geometric framework for simultaneous non-rigid registration and segmentation of surfaces
Alternate Title: Department of Computer and Information Science and Engineering Technical Report ; TR-06-001
Physical Description: Book
Language: English
Creator: Lord, N.
Ho, J.
Vemuri, Baba C.
Publisher: Department of Computer and Information Science and Engineering, University of Florida
Place of Publication: Gainesville, Fla.
Publication Date: April 18, 2006
Copyright Date: 2006
 Record Information
Bibliographic ID: UF00095637
Volume ID: VID00001
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.


This item has the following downloads:

2006367 ( PDF )

Full Text

An Intrinsic Geometric
Framework for Simultaneous
Non-Rigid Registration and
Segmentation of Surfaces *

N. Lord, J. Ho, and B.C. Vemuri
e-mail: [nal, jho, vemuri]

UF-CISE Technical Report TR-06-001
Computer and Information Sciences Department
University of Florida
Bldg. CSE, Room 301
Gainesville, Florida 32611
April 18, 2006

* This research was supported in part by the grant NIH-R01-NS04812 to BCV and the
UF Alumni Fellowship to NAL. Data sets were provided by Dr. Christiana Leonard
of the UF McKnight Brain Institute.

Abstract. In clinical applications where structural asymmetries between
homologous shapes have been correlated with pathology, the questions
of definition and quantification of 'asymmetry' arise naturally. When not
only the degree but the position of deformity is thought relevant, asym-
metry localization must also be addressed. Asymmetries between paired
shapes can and have been formulated in the literature in terms of (non-
rigid) diffeomorphisms between the shapes. For the infinity of such maps
possible for a given pair, we define optimality as the minimization of
total distortion, where 'distortion' is in turn defend as deviation from
isometry. We thus propose a novel variational formulation for segmenting
asymmetric regions from surface pairs based on the minimization of a
functional of both the deformation map and the segmentation boundary,
which controls gradient discontinuity of the map. This minimization is
achieved via a quasi-simultaneous evolution of the map and curve. Our
formulation is inherently intrinsic and parameterization-independent. We
present examples using both synthetic data and pairs of left and right
hippocampal structures, hippocampus malformation being linked with
such neurological disorders as !il and schizophrenia.

1 Introduction

Many neurological disorders have been clinically correlated with shape abnormal-
ities in specific brain structures. For example, hippocampal shrinkage is associ-
ated with epilepsy and schizophrenia, among other conditions. A precise analysis
of such abnormalities, which are usually considered to be failures of symmetry
of the left and right halves of the structure in question, could unlock the abil-
ity not only to track the progression of existent disease, but to identify at-risk
individuals for preventative treatment. However, devising a framework in which
to conduct the necessary shape comparison is a nontrivial matter. The short-
comings of raw volumetric comparison are obvious. Landmark-based methods
are inapplicable to problems involving many anatomical structures (including
the hippocampus), because many biological structures do not exhibit the sorts
of readily and consistently identifiable local features on which such methods de-
pend. The medial surface representation of Gerig et al. [1][2] represents a more
sophisticated take on volumetric analysis that includes some localization of vol-
ume disparity, but the economy of the method leads both to noise sensitivity and
an acknowledged inability to distinguish between different sorts of deformation
(e.g. certain forms of elongation vs. bending). This is due to the expression of
complicated behaviors in terms of simple distances which discard some amount
of directionality.
More recently, Wang et al.[3][4][5] presented a technique involving the align-
ment of surfaces based on the formation of a 2D parametricc) diffeomorphic
map from Riemannian surface structure information, with this map representing
the non-rigid registration between the surfaces. This approach arguably makes
use of what surface I. ,ilre" information is there to be had from relatively
nondescript shapes, without imposing an ill-defined demand for landmark iden-
tification. However, in forming the correspondence by maximizing the mutual

information of the chosen surface characteristics over the entirety of both sur-
faces, these methods ignore key considerations. For one, while the conceptual
division of the body into anatomical structures is hardly arbitrary in a general
sense, for the purpose of deformation analysis, it may very well be. That is to
say, we may wish to identify a selected portion of the hippocampus as struc-
turally abnormal relative to its correspondent, rather than regarding the entire
hippocampus as deformed. Furthermore, the MI of the domain and target surface
features does not in and of itself restrict the differential properties of the map,
only what properties should be found at given locations on the target surface.
Therefore two regions of the target surface that exhibit the same Riemannian
characteristics cannot be differentiated between, despite the fact that one choice
may be physically unreasonable as a map target. Thus this framework demands
a separate regularizing constraint on the map evolution.

To address these concerns, we present a novel approach to the non-rigid reg-
istration of anatomically correspondent surfaces which rests on the fundamental
assumptions that (a) the relative deformity of anatomical structures for which
symmetry is expected can be intuitively and precisely quantified as the devia-
tion from isometry of the deformation map between their surfaces, and (b) since
deformities may well be locally confined to certain regions of the given surfaces,
the evolution of the global correspondence must allow for a partial disconnect
between normal and abnormal regions, as these are by definition not expected
to exhibit the same deformation patterns. As such, surface segmentation rests
at the heart of our approach, which can be viewed not only as a registration
tool but also as an asymmetry localizer and quantifier. Our definition of system
energy in terms of isometry inherently balances the consideration that the map
should accurately match Riemannian surface characteristics with the idea that
the map should incur as little deformation as possible in doing so. This gives
way to a simple, elegant formulation relatively free of heuristics. Since we here
restrict our consideration to surfaces that are topologically cylindrical (appropri-
ate for hippocampi among many other structures), we are able to parameterize
with a single patch (periodic at a pair of ends) for each surface in the pair, and
as the formulation is intrinsic, we can confine evolution to these 2D parametric
domains, retaining elegance. The evolution is driven quasi-simultaneously by the
energies of both the segmenting curve and the deformation field, in the spirit of
the works of Yezzi et al. and Unal et al.[6] [7], where the former group restricted
its work to rigid image registration and the latter generalized the idea to the
non-rigid deformation case. In [8], Vemuri et al. presented a registration-assisted
image segmentation technique. Of the many proposals in the literature which
followed, however, none addressed the issue of joint manifold segmentation and
registration, which is our focus. Our work differs notably in that we evolve the
segmenting curve and the registration map on a general 2D manifold (rather than
a flat image plane), and rather than registering and segmenting image data, we
register and segment the shape characteristics of the very surface on which the
processes are operating.

The rest of the paper is organized as follows: Section 2 contains the mathe-
matical formulation and numerical algorithm used to solve it. Section 3 contains
experimental results on both synthetic and real data sets. In section 4, we draw

2 Simultaneous Segmentation and Registration

Let S1 and S2 be two surfaces in IR3. The Euclidean metric in IR3 induces Rie-
mannian metrics gl and g2 on S1 and S2, respectively. The goal of the algorithm
is to segment regions in S1 and their corresponding regions in S2 such that the
metric structures of the corresponding regions match optimally. More specifi-
cally, the algorithm computes a homeomorphism f between Si and S2 and a
set of closed curves 71 on Si. Let 72 denote the closed curve in S2 that is the
image of i7 under f. The restriction of f to the complement Sl\7y is a diffeo-
morphism between S1\7y and S2\y2. If U1, -, U, and V1, - V, denote the
collections of open components of S1\7i and S2\72, respectively, then f maps
Ui diffeomorphically onto Vi for each i and f matches the Riemannian structures
between that of Ui and of Vi = f(U). See Figure 1. We solve the simultaneous

St -. s2

Is ii

Fig. 1. Left: The two given surfaces, S1 and S2. Right: The desired result of simul-
taneous segmentation and registration of these two surfaces. The map f is a homeo-
morphism between S1 and S2. f maps the segmenting curve 71 on Si to a curve 72
on S2, and f establishes correspondences between connected components of Si\7i and
S2\72. U1 and U2 are two disjoint open sets in Si\7i, and V1 and V2 are two disjoint
open sets in S2\72. f maps Ui diffeomorphically onto V for each i. The restriction of
f to Ui optimally matches the Riemannian structures of Ui and V.

segmentation and registration problem outlined above using a variational frame-
work. The energy functional E is defined as a functional of a pair (f, y7), where
71 is a set of closed curves on S1 and f : S1 S2 is a homeomorphism which is
C on Sl\yi. Let f*g2 denote the pull-back metric (the first fundamental form)

on gi, and S(f, i7) is defined as

'(f, ~) I |f*g2 -g|2dA+ df(lQ(t))/ /' ./ (1)
JS1\7i .
The second integral computes the length of the segmenting curve 72 on S2, while
the first integral computes the matching (registration) cost between the two Rie-
mannian structures on the complement, Si\7i. By computing the length of the
segmenting curve on S2 instead of Si, the energy functional tightly couples the
two somewhat disparate processes, segmentation and registration. The integrand
of the first term above, If*g2 gfl2, is the usual L-2norm between the two ten-
sors g1 and f*g.9 In local coordinates, it is given as the usual L2-norm between
the two matrices Gi and JtG2J:

IJtG2J Gi 9, (2)

where J is the Jacobian of f expressed in the local coordinates system, and Gi
and G2 are the usual 2 x 2 positive definite matrices expressing the metric tensors
gi and g2 in the given local coordinates on S1 and the induced one (using f) on
At this point, a comparison between our approach and the usual Mumford-
Shah framework [9] is called for. Given an image I, the energy functional for
Mumford-Shah is

MS(I' ) j -1I' 12dA+ VII|dA+ [ (t)Idt.

In the above, D denotes the image and 7 the segmenting curve on D. The third
term, which measures the length of 7, corresponds directly to the second term
of C in 1. The first two terms of Mumford-Shah, representing data fidelity and
smoothness, respectively, may seem to lack analogous entries in E. However,
we argue below that the first integral in the definition of C serves both as the
fidelity and smoothing term for I, and we can thus consider C as a generalization
of the usual Mumford-Shah functional to the problem of simultaneous shape
registration and segmentation.
For our problem, we want to compare two surfaces B1 and B2 via their
respective Riemannian structures gi and g2, and to extract a meaningful seg-
mentation from the comparison. Instead of an image I' as in the usual Mumford-
Shah framework, our variable is the map f between Si and S2. In our case, the
data fidelity requirement appears as the penalization of f for failure to match
the Riemannian structures of the surfaces, while the smoothing term appears
as the penalization of f for large first-order derivative magnitudes. A close ex-
amination of Equation 1 reveals that the first integral indeed contains both
requirements. First, the norm between the two tensors f*g9 and gi as defined
in Equation 2 clearly measures the quality of match between the two metric
tensors. Second, the Jacobian J contains all of the first-order derivatives of f,
and the squared-norm in Equation 2 indirectly measures the magnitudes of the

first-order derivatives of f using both metrics gl and g2. In sum, we can interpret
the first term in Equation 1 as the integration of the first two terms in the usual
Mumford-Shah functional MS. In MS, these two terms are necessarily separate
because there is no natural connection between matching pixel intensities and
restricting the gradient VI. In our problem, the connection between matching
and smoothing is natural: since we are matching Riemannian structures, which
are tensors defined on tangent spaces, any matching between these two struc-
tures has to involve some kind of (implicit) matching between tangent spaces,
and the appearance of derivatives in defining the matching is then unavoidable.

2.1 Minimizing the Energy Functional

In implementing the energy minimization as defined above, we rely on alter-
nating iterations between a Chan-Vese level set curve segmenter [10], extended
to the general 2D manifold domain, which produces a step towards the best
estimate of the segmentation of the energy field given the current map, and a
sparse Hessian-exploiting quasi-Newton optimization scheme acting on the non-
rigid map between the parametric domains, which provides a step towards the
minimum-energy map between the surfaces as constrained by the current po-
sition of the segmenting curve. The manner in which the curve separates the
evolution occurring on its exterior from that occurring on its interior relates to
a Neumann condition imposed across the boundary, akin to what is suggested
in the Mumford-Shah curve evolution implementation in [11]. The curve length
minimizing term implementation is straightforward, though curvature computa-
tion must take nonuniform surface length elements into account. Hippocampus
surface data used in this investigation was brought into rigid alignment by ICP
in a preprocessing step, then 2D sliced and manually segmented by a trained
neuro-anatomist into a 40x21 grid periodic in one direction. This grid in effect
forms a parametric domain for the corresponding surface, and the Riemannian
surface characteristics are thus calculated by numerical methods and stored in
said domain.

3 Results

We here present examples of the algorithm in action, validating its behavior on
synthetic cases and demonstrating results obtained on real hippocampus pairs.
As a starting point, we consider the trivial case of a pair of cylinders, where
one member of the pair has had its surface distorted according to an outward
normal vector field of magnitude dictated by a Gaussian distribution. That is
to say, its surface has been bumped. In Fig. 2, we see a comparative case in-
volving a cylinder with two such distortions versus its unaltered counterpart.
The intuition behind our method asserts that we should "blame" the applied
bumps for observed differences between the Riemannian surface characteristics
of the two shapes. Obviously, then, the desired result is the segmentation of the
two clearly visible bumps on the left cylinder, which is what we observe. (Red

'hot spots,' the curve interiors, represent regions of highest deformation energy
under the current map and segmentation. The curve itself is found at the border
between the red and blue regions.) Note the textured quality of the surfaces,
which occurs as a result of the addition of a SD=0.02 zero-mean Gaussian noise
field to the (X,Y,Z) coordinates of the (unit radius) cylinders prior to output,
for demonstration of robustness in the face of noise in surface scanning.
As indicated previously, the case of the cylinder is a trivial one in this context
in that there exists a "natural" parameterization of the surface which results in
a constant metric tensor field over the parametric domain. Thus, we further val-
idate using synthetic alterations to real data exhibiting arbitrary shape, within
the confines of cylindrical topology. Since the hippocampus meets these criteria
and forms part of our clinical motivation, it seems a natural choice. Selecting
a single hippocampal horn, we again apply a small outward normal Gaussian
deformation to its surface, whose location we hope to see the segmenter recover.
However, rather than illustrating invariance to surface point cloud noise, we now
seek to additionally validate invariance of the solution to reparameterizations of
the given surfaces. To pause for conceptual refreshment, recall that in compar-
ing surfaces we are in fact comparing two parametric domains based on the two
tensor-valued functions corresponding to each of those domains, those functions
being the Riemannian surface structure information of each parameterized sur-
face. But just as two surfaces (or their parametric domains, as we do in our
implementation) can be mapped to one another nonrigidly in infinitely many
ways, so too can each parametric domain map to its topologically-equivalent pa-
rameterized surface in infinitely many ways. The Riemannian structure tensor
fields themselves depend upon which of these infinitely many parameterizations
are chosen. When synthetically producing one surface from another, we pos-
sess a trivial alignment between the domains under which the two tensor fields
are identical, excepting the patch corresponding to the synthetic bump. In real
cases, of course, data will not have this sort of property. The validation is thus
incomplete until we can show that under a suitably deformed reparameteriza-
tion of one of the surfaces being compared, a result comparable to the above is
achieved. We can reparameterize simply by randomly generating a distortion of
the parametric domain's regular grid, then interpolating a new point cloud from
the given data based on the new locations of the gridpoints. It is completely
appropriate to confine said distortion to the order of a pixel in each direction:
after all, we have assumed close rigid alignment and sensibly spaced data points
from our preprocessing, and if these constraints are violated, they can always
be reinstated by repeating the algorithms that enforced them in the first place.
The purpose of the nonrigid registration evolution is to further refine the surface
correspondence following bulk rigid alignment of the shapes based on the sur-
face characteristics observed locally, not to register raw point clouds. We cannot
show the parameterization warp used due to space constraints: it is of 0.1 SD in
each coordinate.
It is visually evident that the segmentation succeeds despite the change of
parameterization of one of the surfaces being compared. The deformation map

II ,.

Fig. 2. Left column: Segmentation of prominent distortions (bumps) from cylindrical
surface, based on comparison to homologous cylinder lacking bumps. Surfaces are both
distorted with distinct random noise fields (seen as textured effect). Top: naive initial-
ization; Middle: 3 iterations; Bottom: 10 iterations. Right column: Segmentation of
synthetic bump from surface of real hippocampus, based on comparison to same hip-
pocampus without bump. The two hippocampi are under different parameterizations.
Top: naive initialization; Middle: 5 iterations; Bottom: 20 iterations.

between the parametric domains (not shown) undergoes considerable evolution
through this process. This evolution can be best understood as the method's
attempt to invert the reparameterization, to achieve the optimal alignment of the
identical surface patches external to the segmented bump, while at the same time
trying to find the minimum energy description of the acknowledged deformation
of the patch on the curve's interior.

In presenting the following real example, we issue the caveat that the method's
success lies in the eye of the beholder, because there is of course no ground truth
associated with this sort of problem. For this reason we have selected a case
with an ".i- ,-. i which can easily be agreed upon. Fig. 3 shows the left and
right hippocampi clinically segmented from an MRI scan of an epilepsy patient.
We take it as self-evident that there exists a localizable region of high distor-
tion between the left and right shapes, and assume that one will agree that our
method has in fact converged to a segmentation of that region. Illustrations of
the segmented energy field and parametric warp evolutions have been withheld
due to space limitations.

Fig. 3. Segmentation of distortion between real hippocampi of epileptic. Top left: naive
initialization; Top right: 1 iteration; Bottom left: 5 iterations; Bottom right: 20 itera-

4 Conclusion

We have presented a novel scheme for simultaneous nonrigid registration and seg-
mentation of homologous shapes, wherein the interdependent registration and
segmentation processes are driven by intrinsic geometric characteristics of the
shapes themselves. As such, we are able not only to identify surface pairs rep-
resenting large deformations, but also to specify which subregions of the shapes
appear most likely to be involved in the deformation, to produce an estimate
of the implicit deformation field, and to quantify the deformation energy of the
segmented subregions separately from that of the remaining surface patches.
Given the number of considerations interplaying, the variational principle driv-
ing the evolution is notably compact. We have seen that the method's results
are reasonably invariant to local parameterization changes and added surface
noise, and results on real data (while inherently lacking in ground truth for
comparison) appear promising. Future efforts could include the application of
the method to extensive real datasets that have been clinically preclassified ac-
cording to the presence or absence of pathology, as well as an attempt to correlate
the locations and intensities of the segmented deformed regions with the given
pathology classifications. It should also be noted that the two-region Chan-Vese
segmentation method (which is often used for its strong balance of effective-
ness and relative simplicity of implementation) is merely a good "first draft" of
a segmentation scheme for this sort of problem. A more robust segmentation
method, such as the general Mumford-Shah functional extended to arbitrary
surfaces, would doubtless prove beneficial, as our energy fields will not always
neatly obey the two-mean framework. Finally, it goes without saying that all
evolution/optimization-based methods (apart from the most trivial cases) suffer
from some amount of initialization-dependence. Intelligent initialization schemes
could be expected to further boost performance.


1. Gerig, G., Muller, K.E., Kistner, E.O., Chi, Y.Y., Chakos, M., Styner, M., Lieber-
man, J.A.: Age and treatment related local hippocampal changes in schizophrenia
explained by a novel shape analysis method. In: MICCAI. (2003) 653-660
2. Styner, M., Lieberman, J.A., Gerig, G.: Boundary and medial shape analysis of
the hippocampus in schizophrenia. In: MICCAI. (2003) 464-471
3. Wang, Y., Gu, X., Hayashi, K.M., Chan, T.F., Thompson, P.M., Yau, S.T.: Surface
parameterization using riemann surface structure. In: Proc. Int. Conf. on Computer
Vision. (2005) 1061-1066
4. Wang, Y., ('C I.... M.C., Thompson, P.M.: Mutual information-based 3d surface
matching with applications to face recognition and brain mapping. In: Proc. Int.
Conf. on Computer Vision. (2005) 527-534
5. Wang, Y., C('I .... M.C., Thompson, P.M.: Automated surface matching using
mutual information applied to riemann surface structures. In: MICCAI. (2005)
6. A.Yezzi, Zollei, L., Kapur, T.: A variational framework for joint segmentation and
registration. In: CVPR MMBIA. (2001)

7. Unal, G.B., Slabaugh, G.G.: Coupled pdes for non-rigid registration and segmen-
tation. In: CVPR. (2005) 168-175
8. Wang, F., Vemuri, B.C.: Simultaneous registration and segmentation of anatomical
structures from brain mri. In: MICCAI. (2005) 17-25
9. Mumford, D., Shah, J.: Optimal approximations by piecewise smooth functions
and associated variational problems. Commun. Pure Appl. Math 42(4) (1989)
10. Chan, T.F., Vese, L.A.: Active contours without edges. IEEE Transactions on
Image Processing 10(2) (2001) 266-277
11. Tsai, A., Yezzi, A.J., Willsky, A.S.: Curve evolution implementation of the
mumford-shah functional for image segmentation, denoising, interpolation, and
magnification. IEEE Transactions on Image Processing 10(8) (2001) 1169-1186

University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - - mvs