Corresponding Supine and Prone Colon Visualization Using Eigenfunction Analysis and Fold Modeling
Abstract
We present a method for registration and visualization of corresponding supine and prone virtual colonoscopy scans based on eigenfunction analysis and fold modeling. In virtual colonoscopy, CT scans are acquired with the patient in two positions, and their registration is desirable so that physicians can corroborate findings between scans. Our algorithm performs this registration efficiently through the use of Fiedler vector representation (the second eigenfunction of the LaplaceBeltrami operator). This representation is employed to first perform global registration of the two colon positions. The registration is then locally refined using the haustral folds, which are automatically segmented using the 3D level sets of the Fiedler vector. The use of Fiedler vectors and the segmented folds presents a precise way of visualizing corresponding regions across datasets and visual modalities. We present multiple methods of visualizing the results, including 2D flattened rendering and the corresponding 3D endoluminal views. The precise fold modeling is used to automatically find a suitable cut for the 2D flattening, which provides a less distorted visualization. Our approach is robust, and we demonstrate its efficiency and efficacy by showing matched views on both the 2D flattened colons and in the 3D endoluminal view. We analytically evaluate the results by measuring the distance between features on the registered colons, and we also assess our fold segmentation against 20 manually labeled datasets. We have compared our results analytically to previous methods, and have found our method to achieve superior results. We also prove the hot spots conjecture for modeling cylindrical topology using Fiedler vector representation, which allows our approach to be used for general cylindrical geometry modeling and feature extraction.
0 \vgtccategoryResearch \authorfooter Saad Nadeem, Joseph Marino, Xianfeng Gu, and Arie Kaufman are with the Computer Science Department, Stony Brook University, Stony Brook, NY 117942424. Email: {sanadeem,jmarino,gu,ari}@cs.stonybrook.edu \CCScatlist \CCScatK.6.1Management of Computing and Information SystemsProject and People ManagementLife Cycle; \CCScatK.7.mThe Computing ProfessionMiscellaneousEthics \teaser (a) (b) \vgtcinsertpkg
Introduction
Colorectal cancer (CRC) is the third most frequently diagnosed cancer worldwide, and the fourth leading cause of cancerrelated mortality with 700,000 deaths per year worldwide. Optical colonoscopy (OC) is an uncomfortable invasive technique commonly used to screen for CRC. In contrast, virtual colonoscopy (VC) was introduced as a noninvasive procedure for the mass screening of polyps [12], the precursors to CRC. VC reconstructs a 3D colon model from a computed tomography (CT) scan of a patient’s abdomen and then a physician virtually navigates through the colon looking for polyps. The low radiation dosage, the recent assignment of CPT code 74263 for insurance reimbursement, the lack of side effects from anesthesia, and the quick turnaround time for results makes VC a prime method for CRC screening.
In a typical VC, the patient is scanned twice, supine (facing up) and prone (facing down) or, more recently, in a side position called the lateral recumbent position to avoid artifacts due to the squishing of the colon in the prone position. The two scans improve the sensitivity and specificity of polyp detection and reduce the extent of uninterpretable collapsed or fluidfilled segments. Because the colon changes between the supine and prone (or side) scans, the two colon models need to be registered to match detections. Computeraided detection (CAD) techniques [13, 36] have been introduced to automatically detect polyps, but they have issues. Registration can aid these CAD techniques by reducing false positives (primarily due to stool and colon folds) when matched across scans. Moreover, the recurrence of cancer due to the incomplete removal of polyps can be tracked via registration of scans captured at different periods. This allows physicians to perform an effective followup examination by visualizing the polyp growth over time. Hence, registration provides an effective tool to distinguish between lesions, aid in CAD, and visualize polyp growth over time. In this paper, we focus on registering VC data across multiple patient positions.
The first nontrivial eigenfunction of the LaplaceBeltrami operator is known as the Fiedler vector representation. This representation can be used by itself for complex shape modeling of tubular structures such as the colon. Lai et al. [16] proposed the use of the Fiedler vector for supine and prone colon registration. However, their registration depends on the accurate detection of important landmarks, based on assumptions which often did not hold in our patient colon datasets, leading to a higher registration distance error. Instead, we integrate the Fiedler vector representation with a more accurate feature detection based on haustral fold modeling to achieve better registration and visualization results. We first use the Fiedler vector representation for global registration of the scans. We then use 3D level sets computed based on this representation to segment the folds. Finally, we refine the registration locally via these folds as anatomical references. An overview of our pipeline is given in Figure 1. The accurate modeling of the colon folds and the Fiedler vector representation results in automatic cross visualization of consistent endoluminal views across multiple orientations.
Different flattening approaches have been introduced to visualize the complex geometry of medical organs [3, 10, 11, 13, 21, 20, 27]. The Fiedler vector representation and fold segmentation also allows us to create more accurate 2D flattening visualizations by automatically computing consistent cuts throughout the colon without crossing the folds. These cuts are then used as an input to a quasiconformal mapping [35], which is used to flatten the colon. Previous flattening approaches require manually selecting extrema on the colon and computing geodesic paths to create a cut. This cut, however, crosses folds on sharp bends since this is considered the shortest path, resulting in flattened visualizations with chopped folds, as shown in Figure 7.
There have been other registration approaches which require either using the centerline or flattening approaches which distort texture and require high quality data. Moreover, once the data is flattened to a 2D rectangle, it is difficult to find correspondences between the original data and the flattened data due to the loss of geometry. It is our endeavor here to be able to easily visualize this correspondence and correlate anatomical features across scans in the original and the flattened data. Our registration framework makes use of accurate fold models to segment out folds with high precision. These fold segmentations are then used as anatomical references to register the VC data across scans. Moreover, the fold segmentations allow us to automatically find suitable cuts to map the 3D data to simpler 2D topologies and visualize the resultant data more effectively. Our contributions are as follows:

Integrated colon registration framework using Fiedler vector representation for global registration and 3D level set fold segmentation for local refinements.

Accurate fold modeling which provides stateoftheart detection and segmentation results on colon models.

Accurate correspondence visualization between 3D and 2D flattening visualizations due to the Fiedler vector representation and accurate fold segmentation.

Extraction of a consistent cut along the colon for visually accurate colon flattening.

The theoretical proof of the hot spots conjecture for modeling general cylindrical topology using Fiedler vector representation.
The paper is organized as follows. In the next section, we discuss related work, followed by a brief theoretical background of Fiedler vector representation in our context in Section 2. This provides the basis for our colon registration algorithm, which we outline in Section 3. In Section 4, we highlight the supine and prone correspondence and more effective flattening visualizations based on our algorithms. Finally, we evaluate our algorithm and show additional results in Section 5, followed by the conclusion and future work in Section 6.
1 Related Work
There have been several works in the colon registration domain based on the identification of landmarks and centerline matching across supine/prone/side scans. A basic method applies linear stretching and shrinking operations to the centerline, where local extrema are matched and used to drive the deformations [1, 2]. Correlating individual points along the centerline through the use of dynamic programming has also been suggested [7, 17, 24].
More recently, the taeniae coli (three bands of smooth muscle along the colon surface) have been used as features which can be correlated between the two scans [14]. This relies on a manual identification of one of the three taeniae coli, and then an automatic algorithm repeats the line representing the identified taenia coli at equal distances. Further progress has been made where the haustral folds and the points between them can be automatically detected, and the taeniae coli are identified by connecting these points [30]. However, this method is only feasible on the ascending and transverse portions of the colon.
Deformation fields have also been suggested for use in supineprone registration. Motion vectors can be identified for matched centerline regions, interpolated for nonmatched regions, and then propagated to the entire volume [29]. It has also been proposed to use a freeform deformation grid to model the possible changes in the colon shape from supine to prone [26].
Conformal mapping has been successfully used for many medical applications, including a brain cortex surface morphology study [11] and colonic polyp detection [13]. Quasiconformal mapping was used to map supine and prone colon datasets to a rectangular plane and convert the 3D registration problem into a 2D image matching problem [35]. This, however, required landmarks to be extracted in order to divide the colon into different welldefined segments. Moreover, the whole registration pipeline was dependent on flattening these segments by tenaie coli extraction and then using graph matching. The flattening was computationally expensive, and precise fold segmentation could have improved the accuracy of the graph matching algorithm.
(a) Supine  (b) Side 
There has been minimal use of folds as anatomical references for registration due to the imprecise modeling of folds, thus leading to poor fold detection and segmentation results. Effective segmentation of the colon folds is necessary for detecting polyps on the folds [34], extracting teniae coli [14, 38], and performing supineprone registration [2, 1, 31, 35]. Most of the previous works have focused on colon fold detection [5, 14, 33], disregarding the precise segmentation or delineation of the boundaries of the fold. A fold boundary modeling and segmentation algorithm was introduced, and the segmentedarea ratio () metric was also introduced to measure the accuracy of the segmentation [37]. In this work, we model and segment the folds (using Fiedler vector representation) on the supineprone or supineside VC dataset pairs which in turn are used as anatomical references to register these pairs.
Fiedler vector representation has been used in various applications including mesh processing [8], mesh parameterization [23], and shape characterization [4]. Lai et al. [16] have proposed the use of Fiedler vector representation for supine and prone colon registration. Based on the Fiedler vector computation for supine and prone colons, they detect the hepatic and splenic flexures, and register the supine and prone colons piecewise using these detected landmarks. They assume that the splenic and hepatic flexures are denoted by local maximum zcoordinates nearest to the maximum and minimum Fiedler vector values, respectively. This assumption, however, is frequently invalid in real patient colon datasets which we have tested on and can result in higher distance registration error compared to the registration done solely based on the Fiedler vector without considering the landmarks. Moreover, their algorithm does not take into account the Fiedler vector flips that can occur between the corresponding supine and prone colons which could lead to erroneous registration results (e.g. the rectumsplenic segment being registered to the cecumhepatic segment).
(a) Fiedler vector on colon segment 
(b) Fiedler vector level sets on colon segment 
(c) Segmented folds on colon segment 
(d) Fiedler vector on flattened colon segment 
(e) Fiedler vector level sets and cross curves on flattened segment 
(f) Segmented folds on flattened colon segment 
In contrast, we use the Fiedler vector representation to model the human colon, register the supine/prone/side models in 3D, and segment folds using level sets based on this representation. The Fiedler vector representation is used to register the supine/prone/side datasets globally and then the segmented folds are used as anatomical references to locally refine the registration results. The fold segmentation based on the level sets is efficient and the overall registration process does not depend on the computation of the centerline, the teniae coli, or the flexures. Our use of the Fiedler vector representation for fold segmentation is closely related to spectral image [19] and shape segmentation [28]. The colon registration in 3D, in essence, is similar to shape matching using functional maps [15, 25].
2 Theoretical Background
In this section, we briefly introduce the theoretical background necessary for this work, which allows us to compute extrema and model the complex geometry of the colon. The hot spots conjecture [6] states that the minimum and maximum of the Fiedler vector are the two points with the greatest geodesic distance on the surface. We prove in Lemma 2 that this conjecture holds for a topological cylinder.
(a)  (b)  (c) 
Suppose is a surface embedded in the threedimensional Euclidean space , and the induced Euclidean Riemannian metric is denoted as . One can choose a special type of local coordinates , called isothermal coordinates, such that the metric has a canonical form , where the function is called the conformal factor. The Laplace Beltrami operator induced by the metric is given by . The Gaussian curvature of the surface is given by . An eigenfunction of the LaplaceBeltrami operator is given by , where the eigenvalue is a nonnegative real number. The LaplaceBeltrami operator has infinitely many eigenvalues and eigenfunctions. The eigenvalues are all sorted in ascending order . The eigenfunctions form an orthonormal basis of the functional space of the surface, . The first eigenfunction is a constant function. The second eigenfunction is called Fiedler’s function.
The heat kernel of the surface is . The heat diffusion equation on the surface is
(1) 
and the solution to the heat equation is given by the heat kernel . Another equivalent way to represent the solution to the heat equation is
(2) 
where
When goes to infinity, the right hand side of Equation 1 goes to , which means that the temperature function becomes harmonic.
Lemma 1
Suppose is a harmonic function, then for any interior point ,
(3) 
where is a small circle surrounding .
Proof 1
If we choose isothermal coordinates , then is a harmonic function and . We construct the conjugate harmonic function , such that and , then the complex function is a holomorphic function, where . Note that is constructed using as a real part and as an imaginary part. By Cauchy’s formula, we obtain
where , . Comparing the real and imaginary parts, we obtain Equation 3.
Corollary 1
A metric surface has boundaries . Suppose is a harmonic function on , then the maximal and minimal points of are on .
Proof 2
Equation 3 means that the value of an interior point is the mean of the values of its neighboring points on the small circle . Therefore, the maximal and minimal points cannot be in the interior of and hence, these points must be on the boundaries .
Lemma 2
Suppose the metric surface is a topological cylinder, then the maximal and minimal points of its Fiedler vector are on the boundaries of the surface .
Proof 3
Suppose the solution to Equation 1 is . When goes to infinity, converges to a harmonic function, and therefore the maximal and minimal points of are close to the boundaries .
From Equation 2, when becomes large enough, the high order eigenfunctions go to much faster. Therefore, the behavior of is mainly controlled by the first two terms,
where is constant. Therefore, the maximal and minimal points of are approximated by those of , and thus the maximal and minimal points of approach the surface boundaries .
Level Sets of Fiedler’s Function
Suppose is a diffeomorphic mapping between two metric surfaces. The local coordinates of are and those of are . The Jacobian matrix of the mapping is . The pull back metric induced by the mapping is defined by . The mapping is conformal if there is a function ,
(4) 
According to conformal geometry, if is a topological cylinder, then there is a conformal mapping which maps the surface onto a flat cylinder , . The second eigenfunction of is . The level sets of the second eigenfunction are circles with constant height .
Suppose a topological cylindrical surface is conformally mapped onto a flat cylinder . Furthermore, the mapping is nearisometric, namely, the conformal factor in Equation 4 is close to . Then the level sets of the second eigenfunction of are similar to those on . In our current work, we use the level sets of the second eigenfunction to locate the folds on the colon surface.
3 Computational Algorithm
(a)  (b)  (c)  (d)  (e) 
(f) Flattened colon Fiedler vector representation  
(g) Flattened fold segmentation 
In this section, we explain the computational algorithms for our registration method. We globally register the supine and prone colons based on their respective Fiedler vector representations. This is followed by bounded fold contour detection using 3D level set bundles computed using the Fiedler vector representation, and the eventual segmentation into individual folds based on the 2D normalization of the projected detected fold contour planes. These segmented folds are used as references on the globally registered colons to locally refine the registration and result in the final output of our registration algorithm.
3.1 Fiedler Vector
Suppose the colon surface is , where is the Riemannian metric. The colon surface is a topological cylinder with two boundaries, . We want to compute the Fiedler vector with the Dirichlet boundary condition
where is the LaplaceBeltrami operator [22] induced by .
In the discrete setting, the colon surface is extracted from CT images, and approximated by a triangular mesh, denoted as , where is the set of vertices, the set of edges, and the set of faces. The cotangent edge weight is defined as follows: suppose is an edge, shared by two triangular faces and , then the edge weight is , where is the corner angle at vertex in the face . Suppose a function is defined on the vertex, . The discrete LaplaceBeltrami operator is defined as , where means the vertex is adjacent to . Therefore, the LaplaceBeltrami operator has a matrix representation , where
We compute the eigen decomposition of ; the eigenvalues are and the eigenvectors are . The first eigenvalue of is and the first eigenvector is . The first positive eigenvalue is and the Fiedler vector is the corresponding eigenvector . We scale such that the maximum value of is and the minimum value is . In the following discussion, we treat as a function defined on the vertex set. Figure 3(a) illustrates the Fiedler vector of a colon surface with color encoding, and Figure 3(d) shows the Fiedler vector on the flattened colon surface.
3.2 Global Registration
Two colon surfaces can be registered directly based on their Fielder vectors. Suppose the two colon surfaces are and , then the corresponding Fiedler vectors are and , respectively. Furthermore, the boundaries of the two surfaces are . We first match the boundary curves using their arc lengths. We choose base points and , where the base points are the extrema which are computed using the maximum and minimum Fiedler vector values. To avoid Fiedler vector flips between the corresponding supine and prone data, we manually make these maximum or minimum Fiedler vector values consistent on the rectum or the cecum for both the supine and prone datasets. The arc length is then used to parameterize the boundary curves and normalize the total length to . The arc length parameterizations are denoted as and . The mapping gives the mapping from to .
We compute the gradient field on the surfaces, denoted as and respectively, then we trace the integration curves of the gradient fields. The curve represents an integration curve on , and . Similarly, is an integration curve on starting from . A point is the intersection of a level set of and an integration curve starting from a point on with the arc length parameter . The whole colon surface is parameterized by , and is globally parameterized by . The initial global registration is given by the mapping .
3.3 Fold Detection
Suppose is a colon surface with Fiedler vector . We denote the level sets of as , where , namely . We uniformly sample the level sets on the colon surface and obtain a family of level sets , as shown on the colon surface in Figure 3(b) and the flattened colon surface in Figure 3(e). We also compute the integration curves of the gradient field of and obtain a family of integration curves , , as shown in Figure 4, where the red loops are the level sets and the yellow curves are the integration curves .
We compute the normal curvature of points . If a point is in a fold area, then the normal curvature of is negative; in other flat or convex areas, the normal curvature of is zero or positive. We compute the inflection points of the integration curves, where the normal curvatures are ’s, and the points with minimal normal curvatures. Fixing an integration curvature , suppose , where and are inflection points, is the minimal curvature point, then the level sets , form a level set bundle. In this way, we compute the clusters of level sets. We discard the initial uniformly sampled level sets, which do not belong to any bundles. The level set bundles are shown in Figure 5(d).
We densely sample the level sets within each bundle, and compute the normal curvatures. The points with negative normal curvature of curve are in the fold area. In this way, we can detect the fold contour, as shown in Figure 5(e).
Bowel preparation, done prior to VC, might lead to local underdistention of some colon regions, which we refer to as collapsed regions. In this case, we ignore the extraction of folds in these regions since these folds do not exhibit any meaningful information. We set a threshold on the level set bundle such that if the encompassing circle is below 0.1 in radius, when normalized in 2D, we mark that region as locally underdistended or collapsed, as shown in Figure 8. However, severe cases of local underdistension can lead to a complete collapse of the colon at a particular region, resulting in multiple colon segments during the segmentation process. Our method does not cater to these severe cases, which will be a focus of future research.
3.4 Haustral Fold Segmentation
For each level set , we find a best fit plane as follows. We choose samples on , denoted as , and the center of the samples is given by . We compute the covariance matrix of the samples, , then compute the eigen decomposition of . The eigenvectors are , which form an orthonormal, then the fitting plane goes through the center and is spanned by and , namely, is given by the equation . The fitting planes of the level sets are illustrated in Figure 5(c).
We then project the level set onto the fitting plane . Suppose , then the projection on the plane is given by . We denote the projected level set as . Suppose and are two fitting planes, with centers and , we align them together as follows: we shift to by a translation, then rotate to coincide with by a rotation , such that the rotation angle is the angle between the normals to the two planes and the rotation axis is along the cross product of the two normals. In this way, we align the fitting planes of all level sets within a bundle, and the spatial level set curves become a cluster of planar curves as shown in Figure 4(c). We denote each projected and aligned level set as .
An integration curve on the original surface becomes a planar curve by mapping to . Suppose the level set bundle is , where . We compute the total length of each and find the local minima with respect to . Then each local maximum corresponds to a haustral fold. This procedure produces the haustral fold segmentation, as shown in Figure 5(e) for the 3D view and Figure 5(g) for the flattened view. We compare our automatic haustral fold segmentation with manual segmentation as shown in Figures 6 and 9, where the green contour shows the manual segmentation results and the red contours show the automatic results. One can see that the automatic segmentation has very high accuracy.
3.5 Local Registration Refinement
The folding areas and the segmented haustral folds are used as anatomical references. From the initial registration, we can find the correspondences among the haustral folds on colon surfaces obtained from the supine/prone/side data. We then locally deform the initial mapping in order to align the haustral folds more accurately. We denote the two flattened colon surfaces as and respectively, is the mapping between them. We define the characteristic function , , where equals if is inside a haustral fold, and otherwise. Then we smooth the characteristic functions out by convolving a Gaussian filter. We define an energy for the mapping as follows:
The first term ensures that haustral folds match haustral folds, and the second term measures the smoothness of the mapping. We obtain the EulerLagrange equation, and the following flow minimizes the energy,
By deforming the mapping along the flow, the registration is improved significantly. Figure 10 shows the registration result using the haustral folds as references.
(a) Fiedler vector on colon segment 
(b) Flattened colon segment with a cut via geodesic path 
(c) Flattened colon segment with a cut via our fold segmentation 
4 Visualization
Due to the use of the Fiedler vector representation throughout the entire colon, it becomes possible to easily colocate positions between two scans, either manually or automatically. Using consistent cuts along the colon based on the segmented folds, we achieve more effective flattened visualizations than other works. The Fiedler vector representation also provides easier visualizations to find correspondences between the 3D model and the flattened view while allowing for better polyp localization via labeling of segmented folds.
(a) Supine  (b) Prone 
(c) Supine  (d) Prone 
4.1 Corresponding Supine and Prone Visualization
The Fiedler vector representation can be used to automatically visualize corresponding regions in multiple patient orientations and in respective 3D and 2D mapped views. The consistent endoluminal views in Figures Corresponding Supine and Prone Colon Visualization Using Eigenfunction Analysis and Fold Modeling and 2 highlight the polyps in supineprone and supineside orientations, respectively. In order to visualize both supine and prone in a consistent manner, we map the camera frustum viewing the supine mucosa (inner surface) to the camera frustum viewing the prone mucosa. This is done by computing the means of the level sets based on the Fiedler vector which result in corresponding centerlines on both supine and prone. We then match the camera orthonormal coordinate frames for the two views.
The centerlines are matched based on the eigenfunction bundles following the registration algorithm, outlined in Section 3. Given a point on the discretized supine centerline, we want to find the corresponding point on the discretized prone centerline. We first find the point on the supine centerline, which is the point nearest in . Given this point , the corresponding point on the prone centerline is found at the same index. can then be found as the point on the prone centerline closest in to .
We focus our endoluminal correlation work on situations where the user is looking at something on the colon wall and wants to view the corresponding location in the other scan. Using the correspondence along the skeletons, it is possible to create a correlated automatic navigation, though we found running two flythroughs sidebyside to be distracting and more of a hindrance than help. Orienting rotation around the centerline for two automatic navigation views is possible based on the individual haustral fold registration on the corresponding supine and prone colons.
4.2 Flattened Visualizations
Since the Fiedler vector representation captures extreme points on a given colon dataset, we can use these extreme points along with the fold segmentation results to trace a consistent cut throughout the colon. This cut can then be used as input to a quasiconformal mapping algorithm [35] to flatten the colon to a 2D plane. In essence, the following steps are used to automatically extract a consistent cut: (1) extract the extreme points on the colon dataset, (2) remove the 1ring vertex set corresponding to the extreme vertices, (3) pick two haustral folds (in the same level set bundle) closest to one of the extreme points, (4) connect the mean of the endpoints of these folds to the next level set bundle of haustral fold endpoints, closest to the centroids of the folds in the previous bundle, and (5) repeat the previous step until the other extreme point is reached. The result of this process, compared to a typical cut based on a geodesic path, is illustrated in Figure 7.
To the best of our knowledge, there are no automatic algorithms for extracting a consistent cut throughout the colon. Normally, two extreme points are manually picked on a given colon dataset and a geodesic path is computed from one point to the other and this constitutes the cut. However, the problem with this cut is that it crosses the folds when a sharp bend is incurred, as shown in Figure 7, whereas our method takes into account individual haustral folds before tracing the cut which allows for consistency.
In addition to greater consistency, our cut produces less angle distortion as compared to the geodesic path cut. We quantify the angle distortion metric of the resultant flattened colon by using the signed singular values of the Jacobian of the transformation for each triangle [18]. Small angular distortion is indicated by a distortion value approaching 2. The colon flattened via geodesic path cut in Figure 7(b) produces an overall angle distortion value of 2.181, whereas our method in Figure 7(c) produces a value of 2.126.
5 Experimental Results
We have validated our algorithms using real VC colon data from the publicly available National Institute of Biomedical Imaging and Bioengineering (NIBIB) Image and Clinical Data Repository provided by the National Institute of Health (NIH). We perform electronic colon cleansing incorporating the partial volume effect [32], segmentation with topological simplification [13], and reconstruction of the colon surface as a triangular mesh via surface nets [9] on the original CT images in a preprocessing step. Though the size and resolution of each CT volume varies between clinical datasets, the general data size is approximately voxels with a resolution of approximately mm.
We have developed our algorithms using generic C++ on the Windows 7 platform. The linear systems for the Laplace equation were solved using the Matlab C++ library. All of the experiments are conducted on a workstation with a Core 2 Quad 2.50GHz CPU with 4GB RAM. On average, the Fiedler vector computation takes 3.5 seconds, the haustral fold detection and segmentation takes 3.2 seconds, and the final local refinement takes 1.8 seconds.
Scan ID  # True  Sensitivity  #FNs  #FPs  
1Supine  285  0.94  18  0  0.82 
1Side  272  0.94  16  1  0.92 
2Supine  191  0.86  27  0  0.87 
2Side  169  0.83  28  0  0.86 
3Supine  293  0.99  2  0  0.93 
3Prone  268  0.97  7  1  0.81 
4Supine  274  0.93  18  0  0.91 
4Prone  251  0.90  24  0  0.83 
5Supine  165  0.84  26  1  0.86 
5Prone  152  0.82  27  0  0.76 
6Supine  223  0.94  13  0  0.86 
6Prone  202  0.88  25  1  0.76 
7Supine  196  0.88  24  1  0.91 
7Prone  161  0.94  9  0  0.78 
8Supine  247  0.97  8  0  0.94 
8Prone  210  0.91  19  0  0.86 
9Supine  187  0.82  33  1  0.84 
9Prone  158  0.91  15  0  0.75 
10Supine  173  0.87  22  0  0.82 
10Prone  142  0.92  11  0  0.89 
Average  211.0  0.90  18.6  0.30  0.85 
49.4  0.05  8.36  0.47  0.06 
5.1 Fold Detection and Segmentation Evaluation
To evaluate our fold segmentation algorithm, we manually segmented the folds on 10 supineprone and supineside colon pair datasets. These manual segmentations were approved by a VCtrained radiologist who reviewed them. The accuracy of the proposed segmentation algorithm was measured on a perfold basis, using the metric [37]. Given the manuallyestablished ground truth, a fold was assumed to be successfully detected if more than 50% of its area has been detected. Hence, we can calculate the true positives (TP), false positives (FP), true negatives (TN), false negatives (FN), sensitivity (), and specificity (). We calculate the as follows:
where is the ground truth area of the fold. is the area of the segmented fold, and is the area of the overlap of the above two areas. is defined as the ratio between the area of the intersection and the area of the union of the expertdrawn folds and the automatically segmented folds. In essence, a larger suggests a better segmentation result, and indicates that the segmented fold perfectly matches the ground truth without any over or undersegmentation.
Table 1 lists the quantitative results from the 20 patient scans. About 90.4% of all the ground truth folds were detected with approximately 18 missed per dataset, under the assumption that a fold would be detected if more than 50% of its area has been segmented. The missed folds were mostly shallow and flat shaped structures. The FP detections were very few for each scan. For all the detected folds, the is approximately 84.8%, indicating that the automatically segmented fold boundaries have matched fairly well with those of the manuallydrawn counterparts. A visual comparison of the detection accuracy of our algorithm with a previous fold detection method using heat diffusion and fuzzy Cmeans clustering is provided as part of the Appendix A.2.
(a) Unregistered supine and prone colon section 1 
(b) Registered supine and prone section 1 
(c) Unregistered supine and prone section 2 
(d) Registered supine and prone section 2 
5.2 Analytic Registration Evaluation
To evaluate our registration results, we use the distance measurement between corresponding features located on registered colon surfaces. We use the haustral fold segmentation algorithm, detailed in Section 3, to find corresponding anatomical references on a given supineprone colon pair. Based on these computed references and our manual labels, we use a subset of the segmented features for our registration and measure the distance errors on the remaining features. We generally extract more than 100 features from both pairs. For this registration error evaluation, we used 75 feature points for our registration and measure the distance errors on the remaining features.
A comparison between our method and other methods is performed using our analytic evaluation results in . For those papers that present their distance error, we compare our results with their results in Table 2. Our method produces a registration with significantly smaller distance error between corresponding points than other methods. When comparing to the published algorithm of Lai et al. [16], 3 out of our 10 colon pair datasets had flips and needed to be made consistent manually in order for their algorithm to work. Moreover, in 6 out of our 10 colon pair datasets, the Lai et al. algorithm performed worse than our global registration baseline ( in their case) due to the incorrect detection of the landmarks. Further details about issues we encountered with their algorithm can be found in Appendix A.1.
Method  Dist. Error 

Fiedler Vector Representation Fold matching  5.24 mm 
Quasiconformal mapping [35]  7.85 mm 
Fiedler Vector Representation  11.98 mm 
Centerline registration + statistical analysis [17]  12.66 mm 
Linear stretching/shrinking of centerline [2]  13.20 mm 
Centerline feature match + lumen deformation [29]  13.77 mm 
Fiedler vector + piecewise registration [16]  14.19 mm 
Centerline point correlation [7]  20.00 mm 
Taeniae coli correlation [14]  23.33 mm 
5.3 Visual Registration Verification
We provide a visual registration verification of our algorithm by flattening the entire supine and prone colon models and then using the segmented folds (see Figure 8) as anatomical references to align the supine and prone datasets. In Figure 10, we show the entirety of both the supine and prone colon models, mapped to the plane, both unregistered and registered. The images of the registered segments clearly show very good alignment of the supine and prone colon structures, whereas the unregistered segments show poor alignment. Compared to previous methods, our algorithm provides a visually superior result for viewing corresponding locations on the two colon models.
We have also shown our results to a VCtrained radiologist. Due to the use of the Fiedler vector color mapping on the colon, he was able to easily find the corresponding regions across supine and prone colons along with a strong correlation between the 3D endoluminal views and the flattened views.
6 Conclusions and Future Work
Shape registration is fundamental for shape analysis problems, especially for abnormality detection in medical applications. We have introduced an efficient framework for the registration of supine and prone colons, through the use of Fiedler vector representation, to improve the accuracy of polyp detection. Specifically, we use the Fiedler vector representation to globally register the supine and prone colon models. We then use level sets computed based on the Fiedler vector representation to detect and segment folds which are in turn used as anatomical references to locally refine the global registration results. The use of the Fiedler vector representation can help in easily visualizing the corresponding regions on 3D and 2D mappings as well as across supine and prone models. We have also proven the hot spots conjecture for modeling cylindrical topology using Fiedler vector representation, which allows our approach to be used for general cylindrical modeling and feature extraction. Furthermore, the fold segmentation allows for more consistent cuts along the colon surface which in turn provides more accurate flattened visualizations.
We have provided a thorough evaluation of our fold segmentation results by using the metric on 20 manually labeled datasets (10 multiorientation colon pairs). We have compared our method with other registration algorithms based on the computed registration error metric and have found our method to provide a significant improvement. Finally, we have also provided visual verification of our results on complete supine and prone colon pairs.
In the future, we will leverage the Fiedler vector level set approach for polyp detection in a given dataset and create a more integrated framework for computeraided detection based on our registration results. It has been shown in previous works that the performance and accuracy of the computeraided detection can be increased dramatically if the detection of polyps is done separately on haustral folds and endoluminal walls. We will use our fold segmentation framework to build this dichotomy and to help improve the accuracy of the CAD algorithms. We will also deploy our method in longitudinal analysis to register colons for the same patient across multiple visits, rather than just for single patient visits. Moreover, in the future, we will also cater to the severe underdistended colon cases where the segmentation can result in more than one colon segment per patient dataset. These cases can vary in complexity from a simpler (single supine/prone colon segment being registered to multiple supine/prone colon segments) to a more cumbersome (multiple supine/prone colon segments being registered to multiple supine/prone colon segments) scenario.
Acknowledgements.
The datasets are courtesy of Stony Brook University Hospital (SBUH) and Dr. Richard Choi, Walter Reed Army Medical Center. We would like to thank Dr. Kevin Baker of SBUH for his help in this project. This work has been partially supported by NSF grants DMS1418255, CNS0959979, IIP1069147, CNS1302246, and the Marcus Foundation.(a) Fold detection algorithm using heat diffusion and fuzzy Cmeans clustering [5] 
(b) Our fold detection and segmentation algorithm 
(a) Supine  (b) Prone 
(c) Supine  (d) Prone 
(e) Supine  (f) Side 
Appendix
A.1. Previous Fiedler vector colon registration
We have compared our work against a similar method using Fiedler vector introduced by Lai et al. [16]. Lai et al.’s automatic registration algorithm is dependent on the detection of two landmarks, namely the splenic and hepatic flexures, once the Fiedler vector has been computed on both the supine and prone models. Theoretically, the hepatic flexure forms the topmost point of the ascending colon and the splenic flexure forms the topmost point of the descending colon; this is the assumption of the Lai et al. algorithm as well. However, on real patient colon datasets that we tested, splenic and hepatic flexure detection is a much more complex problem and cannot simply be characterized by the local maximum zcoordinate near one Fiedler vector extrema or the other (as described in Lai et al.’s paper) due to the inflection points and the large variation in the troughs and ridges of the folds at the top. In general, local maximum zcoordinates can span several folds at the top and hence, the flexures detected within an epsilon value (using Lai et al.’s published approach) on the supine and prone colon models can vary considerably.
Due to this issue on some datasets when using Lai et al.’s registration algorithm, rather than improving the results (over our global registration baseline, which is essentially the registration with epsilon value of 0 in their case) once the landmarks (splenic and hepatic flexures) are detected, the results deteriorated in six out of our ten colon pair datasets (with an average registration distance error of 14.19 mm) due to the incorrect detection of the flexures on the corresponding colon pair (supine/prone/side) datasets. The average registration distance error for our complete algorithm, including the fold matching, was 5.24 mm for the same datasets.
Some of the registration problems that we encountered with the Lai et al. algorithm are highlighted in Figure 12 with three of our datasets. In Figure 12(b) for example, two local maxima are detected on the prone model within the value. This case is not handled based on the published Lai et al. approach. To compute the results, we picked the one closer to the maximum Fiedler vector value (as done by the Lai et al. approach for the supine but not for the prone model).
Moreover, the details of Lai et al.’s published algorithm do not deal with the flipping of the Fiedler vector between the supine and prone colon pairs. If the Fiedler vector values are flipped, this can lead to the cecumhepatic segment being registered with the rectumsplenic segment on the corresponding model, which gives an erroneous result. To avoid the flipping in our case, we manually make the Fiedler vector minimum and maximum values consistent between the colon pairs.
A.2. Previous fold segmentation using fuzzy Cmeans
We provide a visual comparison of the detection accuracy of our algorithm with a previous fold detection algorithm using heat diffusion and fuzzy Cmeans clustering [5] in Figure 11. As can be observed from a visual inspection of the results, our method succeeds in segmenting the individual folds, whereas the fuzzy Cmeans clustering often results in multiple folds being segmented together as a single fold.
References
 B. Acar, S. Napel, D. Paik, P. Li, J. Yee, C. Beaulieu, and R. Jeffrey. Registration of supine and prone CT colonography data: Method and evaluation. Radiology, 221(332):332, 2001.
 B. Acar, S. Napel, D. S. Paik, P. Li, J. Yee, R. B. Jeffrey, Jr., and C. Beaulieu. Medial axis registration of supine and prone CT colonography data. Proc. of Engineering in Medicine and Biology Society (EMBS), pp. 2433–2436, Oct. 2001.
 A. V. Bartrolí, R. Wegenkittl, A. König, and E. Gröller. Nonlinear virtual colon unfolding. IEEE Visualization, pp. 411–420, 2001.
 D. Boscaini, D. Eynard, D. Kourounis, and M. M. Bronstein. Shapefromoperator: Recovering shapes from intrinsic operators. Computer Graphics Forum, 34(2):265–274, 2015.
 A. S. Chowdhury, S. Tan, J. Yao, and R. M. Summers. Colonic fold detection from computed tomographic colonography images using diffusionFCM and level sets. Pattern Recognition Letters, 31(9):876–883, 2010.
 M. K. Chung, S. Seo, N. Adluru, and H. K. Vorperian. Hot spots conjecture and its application to modeling tubular structures. Machine Learning in Medical Imaging, pp. 225–232, 2011.
 A. H. de Vries, R. Truyen, J. van der Peijl, J. Florie, R. E. van Gelder, F. Gerritsen, and J. Stoker. Feasibility of automated matching of supine and prone CTcolonography examinations. British Journal of Radiology, 79:740–744, Sept. 2006.
 K. Gebal, J. A. Bærentzen, H. Aanæs, and R. Larsen. Shape analysis using the auto diffusion function. Computer Graphics Forum, 28(5):1405–1413, 2009.
 S. F. F. Gibson. Constrained elastic surface nets: Generating smooth surfaces from binary segmented data. Proc. of MICCAI, pp. 888–898, 1998.
 K. C. Gurijala, R. Shi, W. Zeng, X. Gu, and A. Kaufman. Colon flattening using heat diffusion riemannian metric. IEEE Transactions on Visualization and Computer Graphics, 19(12):2848–2857, 2013.
 S. Haker, S. Angenent, A. Tannenbaum, R. Kikinis, G. Sapiro, and M. Halle. Conformal surface parameterization for texture mapping. IEEE Transactions on Visualization and Computer Graphics, (2):181–189, 2000.
 L. Hong, S. Muraki, A. Kaufman, D. Bartz, and T. He. Virtual voyage: Interactive navigation in the human colon. Proc. of SIGGRAPH, pp. 27–34, 1997.
 W. Hong, F. Qiu, and A. Kaufman. A pipeline for computer aided polyp detection. IEEE Transactions on Visualization and Computer Graphics, 12(5):861–868, Sept. 2006.
 A. Huang, D. Roy, M. Franaszek, and R. M. Summers. Teniae coli guided navigation and registration for virtual colonoscopy. IEEE Visualization, pp. 279–285, Oct. 2005.
 A. Kovnatsky, M. M. Bronstein, A. M. Bronstein, K. Glashoff, and R. Kimmel. Coupled quasiharmonic bases. Computer Graphics Forum, 32(24):439–448, 2013.
 Z. Lai, J. Hu, C. Liu, V. Taimouri, D. Pai, J. Zhu, J. Xu, and J. Hua. Intrapatient supineprone colon registration in CT colonography using shape spectrum. Medical Image Computing and ComputerAssisted Intervention, pp. 332–339, 2010.
 P. Li, S. Napel, B. Acar, D. S. Paik, R. B. Jeffrey, Jr., and C. F. Beaulieu. Registration of central paths and colonic polyps between supine and prone scans in computed tomography colonography: Pilot study. Medical Physics, 31(10):2912–23, Oct. 2004.
 L. Liu, L. Zhang, Y. Xu, C. Gotsman, and S. J. Gortler. A local/global approach to mesh parameterization. Computer Graphics Forum, 27(5):1495–1504, 2008.
 R. Liu and H. Zhang. Segmentation of 3D meshes through spectral clustering. Pacific Conference on Computer Graphics and Applications, pp. 298–305, 2004.
 J. Marino and A. Kaufman. Planar visualization of treelike structures. IEEE Transactions on Visualization and Computer Graphics, 22(1):906–915, 2016.
 J. Marino, W. Zeng, X. Gu, and A. Kaufman. Context preserving maps of tubular structures. IEEE Transactions on Visualization and Computer Graphics, 17(12):1997–2004, 2011.
 M. Meyer, M. Desbrun, P. Schröder, and A. H. Barr. Discrete differentialgeometry operators for triangulated 2manifolds. Visualization and mathematics III, pp. 35–57, 2003.
 S. Nadeem, Z. Su, W. Zeng, A. Kaufman, and X. Gu. Spherical parameterization balancing angle and area distortions. IEEE Transactions on Visualization and Computer Graphics, PP(99):1–14, 2016.
 D. Nain, S. Haker, W. E. L. Grimson, E. Cosman Jr., W. W. Wells, H. Ji, R. Kikinis, and C. F. Westin. Intrapatient prone to supine colon registration for synchronized virtual colonoscopy. Proc. of MICCAI, pp. 573–580, 2002.
 M. Ovsjanikov, M. BenChen, J. Solomon, A. Butscher, and L. Guibas. Functional maps: a flexible representation of maps between shapes. ACM Transactions on Graphics (TOG), 31(4):30, 2012.
 W. Plishker and R. Shekhar. Virtual colonoscopy registration regularization with global chainmil. Proc. of MICCAI Workshop on Virtual Colonoscopy, pp. 116–121, 2008.
 E. L. Schwartz and B. Merker. Computeraided neuroanatomy: Differential geometry of cortical surfaces and an optimal flattening algorithm. IEEE Computer Graphics and Applications, 6(3):36–44, 1986.
 J. Shi and J. Malik. Normalized cuts and image segmentation. IEEE Transactions on Pattern Analysis and Machine Intelligence, 22(8):888–905, 2000.
 J. W. Suh and C. L. Wyatt. Deformable registration of supine and prone colons for computed tomographic colonography. Journal of Computer Assisted Tomography, 33(6):902–911, Nov. 2009.
 Y. Umemoto, M. Oda, T. Kitasaka, K. Mori, Y. Hayashi, Y. Suenaga, T. Takayama, and H. Natori. Extraction of teniae coli from CT volumes for assisting virtual colonoscopy. Proc. of SPIE Medical Imaging, 6916:69160–D–1–69160–7, 2008.
 S. Wang, J. Yao, J. Liu, N. Petrick, R. L. Van Uitert, S. Periaswamy, and R. M. Summers. Registration of prone and supine CT colonography scans using correlation optimized warping and canonical correlation analysis. Medical Physics, 36(12):5595–5603, 2009.
 Z. Wang, Z. Liang, L. Li, B. Li, D. Eremina, and H. Lu. An improved electronic colon cleansing method for detection of colonic polyps by virtual colonoscopy. IEEE Transactions on Biomedical Engineering, 53(8):1635–46, Aug. 2006.
 Z. Wei, J. Yao, S. Wang, and R. M. Summers. Teniae coli extraction in human colon for computed tomographic colonography images. Virtual Colonoscopy and Abdominal Imaging. Computational Challenges and Clinical Opportunities, pp. 98–104, 2010.
 J. Yao and R. M. Summers. Detection and segmentation of colonic polyps on haustral folds. IEEE International Symposium on Biomedical Imaging: From Nano to Macro, pp. 900–903, 2007.
 W. Zeng, J. Marino, K. C. Gurijala, X. Gu, and A. Kaufman. Supine and prone colon registration using quasiconformal mapping. IEEE Transactions on Visualization and Computer Graphics, 16(6):1348–1357, 2010.
 L. Zhao, C. P. Botha, J. O. Bescos, R. Truyen, F. M. Vos, and F. H. Post. Lines of curvature for polyp detection in virtual colonoscopy. IEEE Transactions on Visualization and Computer Graphics, 12(5):885–892, Sept. 2006.
 H. Zhu, M. Barish, P. Pickhardt, and Z. Liang. Haustral fold segmentation with curvatureguided level set evolution. IEEE Transactions on Biomedical Engineering, 60(2):321–331, 2013.
 H. Zhu, L. Li, Y. Fan, Q. Lin, H. Lu, X. Gu, and Z. Liang. Automatic teniae coli detection for computed tomography colonography. SPIE Medical Imaging, 7963:79632N–1–79632N–7, 2011.