Template-Based Active Contours
We develop a generalized active contour formalism for image segmentation based on shape templates. The shape template is subjected to a restricted affine transformation (RAT) in order to segment the object of interest. RAT allows for translation, rotation, and scaling, which give a total of five degrees of freedom. The proposed active contour comprises an inner and outer contour pair, which are closed and concentric. The active contour energy is a contrast function defined based on the intensities of pixels that lie inside the inner contour and those that lie in the annulus between the inner and outer contours. We show that the contrast energy functional is optimal under certain conditions. The optimal RAT parameters are computed by maximizing the contrast function using a gradient descent optimizer. We show that the calculations are made efficient through use of Green’s theorem. The proposed formalism is capable of handling a variety of shapes because for a chosen template, optimization is carried with respect to the RAT parameters only. The proposed formalism is validated on multiple images to show robustness to Gaussian and Poisson noise, to initialization, and to partial loss of structure in the object to be segmented.
Active Contours, Snakes, Affine matching, Contrast function, Shape constraint, Image segmentation.
The problem of segmentation in computer-aided image analysis relates to outlining the boundaries of the object of interest. The optimal boundary/contour cannot often be solved in closed form. Instead, it is solved for iteratively using optimization techniques. The contour evolves from a user-specified initialization, and converges to the object boundaries. Such contours are referred as active contours or snakes, a term coined by Kass et al. .
Snakes are evolutionary contours that move under the influence of certain forces: gradients or functions derived therefrom, which direct the contour to the boundary of interest — internal forces, which ensure certain degree of smoothness or regularity of the curve, these forces are referred as Snake energies. The novelty in snake design lies in the specification of suitable energies and in representation of the curve [5, 2, 6, 3, 7, 4]. The key challenges in segmentation are low signal-to-noise ratio, occlusion, non-uniformity of regional intensities, broken or diffuse boundaries or a combination of them. To mitigate these problems, a common approach is to incorporate prior shape information relating to the object of interest in the optimization framework [9, 8, 10, 11, 12]. The effectiveness of incorporating shape priors finds application in medical image segmentation [9, 8, 13, 14, 15, 16], where the objects of interest exhibit a certain degree of shape specificity.
1.1 Prior Art
Cootes et al.  proposed a method to incorporate shape-prior information in the form of a model derived from a set of training samples, which is assumed to reliably capture variations in the shape of the object of interest. However, this technique relies on image training data, which may be expensive, difficult, and in some cases not reliable. In the absence of such data, segmentation using a single shape prior is an alternate option. We briefly review previous literature on segmentation using a single shape prior, follow by our approach.
Techniques on single shape prior based segmentation employ regularization energies in addition to image-derived energies for optimization [10, 11, 12, 8, 18, 17, 13, 14]. Regularization energies penalize segmenting curve for deviating from the geometric transformation space (GTS) of the reference curve, thereby imposing a shape constraint. The regularization functional is a shape distance metric (SDM) that measures the deviation in shape between the segmenting curve and the reference. Optimization within this regularization framework is carried out iteratively and in two stages [10, 11, 12, 18, 17]. In the first stage, the curve is evolved to optimize the image derived energies, while in the second stage, the curve is modified to satisfy the shape constraint by minimizing the SDM between the evolving curve and the reference in the GTS space governing the reference curve. As gradient descent is employed for minimization, the local optimality of the solution outputted by the algorithm may not satisfy the shape constraint.
Thévenaz et al. proposed a technique that mitigates the two stage approach with in the regularization framework, for segmenting circular  and elliptical blobs . The active contour comprises two concentric circles (snakuscule) or ellipses (ovuscule) as a shape template for segmentation. In both cases, the snake energy is defined as the contrast between the pixel intensities in the inner contour and those lying in the annular region between the inner and outer contours. Snakuscules are parametrized by a pair of diametrically opposite points, which are constrained to always possess the same ordinate by means of a regularization functional. Ovusucules are specified by three points on the inner ellipse. However, the shape specific parameterization and regularization employed is not easily generalizable to any shape. A regularization-free approach was proposed by Chen et al. . They combined image registration techniques with segmentation. This algorithm exhaustively search for optimal similarity transform (ST) parameters globally, to match reference shape and object of interest. The technique is iterative and each iteration consists of two steps. In the first step, the translation parameters are updated following a global search in the image pixel coordinates.
In the second step, rotation and scaling parameters are updated using similar global search in the log polar representation of the image. The computations for the search are made efficient by the use of fast Fourier transform (FFT).
In , we extended the formalism of Thévenaz et al. [13, 14] to concentric rectangular templates for segmentation of gel electrophoresis images. In , we proposed a unified approach for designing snakuscules and ovuscules from a pair of concentric circles and optimizing for a restricted affine transformation (RAT). RAT has five degrees of freedom allowing for different scaling along horizontal and vertical directions, rotation, and translation. This unified approach is simpler to formulate and optimize, while retaining the shape specificity.
1.2 This Paper
In this paper, we consider a template based active contour formalism for segmentation. Specifically, we extend the formalism in [13, 14, 15, 16] to a generic shape template, which may even be non-convex. The salient features of our formulation are:
The mathematical formalism is shape-independent, devoid of explicit regularization energies and requires calculation of a fixed five parameters irrespective of the shape.
Using Green’s theorem, we convert region integrals into contour integrals. As contour integrals require less number of computations compared to region integrals, this step will result in computational savings.
We will show through experiments that proposed active contour formalism is fast enough to enable real time deployment.
We retained the nested structure in [13, 14] for specifying the shape template. The inner contour of the snake lies completely inside the outer contour. The user-defined shape template serves as the inner contour while the outer contour is obtained either by scaling the inner contour (for convex shapes) or by a manual specification. The contours are parametrized using B-splines; linear B-splines are used for representing polygonal shapes and cubic B-splines are used for designing smooth closed contours. The snake energy is a normalized contrast function, which measures the difference in average intensities between the inner contour and the annular region between the inner and outer contours. We modify the energy proposed by Thévenaz et al. [13, 14], by adding a normalization using the areas enclosed by the inner and outer contours to ensure stability of fit and to eliminate explicit regularization. Upon optimization, the inner contour locks on to the object of interest. The outer contour guides the computation of the local contrast metric. The optimization of the snake energy is carried out with respect to the RAT parameters of the shape template. However, the proposed formalism can also be adapted to a generic transformation.
The rest of the paper is organized as follows. In Section 2.1, we briefly review some key aspects of B-spline kernel based curve parameterization. We also give a formulation of the active contour and introduce the snake energy. In Section 3, we address the active contour optimization problem. The computational aspects are discussed in Section 4. In Section 5, we present experimental results on synthesized images, medical, and biological images. Robustness to Gaussian, Poisson noise conditions, and to initialization are also reported.
2 Active Contour Formulation
Active contour design involves two major steps: design of the contour, and design of the energy functional. We design the contour using a B-spline based parameterization which have some interesting properties such as non-negativity, compact support, minimum curvature. Shifted cubic B-splines also form Riesz basis. We briefly review some of those properties to justify their usage in the proposed formalism.
Usage of B-splines will lead to a case where we can use the parameters of the curve (referred as spline knots) for optimizing the RAT parameters, instead of the curve itself. This will be shown in Section 3, when we derive the derivatives of RAT parameters. As we are going to work on the spline coefficients, we review the Riesz basis property satisfied by B-splines, which allows us to operate on the discrete spline coefficients.
2.1 Riesz basis property
Consider a signal in the space spanned by the shifted cubic B-spline. We can write , where is a cubic-B-spline. As the shifted cubic B-splines are linearly independent, they form Riesz basis and satisfy the following Riesz basis conditions:
This property ensures that a significant change in the value of coefficients () leads to a significant change in the shape of the curve and vice versa.
2.2 B-spline-based template parametrization
Consider a pair of inner and outer contours parametrized as and , respectively, where is the independent variable. The subscripts and indicate the inner and outer contours, respectively. The functions , , , and , are expressed in terms of a basis function and its integer translates. The structure of the parametrization is along the lines of [4, 22], and results in a vector representation for the shape template:
where , are the coefficient vectors, which are referred to as spline coefficients in spline literature. The Riesz basis property (1) of the B-splines allows us to operate on the knots to control the contour. Closed curves result in periodic and , for which we have the equivalent representation:
where is -periodic; , , , are the periodized coefficient sequences with period ( being the number of knots). In Figure 1, we show a shape template designed using cubic B-splines.
For convex shapes centered at the origin, the outer contour () may be chosen as a scaled version of the inner contour () that is, , is a scalar greater than unity. In this case, the coefficients of the outer contour are times those of the inner contour. The parameter gives direct control over the width of the annular region and determines the region over which the local contrast function is computed. For non-convex shapes, the outer contour must be parametrized independently of the inner contour.
2.3 From shape templates to active contours
Template based active contours are derived from the shape template in (2) according to the transformation:
where , and and are the inner and outer contours of the template based active contour respectively. and represent the scale parameters along axial and perpendicular directions, respectively; represents the angle of rotation (clockwise), and are the translation parameters. The nested contour design is general as the outer and inner contours can be chosen differently.
2.4 Snake energy specification
In this section, we first provide a short structural overview of snake energies in general and discuss the snake energies that are relevant to our work. We then propose a new snake energy and establish the connection shared by the proposed formulation with other snake energies.
Most snake energies [1, 5, 2, 6, 3, 7, 4, 11, 12, 13, 14, 15, 16] can be represented by (5). The snake energy () comprises image and contour derived functions. The contour related functions () are aimed at ensuring stability and smoothness of the evolving contour, while, the image related functions are responsible for driving the snake towards the object boundaries:
The image related energies are broadly classified into region () and gradient energies (). Region energies rely on image statistical measures for segmentation, whereas gradient energies finely segment object boundaries by capturing steep transitions in pixel intensity using derivative operators. Typically, a linear combination of the energies is used for segmentation. Assigning the weight factors () for the combination is not straightforward and may require tuning on an image by image basis.
Recently, there is growing interest to use only region based energies as the image derived component of the snake energy [23, 24, 25, 26]. The motivation is due to the fact that gradient energies highlight noise and minor artifacts in the image, thereby degrading the segmentation performance of the snake. Also, it is known that region energy-based snakes are more robust to noise and initialization than gradient energy snakes as they rely on non-local measures. In view of these observations, we favor the use of region based energies in our formulation. We briefly review a few popular region based energies and gradually introduce our snake energy.
For an active contour initialized on an image , Chan and Vese  proposed the following energy model, which is a special case of Mumfard-Shah functional:
comprises two terms and , where represents curve regularization energies: length and area of the evolving contour C. The region within the evolving contour (C) is considered as foreground and outside the contour as background. penalizes the variances of foreground and background regions. and are the mean intensity values of the foreground and background respectively. Segmentation is achieved based on the assumption that the foreground and background regions are distinguishable based on and . Yezzi et al.  proposed the mean separation energy () shown in (7) for segmentation. Instead of approximating the foreground and background regions with constant pixel intensities, the formulation aims to segment the image by maximally separating their mean intensity values ( and ). The formulation includes an arc length regularizer, penalizing the active contour for encircling small noisy regions:
Thévenaz et al. proposed active contours to segment locally bright blob like structures. They considered concentric circles and ellipses. They used a local contrast function that is similar to the mean separation energy.
In (2.4), is region bound by the inner contour and , by the outer contour, with . The contrast term in (2.4) forces the snake to maximize the difference in average intensity between the annular region and inner contour. The parameterization adopted by Thévenaz et al. yields non-unique solutions, which are resolved with the help of the regularization term (). The shape specific nature of the parameterization leads to a regularization functional that is also shape dependent. The contrast function considers the area bounded by the inner contour as foreground and the annular region as background. The area occupied by these regions on the image varies as the snake evolves. This is in contrast to the techniques in [24, 26] where the union of the foreground and background region is the constant user specified bounding box on the image.
Motivated by the definition of local contrast proposed by Thévenaz et al., we propose a modified version of local contrast. Consider the nested snake design introduced in Section 2.3. For a snake initialized on an image , let and be the regions enclosed by and , respectively (as shown in Figure 2). We define our snake energy as follows:
where and are the areas enclosed by the inner and outer contours of the shape template respectively, is the area of the annular region of the shape template. Minimizing enables the active contour to lock on to objects that are brighter than their immediate surroundings, and maximizing it would have the opposite effect. The normalization term in the denominator removes the ambiguity suffered by active contours. To explain the benefit of normalization, consider the example shown in Figure 3, which consists of a dark square in a bright background. In order to capture the dark square, needs to be maximized. Optimizing with the help of the normalization term results in the outline shown in Figure 3(a), whereas (b) and (c) show the results obtained without normalization. For all three cases, the value of is the same, but the fit is optimal when is normalized by , and suboptimal otherwise. It is also desirable for the snake to remain stationary when placed on a constant region. We see that if (a constant), then , which means that the snake is inert on constant regions. Due to shape-specific nature of the proposed snake and one-one relationship between shape and parameters, our formulation does not require additional regularization energies. For symmetric shapes this one-one relationship does not hold. However, there are finite number of configurations of RAT for the same shape. All these local minima are strict local minima and the gradient descent technique will converge to one of the shapes without oscillating. The only exception is circular active contour, for which uncountably-infinite configurations of correspond to the same shape; we resort to the method proposed by us in  for segmenting circular shapes.
2.5 Fischer ratio interpretation of the snake energy
We next show that optimizing in (2.4) is analogous to maximizing the Fischer ratio  under certain conditions. Consider an image that contains a homogeneous bright object to be segmented in the presence of a relatively darker background. Assume that the image is corrupted with Gaussian noise with mean and variance . The distributions of within the bright object and outside of it are both Gaussian, with same variances, but with different means. The problem of outlining is therefore equivalent to a Gauss-Gauss detection problem , which can be solved using a Fischer ratio approach.
The energy function in (2.4) can be modified as
For a snake initialized on an image, if , minimizing the snake energy (2.5) is equivalent to maximizing , where are the mean intensities of in the annular region and inner contour, respectively. The Fischer discriminant for this Gauss-Gauss detection problem is , where , are variances of in the annular region and the inner contour, respectively. Since under uniform Gaussian noise conditions, the proposed snake energy optimization is directly related to Fischer ratio maximization.
2.6 Advantages of contrast energy over squared contrast energy
We assumed that the object is locally brighter or darker than its immediate surroundings and correspondingly selected the sign of the energy in (2.4) for minimization. By squaring the energy and adding a negative sign in (2.4), we would obtain the sign-agnostic version of the contrast energy, called squared contrast energy (), similar to mean separation energy in (7). Though would overcome the problem of choosing the sign of the energy, we show that the energy in (2.4) provides certain advantages over . We perform two experiments to compare local contrast function () in (2.4) and local squared contrast energy (). We show that local contrast function is preferable to local squared contrast energy.
In the first experiment, we consider a bright ring shaped structure against a dark background as shown in Figure 4(a) having an inner and outer radius of 100 and 141 ( 100) pixels, respectively. We choose a circular template (ratio of inner to outer radius of template fixed at 1:) for segmentation and initialize the snake as shown in Figure 4(a). We plot the energies and as the template is uniformly scaled (i.e. ) between (0, 400) about the center of the bright ring. The plot is shown in Figure 5, where the green curve corresponds to and the red curve to . The energies are plotted against the radius of the inner contour of the shape template on the horizontal axis. From the plot we observe that for a snake initialized as shown in Figure 4(a) and having an initial inner contour radius , minimization of causes the snake to converge to the result shown in Figure 4(b). Maximization causes the snake to converge to the result shown in Figure 4(c). However, for snake energy, the final result is dependent on the initialization. Specifically, the result in Figure 4(c) requires an initialization between 100 and 120 pixels while the result in Figure 4(b) requires an initialization beyond 120 pixels, shown as ‘b’ in Figure 5.
In second experiment, we consider an image containing two overlapping circular objects as shown in Figure 4(d). The pixel intensity value of the two objects are chosen differently while the region of overlap and the background are set to the average value. We consider a snake (derived from a circular template) initialized on the intersection of the two circular objects as shown in Figure 4(d). Minimizing gives the result in Figure 4(e), while maximizing gives the result in Figure 4(f). However, the relative area of overlap with each object at the time of initialization determined the final output for optimizing snake. In both experiments, had a larger basin of attraction as compared to , which resulted in better control of the final segmentation output of the snake.
We further substantiate that snakes will have larger basin of attraction compared to snakes by showing that some spurious local minima are present in snakes which are not local minima for snakes. Consider the local maxima for snakes (for 0). These local maxima will become local minima for as (). These additional minima locations correspond to sub-optimal solutions and decrease the basin of attraction. Also, the snakes optimizing the local contrast carry more information (sign of energy function) and hence, their segmentation performance is superior to snakes. To conclude, local square contrast function is more sensitive to initialization and less controllable compared to local squared contrast function. We use local contrast () as energy functional in the rest of the paper.
3 Active Contour Optimization
We next proceed with the optimization of the snake energy with respect to the RAT parameters and . We need the partial derivatives of with respect to the parameters to enable gradient-descent optimization. The energy involves double integrals, and the parameters define the boundaries. These, in turn, define the regions and over which the integrals are evaluated. Direct optimization seems to be a difficult task; however, thanks to Green’s theorem the surface integral can be expressed conveniently as a contour integral. This methodology simplifies the calculations of the partial derivatives and consequently, we have to deal with contour integrals.
We first perform a coordinate-axes transformation from to such that
which converts the image function to
We could have preserved the image as it is, and performed the axes transformation on the inner and outer contours, but we found that transforming the image to the new coordinate system simplifies the calculations and also results in brevity of notation. The resulting energies and take the form:
For brevity of notation, we have dropped the parameter and denoted as , respectively. We next invoke Green’s theorem for computation of and . We show the calculations only for , and those for follow analogously. Applying Green’s theorem to leads to the following expression:
is a function of , which are in turn functions of the parameters of the contours. The partial derivative of with respect to is given by the chain rule as follows:
Similarly, the partial derivatives of with respect to the parameters and are given as,
The generalized expressions for the partial derivatives of in (2.4) with respect to the RAT parameters are given below:
We simplified (15) by substituting the values of partial derivatives of and with respect to RAT parameters. Expressing them in terms of , , , and the image , we have
4 Computational complexity
where , and = denotes the number of knots. and
The main steps involved in the computation of partial derivatives are given below.
The continuous-domain expressions for and in (16) have to be discretized for the purpose of computation. This is achieved by sampling the contour finely, which leads to the approximation:
where is the number of knots, is the number of discretization steps over , leading to an overall number of points along the closed contour. The terms involving are spline-dependent and can be precomputed. Also, is common in the partial derivatives with respect to , , and ; is common in the partial derivatives with respect to and .
The use of gradient descent for optimization warrants a stopping criterion, for which we chose to monitor the snake energy . Let be the value of in the iteration during snake optimization, and the value at iteration. Denote . If the number of zero-crossings of crosses a predefined threshold, then we stop the gradient descent. Zero crossings is a good indicator as RAT parameters oscillate about the optimal solution, a behavior typical to gradient descent.
5 Experimental Results
In this section, we examine the performance of the proposed technique on multiple images, and the effect of initialization, noise, and occlusion.
5.1 Effect of initialization
Consider a Shepp-Logan phantom image  shown in Figure 6, which is widely used for validation of computed tomography reconstruction algorithms. The initializations shown in Figure 6(a)–(d) converge to the optimal solution shown in (e), whereas a poor initialization (as shown in Figure 6(f), where the outer contour overlaps with the dark boundary) leads to slightly suboptimal solution shown in (g). Also, as expected, a good initialization (Figure 6(b)) results in faster convergence. The initializations shown in (a)–(d) converged to (e) in 121, 71, 74, and 77 milliseconds, respectively. The execution times are given with respect to an ImageJ implementation of the algorithm running on Mac OSX 2.66 GHz, Intel Pentium Dual CPU T3400 @2.16GHz. A similar set of experiments for a non-convex shape was carried out using MR images shown in the second row of Figure 7. Since the desired object of interest in this case is dark and the surrounding background is bright, we negate the energy function . The convergence times in milli-seconds for the initializations shown in (a)-(d) are 110, 78, 95, and 83, respectively. We attribute the fast convergence to the local nature of snake energy and direct computation of line integrals in (16) and (17).
5.2 Robustness to Poisson and Gaussian noise
We validate the performance of the proposed technique on noisy images, by considering two types of noise distributions — Gaussian noise and Poisson noise. Gaussian noise is independent of the image data, additive and stationary, whereas Poisson noise is dependent on the image data and is multiplicative. Images in the top row of Figure 7 (Shepp-Logan phantom ) are corrupted with Gaussian noise and the images in the bottom row with Poisson noise. Standard deviations of the Gaussian noise added to the four images are 25, 50, 75, and 100 respectively. Peak signal-to-noise ratios (PSNRs) computed therefrom are indicated in the figure captions. For Gaussian noise images, we used the standard PSNR definition given in , and for Poisson noise images, the definition given in . Section 2.5 gives theoretical guarantees for gaussian noise corrupted images, but from the experimental results, we observe that performance of active contour technique is robust even to Poisson noise.
5.3 Elliptical templates
Many of the bio-medical cells can be approximated with an ellipse . We present results using two elliptical shaped templates: a spline synthesized ellipse for which we use the partial derivatives given in equations (16) and (17) directly, an exact ellipse for which we use the parameterization given in  and use the equation (15). First row of Figure 8 corresponds to a hand drawn ellipse and the second row correspond to a parametric ellipse. We considered these two parameterizations to show that the algorithm performance is not critically dependent on the parameterization. These are fundus images taken from . The presence of veins in the fundus anatomy obscures snakes whose energies rely on image derivative functions. From the results shown in Figure 8, we infer that the proposed technique is less affected by the vein structures or the parameterization and captures the near elliptical shape of the fundus outline.
We next conducted an experiment on a mammogram taken from the mini-MIAS database , which contains some benignant and malignant tumor images. The images (filename: mdb028) containing lumps as shown in Figure 9, were specifically classified in the database as circular in shape with a certain center and radius. We employed the elliptical shape template to localize the mass and measured the resulting dimensions of the fit. For the top row, the converged ellipse fit was found to have a semi-major and semi-minor axes of 58 and 56 pixels, respectively. The radius of the approximate circumscribing circle was given in the database as 56 pixels. The center of converged contour matches the one given in the database. For the bottom row, semi-major and semi-minor axes are 64 and 42 pixels, whereas the true radius is 68. The center of the converged snake differs by 4 pixels in both horizontal and vertical directions from that given in the database.