- Email: [email protected]
NURBS snakes
Image and Vision Computing 21 (2003) 551–562 www.elsevier.com/locate/imavis
NURBS snakes Ravinda G.N. Meegama, Jagath C. Rajapakse* School of Computer Engineering, Nanyang Technological University, Nanyang Avenue 639798, Singapore Received 12 November 2002; received in revised form 7 March 2003; accepted 18 March 2003
Abstract A new formulation of energy minimising splines, based on Non-Uniform Rational B-Spline (NURBS) curves, is presented. The existing snake models do not facilitate the local control of their shapes without increasing the number of control points. The introduction of NURBS renders more flexibility into the snake model due to the weighting parameters that are made to adapt according to the curvature of the contour. The proposed snake model was tested on both synthetic and real images; the new model demonstrated more local flexibility than the previous models and is capable of attracting to complex boundaries of the objects with increased accuracy. q 2003 Elsevier Science B.V. All rights reserved. Keywords: Active contours; Non-Uniform Rational B-Spline; Energy minimisation; Boundary extraction; Snakes
1. Introduction Extracting object boundaries accurately from images has become an increasingly important and a challenging task in the arena of image analysis research over the past two decades. This task is extremely difficult when the shape of the object exhibits complex features such as high curvature regions. To address this issue, Kaas, Witkin and Terzopoulos introduced a model-based active shape contour referred to as a snake [1]. A snake is initialised approximately around the object of interest by manually placing a set of discrete points, referred to as control points, that govern the movement of the contour. These points evolve in time according to internal forces (elastic forces) determined by the shape of the curve and external forces (image forces) obtained by local intensity statistics of the image. Since their inception, the behaviour of snakes near a highly convoluted and a complex boundary, such as the surface of the human brain, has been of highly interest [2 – 4]. Snakes have also been used extensively for segmentation of medical images such as measuring the area of leg ulcers [5], locating the boundaries of the head [6], bladder [7], heart [8] and liver [9]. A majority of these research focused on generalising the snake model and improving the convergence and stability of its evolution * Corresponding author. Tel.: þ 65-67905802; fax: þ 65-67926559. E-mail addresses: [email protected] (J.C. Rajapakse); [email protected] (R.G. Meegama).
[10 –14]. Several new snake models have emerged as a result of this, for example, topological adaptive snakes [15] and B-spline snakes [16]. However, these snake models are incapable of increasing the local flexibility of the contour without adding more control points. To attract a snake onto a boundary having a high curvature, the changes made to the positions of the control points near a particular segment of the snake should affect the shape of only the immediate neighbouring segments of the snake. The local control of an active contour is vital in cases where the shape to be extracted exhibits both high and low curvature regions. The discrete snakes [1] cannot adapt to complex shapes of object boundaries because the number of control points remains constant during the evolution. Although a discrete snake is simple to implement, the accuracy of the extracted boundary is limited because the snake is handicapped by the fixed number of control points. A deformation scheme that extends the geometrical flexibility of the snake to deal with complex shapes, referred to as topologically adaptive snakes (T-snakes), is revealed in [15]. Both discrete and T-snake suffer from the drawback that the movement of a single control point affects the entire length of the contour; a small change in the position of a control point tends to strongly propagate throughout the entire snake. This can be described as a global propagation of change of a particular control point. The B-spline snake [16], in which the contour is constructed from a parametric B-spline curve, has improved
0262-8856/03/$ - see front matter q 2003 Elsevier Science B.V. All rights reserved. doi:10.1016/S0262-8856(03)00066-0
552
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
local control of the contour. The curve is approximated by piecewise polynomial functions that define specific segments on the curve. Although B-spline snakes perform much better than the former models with respect to accuracy, the distinct control points may have to be duplicated to detect high curvature regions, resulting in increased redundant computations. NURBS-based deformable models, called D-NURBS (or Dynamic NURBS), have been investigated in [17] for computer aided geometric design. In this technique, the weights near a desired region have to be interactively adjusted by the designer in order to influence the shape of the curve. A NURBS-based active contour is outlined in [18] to segment objects in blurred digital images. This method, however, does not make use of the full potential of the NURBS representation of curves, such as shape modification with respect to the weights, to extract the complex boundaries. The main goal of our work is to address the problems encountered in the above methods by introducing a new snake model that is equipped with an additional parameter that can be automatically adjusted to change the local shape of the curve. Our proposed snake model, referred to as a NURBS snake, uses the NURBS representation of curves [19]; the weighting parameters in NURBS are automatically adjusted to control the flexibility of the curve at each control point. During the evolution, the weight associated with each distinct control point is modified depending on the curvature of the points on the relevant spline segment. Such a modification of weights influences the snake to attract towards high curvature regions without adding more control points. The paper is structured as follows: Section 2 presents major definitions in the traditional snake model. Section 3 outlines the principles of basis functions and B-splines. The NURBS snakes model is introduced in Section 4 followed by the implementation details in Section 5. Experimental details are presented in Section 6 and finally, Section 7 discusses in detail the performance of the NURBS model.
2. Snakes definitions This section gives the important definitions pertaining to active contour models and the associated energy terms. In this paper, we assume two-dimensional (2D) images.
[0,1] £ [0,1) ! V that parameterises the image coordinates ðx; yÞ by the arc parameter s and a time parameter t [ ½0; 1 such that vðs; tÞ ¼ ðxðs; tÞ; yðs; tÞÞT In what follows, in a curve vðs; tÞ; we sometimes omit either the time parameter t or the arc parameter s; or both, to simplify the notations according to the context. Definition 3. Internal Energy. The internal energy Eint at a point vðsÞ on a curve v is defined as Eint ðvðsÞÞ ¼ alv0 ðsÞl2 þ blv00 ðsÞl2
ð1Þ
where a; b [ {0; R }; v ðsÞ ¼ ›vðsÞ=›s and v ðsÞ ¼ › vðsÞ= ›s2 : þ
0
00
2
Definition 4. External Energy. The external energy Eext at a point vðsÞ on a curve v is defined as Eext ðvðsÞÞ ¼ 2gl7IðvðsÞÞl2
ð2Þ
where g [ Rþ ; IðvðsÞÞ denotes the intensity at the point vðsÞ and l·l denotes the Euclidean norm. Definition 5. Balloon Force. The balloon force applied on a point vðsÞ is defined as [20] ! l7IðvðsÞÞl2 bðvðsÞÞ ¼ d 1 2 nðvðsÞÞ ð3Þ k7IðvðsÞÞk2 where d [ R; nðvðsÞÞ is the unit normal vector at vðsÞ and k·k denotes the maximum norm of a vector such that k7IðvðsÞÞk ¼ max l7IðvðsÞÞl: s
Definition 6. Total Energy. The total energy E of a curve point vðsÞ is defined as EðvðsÞÞ ¼ Eint ðvðsÞÞ þ Eext ðvðsÞÞ:
ð4Þ
Definition 7. Snake. A snake is a dynamic curve that evolves over time attempting to minimise the total energy of its contour. If vðs; tÞ represents a point on a dynamic curve at a discrete time instance t; the point at time t þ 1; vðs; t þ 1Þ; is given by vðs; t þ 1Þ ¼ arg min Eðvðs; tÞÞ þ bðvðs; tÞÞ vðs;tÞ[n
where n is the set of all possible curve points vðs; tÞ: Definition 1. Curve. If V denotes a 2D image domain, a curve is a one-parameter vector-valued mathematical function v: [0,1] ! V such that an arbitrary point on the curve v is parametrically expressed as vðsÞ ¼ ðxðsÞ; yðsÞÞT ; s [ ½0; 1: Definition 2. Dynamic Curve. A dynamic curve is a twoparameter vector-valued mathematical function v:
The energy of the snake is defined in such a way that starting from a user-specified configuration near an object boundary, it deforms itself to conform to the nearest salient boundary by minimising the combined internal and external energies. The internal energy of the snake, Eint ; models the elasticity of the contour and makes it behave as a membrane and/or thin plate in which the amount of stretching and
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
flexibility is characterised by the parameters a and b; respectively. The external energy, Eext ; is determined by the intensity statistics of the image and leads the snake toward the required edges of the image. The balloon force, b, is applied to an energy minimised point on the contour and it forces the snake to overcome noise and other spurious features on the image [21 – 23].
3. Basis function representation A parametric curve vðsÞ can be defined as a linear combination of a given set of points, referred to as control points, and a set of polynomial functions fi : ½0; 1 ! R; i ¼ 0; 1; …; n 2 1; such that vðsÞ ¼
nX 21
fi ðsÞpi
ð5Þ
i¼0
When the polynomial functions form a basis for the vector space of all the polynomials of the desired degree and continuity, the polynomials are referred to as basis functions. The basis function fi can be a piecewise polynomial function defined on an interval, say ½a; b; with specified continuity constraints such that when the interval is subdivided by a sequence of real values, si where a ¼ s0 , … , sn21 ¼ b; fi is a polynomial within each subinterval ½sj ; sjþ1 Þ: In such a case, fi is referred to as a spline basis function. Spline basis functions provide local control over the shape of a curve if each basis function, fi ; is made non-zero only within a limited number of subintervals of s: As each pi is multiplied by the corresponding basis function fi (Eq. 5), a movement of the control point pi affects only the segment of the curve lying in the subinterval where fi ðsÞ – 0: Therefore, such a spline basis function representation can be used to model snakes having local flexibility. Definition 8. B-spline Curve. A B-spline curve v : ½0; 1 ! V of degree k [ Zþ is defined by vðsÞ ¼
nX 21
Ni;k ðsÞpi
knot, si ; relates the parameter s to each control point pi such that a translation of pi affects only the curve segment vðsÞ; s [ ½si ; siþkþ1 Þ: This local control comes from the Bspline basis functions in which Ni;k ðsÞ ¼ 0 for s ½si ; siþkþ1 Þ: B-spline P basis functions hold the property of partition of unity; i Ni;k ðsÞ ¼ 1; making vðsÞ an affine combination of its control points. Moreover, since ;i ¼ 0; 1; …; n 2 1; Ni;k ðsÞ $ 0; the B-spline curve vðsÞ is a convex combination of its control points, which makes a snake formulated using a B-spline curve always lie inside the convex hull formed by the control points. The points on a uniform B-spline curve have equally spaced parameterisation. Conversely, the parameterisation of points on a ‘Non-Uniform’ B-spline curve may not necessarily be equally spaced. A non-rational B-spline function provides only an approximation to conics whereas a ‘Rational’ Bspline consists of an extra parameter, wi [ Rþ ; called a weight associated with each control point pi that gives not only an extra degree of freedom for shaping the curve but also defines precisely any of the conics. Incidently, NURBS is the acronym for Non-Uniform Rational Bsplines [19].
4. NURBS snake This section gives the definitions related to the proposed NURBS snakes model. Definition 9. Dynamic NURBS Curve. If pti and wti denote the ith control point and its weight at time t; respectively, a dynamic NURBS curve v of degree k is defined as Xn21 t wi Ni;k ðsÞpti vðs; tÞ ¼ Xi¼0 n21 t wi Ni;k ðsÞ i¼0 It can also be written as vðs; tÞ ¼
i¼0
where the ith B-spline basis function of degree k; Ni;k ; is recursively defined as ( 1; if s [ ½si ; siþ1 Þ Ni;0 ðsÞ ¼ 0; otherwise
553
n21 X
Rti ðsÞpti
ð6Þ
i¼0
where Rti ðsÞ ¼
wti Ni;k ðsÞ n21 X t wj Nj;k ðsÞ j¼0
and ðs 2 si ÞNi;k21 ðsÞ ðs 2 sÞNiþ1;k21 ðsÞ þ iþkþ1 Ni;k ðsÞ ¼ siþk 2 si siþkþ1 2 siþ1 where si ; i ¼ 0; 1; …; n þ k; are positive real numbers referred to as knots. The vector (s0 ; s1 ; …; snþk ÞT is called the knot vector of the B-spline curve v. The value of the ith
The NURBS curve in Eq. (6) can be expressed in matrix from as vðsÞ ¼ pT RðsÞ where p ¼ ðpt0 ; pt1 ; …; pTn21 ÞT and RðsÞ ¼ ðR0 ðsÞ; R1ðsÞ ; …; Rn21 ðsÞÞT :The internal energy, external energy and
554
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
the balloon force for a NURBS curve, computed using the definitions given in Section 2, are given, respectively, by Eint ðvðsÞÞ ¼ alpT R0 ðsÞl2 þ blpT R00 ðsÞl; 2 Eext ðvðsÞÞÞ ¼ 2g7IðpT RðsÞÞ and bðvðsÞÞ ¼ d 1 2
l7IðpT RðsÞÞl2 k7IðpT RðsÞÞk2
ð7Þ
! nðpT RðsÞÞ:
ð8Þ
If pi ¼ ðxi ; yi ÞT ; the gradient operation 7I is defined as s ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 › › T 7Iðp RðsÞÞ ¼ wpI þ wpI ›xT RðsÞ ›yT RðsÞ where x ¼ ðx0 ; x1 ; …; xn21 ÞT ; y ¼ ðy0 ; y1 ; …; yn21 ÞT ; w is a smoothing kernel and p denotes the convolution operator. Definition 10. NURBS Snake. A NURBS snake is a dynamic NURBS curve vðs; tÞ whose energy is minimised with respect to the points on the curve such that a curve point at time t þ 1; vðs; t þ 1Þ is given by vðs; t þ 1Þ ¼ arg min ðEint ðvðs; tÞÞ þ Eext ðvðs; tÞÞÞ þ bðvðs; tÞÞ vðs;tÞ
Fig. 1. Movement of a curve point vðsÞ towards the desired control point when the associated weight is increased.
of the control point ptþ1 is updated such that weight wtþ1 i i wtþ1 ¼ wti þ h i
Kt kkðs; tÞk
ð11Þ
where the parameter h [ R controls the amount of attraction of the curve towards the control point and kkðs; tÞk is the maximum curvature of the closed NURBS curve. Note that in Eq. (11), the curvature Kt is divided by the maximum curvature in order to normalise the weight modification.
4.1. Weight-based control 5. Implementation This section presents how the weights of a NURBS snake, without moving or duplicating the relevant control points, effectively control the local shape of the contour based on its curvature properties. Let vðs; wi Þ be the family of curves obtained by changing the weight wi at a control point pi : Also, let vðs; wi ¼ 0Þ be the set of curves when wi is zero. Then, v can be expressed as [19] vðs; wi Þ ¼ vðs; wi ¼ 0Þ þ tðpi 2 vðs; wi ¼ 0ÞÞ
ð9Þ
where
t¼
Ni;k ðsÞwi nX 21
ð10Þ
Nj;k ðsÞwj
j¼0
From Eqs. (9) and (10), it is clear that as the weight wi associated with the control point pi is increased, t is increased, and hence, the point vðs; wi Þ moves toward pi along the vector pi 2 vðs; wi ¼ 0Þ as shown in Fig. 1. Similarly, vðs; wi Þ moves away from pi if t is decreased. Let pti denote the position of a control point with weight t wi at time t and kðvðs; tÞÞ; the curvature at a point vðs; tÞ: If the maximum curvature within the spline segment vðsÞ; s [ ½si ; siþpþ1 Þ; at time t; Kt ; exceeds a pre-defined value (say, the average curvature along the NURBS curve), the new
5.1. Evolution of points In the beginning, the control points of the snake are selected near the object boundaries to ensure that the snake does not converge to a local minimum. The selection of these initial points is based on a priori information or an interactive manual placement. These points can lie inside (or outside) the object; the balloon force is able to push (or pull) the snake accordingly. At a given iteration, new control points are dynamically added to (or removed from) the control polygon based on both the Euclidean distance between two consecutive control points and the curvature properties within the relevant spline segment. The removal of control points is due to the fact that only a fewer points are needed for the snake to expand at low curvature regions (in B-spline snakes, however, new points must be added at high curvature regions). The positions of the control points are now moved according to the minimal energy configuration of the spline and a NURBS curve is constructed based on these new locations of the control points. The major components of our algorithm is illustrated in Fig. 2. We used dynamic programming, as highlighted in [24], to evolve the control points. Let d ¼ lpiþ1 2 pi l and kðsl Þ be the curvature at a sampled point vðsl Þ on the spline segment. Also, let Td and Tk be the average Euclidean distance
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
555
Algorithm 2. : Pseudo code of NURBS snakes. procedure NURBS_snakes initialize vðsÞ repeat Ð { minimise s EðvðsÞÞds apply balloon forces update locations of {pi }n21 do weight_adjustment construct NURBS curve based on {pi }n21 }until (snake converge)
Fig. 2. The block diagram of the NURBS snake algorithm.
between the control points and the average curvature of the NURBS contour, respectively. These values are used as the thresholding limits for the addition (or deletion) of the control points during the iteration process. The implementation of weight adjustments is given in Algorithm 1 while the pseudo code for the NURBS snake is in Algorithm 2. At each iteration, the energy function is calculated among the possible shifts of the control point within a neighbourhood of N pixels. Given a control point pi ; it can only be moved within a specific neighbourhood around it; in our algorithm, we consider a 3 £ 3 neighbourhood ðN ¼ 9Þ: Each control point is located in this window corresponding to the minimum energy configuration of the NURBS snake [24].
Algorithm 1. : Control point addition/deletion and weight adjustment. procedure weight_adjustment for i ¼ 0 to n þ k{ begin spline segment vðsÞ; s [ ½si ; siþpþ1 Þ { K ¼ maxl {vðsl Þ} if (d $ Td and K . Tk ) add new control point modify weight wi elseif (d $ Td and K . Tk ) piþ1 ! pi } end spline segment } endfor
After investigating the converging behaviour of the snake, we select a convergence criteria based on the varying characteristics of the total energy of the snake for a userspecified number of iteration steps. At the end of each iteration step, we check for lEðvðt þ 1ÞÞ 2 EðvðtÞÞl # e ; t # t0 where e is a given tolerance and t0 is the user-specified maximum number of iteration steps. If this condition is satisfied, the snake is considered to have reached equilibrium and hence, the iteration process is stopped.
6. Experiments Experiments with both synthetic images and real images obtained from MR, satellite, clinical and aerial sources are utilised to demonstrate the performance of NURBS snakes. The results are compared with previous snake models, namely, discrete snakes, T-snakes and B-spline snakes. The accuracy of the NURBS-based and other snake algorithm in extracting object boundaries is validated using the average distance error [25]. It is the average distance between each point on the deformed snake v and the closest point on the actual boundary vˆ. The average distance error e between v and vˆ is given by 1 ð eðv; v^ Þ ¼ min lvðsÞ 2 v^ ðs0 Þlds L s s0 [v^ where L is the length of v. The actual boundary was extracted by manually outlining the outer perimeter of the object of interest. The initial contour of the snake was manually placed close to the actual boundary and for all the snake models, the same initial contour was used for a given object in the image. 6.1. Synthetic images A synthetic image of size 256 £ 256 pixels was constructed using a brain mask as the foreground. The extraction process of the boundary of the brain was nontrivial as the mask had both convex and concave parts. Gaussian noise was added to have a signal-to-noise ratio (SNR) of 10 dB. The final boundaries were extracted using various snake models on the image and depicted in Fig. 3.
556
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
Fig. 3. Illustrations of snake models converging to the boundary of the foreground object in the synthetic image having an SNR ratio of 10 dB: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
As seen in this figure, the NURBS snake model gave the nearest actual boundary of the object compared with the other three. The model parameters used; g ¼ 0:02; d ¼ 0:42 and h ¼ 0:10: The evolutions of the four snake models are depicted in Fig. 4 which shows the change in number of control points in each snake model. Both B-spline and NURBS snakes added more control points to the curve. However, the NURBS snake had lesser number of control points at convergence than the former because B-spline snakes required more points to be added near high curvature regions to increase flexibility. The important feature in this situation is that the NURBS snake, even with lesser number of control points than B-spline snake, arrived at the energy minimising state (, 55 iterations) before the latter (, 65 iterations). The synthetic image is also used to analyse the effect of g on the average distance error e for the NURBS snake as depicted in Fig. 5. This weighting parameter, that governs the external energy term, made a significant impact on the accuracy of the segmentation due to large values of g disabling the snake being attracted to strong edges and making the snake wander away from the foreground object. Hence, a small positive value for g (< 0.02) is selected to increase the external energy near the strong edges.
As the direction of the balloon force can be changed by the sign of d (see Eq. 3), negative values of d inflates the snake so much that it moves away from the region of interest which results in a large distance error as seen in Fig. 6. Also,
Fig. 4. Behaviour of different snake models on the synthetic image with respect to the number of control points. The image had an SNR ratio of 10 dB.
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
557
Fig. 5. Effect of g on the performance of NURBS snakes for the synthetic image.
Fig. 7. Effect of parameter h on the performance of NURBS snakes for the synthetic image with SNR ¼ 10 dB.
a large value of þ d pulls the snake over the required boundary with a strong deflation force; it does not stop at the strong edge between the background and the brain mask. We explored the effect of h on the performance of the NURBS snake when applied to the synthetic image as shown in Fig. 7. It can be seen that e rapidly decreases (hence, accuracy increases) until , 40 iterations for the given h values. However, it is less sensitive to h when the iterations are continued for a long time. We also experimented the robustness of the snake models on the synthetic image with different SNR ratios where the values obtained for e are plotted against the SNR ratios for the four snake models as in Fig. 8. As observed from the graph, the performance indices exhibit approximately a linear
relationship with SNR values. The NURBS snake gave the lowest value for e in all the cases which confirms the superior performance of the proposed model even under low SNR ratios compared to other snake models. The parameter d is selected in proportion to the amount of noise in the image such that d ¼ 0:5e20:01SNR :
A coronal slice from a MR head scan was used to demonstrate the capabilities of different snake models in extracting a real object with complex boundaries: the brain contour. Fig. 9 illustrates the results of applying the snake algorithms to a representative MR head scan. Although
Fig. 6. Effect of parameters d on the performance of NURBS snakes for the synthetic image with SNR ¼ 10 dB.
Fig. 8. Performance of the four snake models at different SNR ratios for the synthetic image.
6.2. MR head scans
558
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
Fig. 9. Snake models converging to the brain surface on MR head scans: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
the B-spline snakes managed to converge to the required boundary satisfactorily, the proposed NURBS snake performed best by penetrating into deep sulci while avoiding inherent noise present in MR head scans. In this application, the model parameters: g ¼ 0:02; d ¼ 0:13 and h ¼ 0:10: It is evident from Table 1 that the NURBS snake model extracted the contour that is closest to the actual boundary with a greater accuracy than the other models. 6.3. Solar images The performance of the NURBS snakes was tested on representative solar images to segment sunspots as in Fig. 10. Such images are useful in testing the proposed algorithm because these spots have complex shapes coupled with noise. These satellite images are of vital significance in astronomical studies since the area of a sunspot gives an estimate of the changes in the magnetic field strength on the sun [26]. The values obtained for e in this application are listed in Table 1. In this case too, the NURBS snake managed to converge to the nearest actual boundary of the sunspot. The model parameters: g ¼ 0:02; d ¼ 0:21 and h ¼ 0:10:
6.4. Surgical video images The flexibility of the proposed NURBS snakes was further tested on surgical video images, as in Fig. 11, to detect the area of a melanoma; a skin lesion that arises in the cells producing pigments [27]. Measuring the physical dimensions of such a disease is important for a pathologist to track the spread of the infection and record the progress of healing. As seen in Table 1, the NURBS snake captured the boundary of the skin lesion with a higher precision than the other snake models which demonstrates the capability of the NURBS snake to extract a region of interest with smooth Table 1 Comparison of average distance error e for different snake models applied on real images Snake models
Discrete Topological adaptive B-spline NURBS
e values MR images
Satellite images
Video images
Aerial images
4.9 3.1
2.2 2.9
2.9 1.9
3.8 2.6
1.9 0.7
0.8 0.4
0.9 0.2
1.5 0.6
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
559
Fig. 10. Snake models converging to a sunspots in solar images: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
boundaries. The model parameters used were: g ¼ 0:02; d ¼ 0:22 and h ¼ 0:10: 6.5. Aerial images The performance of the snakes was tested on an aerial image of lake, as in Fig. 12 (Aerial Graphics, Michigan, US, http://www.aerialgraphics.com), because the region of interest consists of both simple and high curvature regions surrounded by random textures. In this application, the surface area measurements of lakes are considered vital for water resource management in agricultural land. It is apparent from Table 1 that the proposed NURBS-based active contour give the lowest average distance error in capturing the outer boundary of the lake. The model parameters: g ¼ 0:02; d ¼ 0:09 and h ¼ 0:10:
7. Discussion Due to the additional weighting parameter, the local shape of the NURBS-based snake can be altered without moving or duplicating the relevant control points. This could be considered as the major advantage of the proposed
NURBS snake over the earlier B-spline model; the latter needs duplication of control points to increase local flexibility. In other words, the B-spline formulation, although is a compact representation for curves, requires more control points, even a particular control point to be defined more than once, to attract the curve towards high curvature regions. Due to the spline basis functions representation, this weight-based flexibility of NURBS does not affect other parts of the curve while preserving the smoothness constraints globally. Furthermore, the modification of weights does not involve user-intervention because the curvature is automatically computed from the evolving curve itself. The novelty in our approach to increase local flexibility is attributed to the fact that the proposed model integrates automatic curvature-based shape control by weight modification; a feature unique for NURBS. For the synthetic images in our experiments, all the snakes are manually initialised with 50 control points. The number of control points is maintained at 50 for discrete snakes because it is incapable of multiplying the control points during evolution. One may argue that discrete snakes could have achieved the same performance as that of NURBS snakes if the number of control points of the former
560
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
Fig. 11. Snake models converging to a skin lesion in surgical video images: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
is increased. However, the major disadvantage of discrete snakes is that all the control points must manually be placed around the region of interest which is tedious, consumes time and subjective if the number of required points is very large. In contrast, NURBS snakes manage to increase the number of control points automatically although the initial number of control points, that are manually placed, is small. The Euclidean distance between two control points alone is not a good criteria to consider when dynamically adding (or removing) control points because the two points could be situated at a large distance apart at a low curvature region of the object of interest. Hence, we used both the Euclidean distance d and the curvature kðsÞ to increase (or decrease) the number of control points (notice in Algorithm 1 that a control point is added if both d and K exceed the relevant thresholding limits). In B-spline snakes, however, a new control point must be added even if d , Td to increase local flexibility which leads to increased computations. The estimation of model parameters a and b in Eq. (1) has been an open problem in active contour research due to difficulties in choosing the optimal configuration for a specific application. Even in most recent literature on snakes [28], these parameters are selected empirically to optimise
the convergence. Although some researchers have spelled several techniques to determine these values automatically [7,18], they are very much application specific. A cubic NURBS curve, on the other hand, allows us to overcome the problem of selecting optimal values for these smoothing parameters because cubic splines are optimal for cost functions involving smoothness constraints based on the first and second order derivatives [29]. Hence, the internal energy function of the traditional active contour model can be dropped by using a cubic NURBS curve. For the synthetic image, the parameter d is estimated based on its SNR ratio as given in Section 6.1 Experiments on the sensitivity of this parameter verify that the NURBS model produced satisfactory results when d is within , 0.1– , 0.4. Hence, we empirically selected a value within this range for the real images. A question may also arise as to why we have considered T-snakes for comparison purposes when all the objects in our experiments constitute the same topology. We must emphasize that our research deals more with extracting object boundaries accurately than topological adaptability. Therefore, we compared the accuracy of T-snakes with the proposed NURBS snakes within the context of extracting
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
561
Fig. 12. Snake models converging to the boundary of a lake in an aerial image: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
the required boundaries in objects of interest. Further research studies could be carried out to verify the adaptability of NURBS snakes on objects with varied topologies. The capability of NURBS snakes was demonstrated on a synthetic image and four types of real images. Also, its robustness to noise was demonstrated on the synthetic image with varying SNR ratios. As for the real images (MR head scans, satellite, video and aerial images), the snakes were applied to show the performance of NURBS snakes at different image features. The NURBS snake contains a balloon force that provides the curve a facility to overcome noise in the image. Notice in Eq. (3) that the balloon energy is dependent on the maximum gradient term along the curve; if k7IðvðsÞÞk ¼ 0; it is unlikely that noise pixel is encountered along the curve vðsÞ and hence, the balloon energy could be inhibited. Several aspects in our technique may be further explored in future. A method to select the optimal knot estimation during the addition (or deletion) of the control points needs to be investigated. When a control point is added, a new value has to be computed and inserted to the remaining knot vector. The difficulty in this scenario is that the shape of the curve depends on both the weight adjustment and spacing of
the knots which in turn effect the accuracy of the extracted boundary.
8. Conclusion We have revealed a new formulation for deformable contours based on NURBS curves. Though active contours are useful in many applications, the earlier models are unable to extract complex shapes accurately due to their inability to change the shape locally at high curvature regions. Compared to other active contour models, the NURBS snakes inherits a unique property that facilitate an accurate extraction of objects with complex boundaries; weight-based shape control. In contrast to the B-spline model, each control point of a NURBS curve is associated with a weighting parameter which makes the snake a combination of rational spline basis functions. The weights are made to change dynamically according to the curvature properties of the evolving curve so that the resulting flexibility influences the local shape of the contour to be attracted towards complex shapes without any user-intervention and duplication of control points.
562
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
The rational approach presented in this paper is for the 2D implementation of snakes. However, extension of the proposed model to 3D is straightforward. Hence, the proposed algorithm provides the initial theoretical framework for flexible active contours which can be extended to 3D NURBS-based curvature dependent deformable surfaces.
References [1] M. Kaas, A. Witkin, D. Terzopoulos, Snakes: active contour models, International Journal of Computer Vision 1 (4) (1988) 321–331. [2] S. Atkins, B. Mackiewich, Fully automated segmentation of the brain in MRI, IEEE Transactions on Medical Imaging 17 (1) (1998) 99–107. [3] C. Davatzikos, J. Prince, An active contour model for mapping the cortex, IEEE Transactions on Medical Imaging 14 (1) (1995) 65– 80. [4] D. Xu, J. Hwana, A topology independent active contour tracking, Proceedings of IEEE Signal Processing Society Workshop on Neural Networks for Signal Processing, 1999, pp. 429–438. [5] T. Jones, P. Plassmann, An active contour model for measuring the leg ulcers, IEEE Transactions on Medical Imaging 19 (2) (2000) 1202–1210. [6] K. Lam, H. Yam, Locating head boundaries by snakes, International Symposium on Speech, Image Processing and Neural Networks (1994) 17–20. [7] S. Fenster, J. Kender, Sectored snakes: evaluating learned-energy segmentation, Proceedings of the International Conference on Computer Vision, 1998, pp. 420 –426. [8] A. Yezzi, S. Kichenessamy, P. Olver, A. Tannenbaum, A geometric snake model for segmentation of medical imagery, IEEE Transactions on Medical Imaging 16 (2) (1997) 199–209. [9] A.K.J. Gao, A. Kaka, A deformable model for human organ extraction, Proceedings of the International Conference on Image Processing, 3, 1998, pp. 323–327. [10] Y. Matsuzawa, T. Abe, Region extraction using competition of multiple active contour models, Proceedings of the International Conference on Image Processing, vol. 3, 1999, pp. 198–202, pp. 198–202. [11] D. Kang, Stable snakes algorithm for convex tracking of MRI sequences, IEEE Electronic Lettters 35 (1) (1999) 1070–1071. [12] K. Lai, R. Chin, Deformable contour modeling and extraction, IEEE Transactions on Pattern Analysis and Machine Intelligence 17 (11) (1995) 1084–1090.
[13] L. Ji, H. Yan, An attractable snake model for contour extraction in MRI images, International Conference on IEEE Engineering in Medicine and Biology Society 20 (2) (1998) 609– 612. [14] W. Neuenschwander, P. Fua, G. Szekely, O. Kubler, Making snakes coverage from minimal initialization, Proceedings of the International Conference on Pattern Recognition, vol. 1, 1994, pp. 613 –615. [15] T. McInerny, D. Terzopoulos, T-snakes: Topology adaptive snakes, Proceedings of the International Conference on Computer Vision, 1995, pp. 840 –845. [16] S. Menet, P. Marc, G. Medioni, Implementation and application to stereo, Image Understanding Workshop, 1990, pp. 720 –726. [17] H. Qin, D. Terzopoulos, D-NURBS: a physics-based framework for geometric design, IEEE Transactions on Visualization and Computer Graphics vol. 2 (1) (1996) 85–95. [18] F. Labroose, P. Willis, Towards continuous image representations, Proceedings of the International Conference in Central Europe on Computer Graphics, Visualization and Computer Vision (WSCG), 2001, pp. 206 –213. [19] L. Piegl, W. Tiller, The NURBS Book, Springer, Berlin, 1995. [20] L. Cohen, On active contour models and balloons, Computer Vision, Graphics and Image Processing: Image Understanding 53 (2) (1991) 211 –218. [21] T. McInerny, D. Terzopoulos, Deformable models in medical image analysis, Proceedings of the Workshop on Mathematical Methods in Biomedical Image Analysis, 1996, pp. 171– 180. [22] X. Li, J. Wang, Adaptive balloon models, International Conference on Computer Vision and Pattern Recognition 2 (1999) 434–439. [23] L. Cohen, I. Cohen, Finite element methods for active contour models and balloons for 2D and 3D images, IEEE Transactions on PAMI 15 (1) (1993) 1131–1147. [24] A. Klein, F. Lee, A. Amini, Quantitative coronary angiography with deformable spline models, IEEE Transactions on Medical Imaging 16 (5) (1997) 468–482. [25] L. Staib, J. Duncan, Model-based deformable surface finding for medical images, IEEE Transactions on Medical Imaging 15 (6) (1996) 859 –870. [26] E. Priest, Solar-Magneto Hydrodynamics, Dordrecht, Holland, 1984. [27] M. Brygel, Video Book of Surgery, North-East Valley Division of General Practice, 2000. [28] X. Pardo, M. Carreira, A. Mosquera, D. Cabello, A snake for CT image segmentation intergrating region and edge information, Image and Vision Computing 19 (1) (2001) 461–475. [29] P. Brigger, J. Hoeg, M. Unser, B-spline snakes: A flexible tool for parametric contour detection, IEEE Transactions on Image Processing 9 (9) (2000) 1484–1496.
NURBS snakes Ravinda G.N. Meegama, Jagath C. Rajapakse* School of Computer Engineering, Nanyang Technological University, Nanyang Avenue 639798, Singapore Received 12 November 2002; received in revised form 7 March 2003; accepted 18 March 2003
Abstract A new formulation of energy minimising splines, based on Non-Uniform Rational B-Spline (NURBS) curves, is presented. The existing snake models do not facilitate the local control of their shapes without increasing the number of control points. The introduction of NURBS renders more flexibility into the snake model due to the weighting parameters that are made to adapt according to the curvature of the contour. The proposed snake model was tested on both synthetic and real images; the new model demonstrated more local flexibility than the previous models and is capable of attracting to complex boundaries of the objects with increased accuracy. q 2003 Elsevier Science B.V. All rights reserved. Keywords: Active contours; Non-Uniform Rational B-Spline; Energy minimisation; Boundary extraction; Snakes
1. Introduction Extracting object boundaries accurately from images has become an increasingly important and a challenging task in the arena of image analysis research over the past two decades. This task is extremely difficult when the shape of the object exhibits complex features such as high curvature regions. To address this issue, Kaas, Witkin and Terzopoulos introduced a model-based active shape contour referred to as a snake [1]. A snake is initialised approximately around the object of interest by manually placing a set of discrete points, referred to as control points, that govern the movement of the contour. These points evolve in time according to internal forces (elastic forces) determined by the shape of the curve and external forces (image forces) obtained by local intensity statistics of the image. Since their inception, the behaviour of snakes near a highly convoluted and a complex boundary, such as the surface of the human brain, has been of highly interest [2 – 4]. Snakes have also been used extensively for segmentation of medical images such as measuring the area of leg ulcers [5], locating the boundaries of the head [6], bladder [7], heart [8] and liver [9]. A majority of these research focused on generalising the snake model and improving the convergence and stability of its evolution * Corresponding author. Tel.: þ 65-67905802; fax: þ 65-67926559. E-mail addresses: [email protected] (J.C. Rajapakse); [email protected] (R.G. Meegama).
[10 –14]. Several new snake models have emerged as a result of this, for example, topological adaptive snakes [15] and B-spline snakes [16]. However, these snake models are incapable of increasing the local flexibility of the contour without adding more control points. To attract a snake onto a boundary having a high curvature, the changes made to the positions of the control points near a particular segment of the snake should affect the shape of only the immediate neighbouring segments of the snake. The local control of an active contour is vital in cases where the shape to be extracted exhibits both high and low curvature regions. The discrete snakes [1] cannot adapt to complex shapes of object boundaries because the number of control points remains constant during the evolution. Although a discrete snake is simple to implement, the accuracy of the extracted boundary is limited because the snake is handicapped by the fixed number of control points. A deformation scheme that extends the geometrical flexibility of the snake to deal with complex shapes, referred to as topologically adaptive snakes (T-snakes), is revealed in [15]. Both discrete and T-snake suffer from the drawback that the movement of a single control point affects the entire length of the contour; a small change in the position of a control point tends to strongly propagate throughout the entire snake. This can be described as a global propagation of change of a particular control point. The B-spline snake [16], in which the contour is constructed from a parametric B-spline curve, has improved
0262-8856/03/$ - see front matter q 2003 Elsevier Science B.V. All rights reserved. doi:10.1016/S0262-8856(03)00066-0
552
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
local control of the contour. The curve is approximated by piecewise polynomial functions that define specific segments on the curve. Although B-spline snakes perform much better than the former models with respect to accuracy, the distinct control points may have to be duplicated to detect high curvature regions, resulting in increased redundant computations. NURBS-based deformable models, called D-NURBS (or Dynamic NURBS), have been investigated in [17] for computer aided geometric design. In this technique, the weights near a desired region have to be interactively adjusted by the designer in order to influence the shape of the curve. A NURBS-based active contour is outlined in [18] to segment objects in blurred digital images. This method, however, does not make use of the full potential of the NURBS representation of curves, such as shape modification with respect to the weights, to extract the complex boundaries. The main goal of our work is to address the problems encountered in the above methods by introducing a new snake model that is equipped with an additional parameter that can be automatically adjusted to change the local shape of the curve. Our proposed snake model, referred to as a NURBS snake, uses the NURBS representation of curves [19]; the weighting parameters in NURBS are automatically adjusted to control the flexibility of the curve at each control point. During the evolution, the weight associated with each distinct control point is modified depending on the curvature of the points on the relevant spline segment. Such a modification of weights influences the snake to attract towards high curvature regions without adding more control points. The paper is structured as follows: Section 2 presents major definitions in the traditional snake model. Section 3 outlines the principles of basis functions and B-splines. The NURBS snakes model is introduced in Section 4 followed by the implementation details in Section 5. Experimental details are presented in Section 6 and finally, Section 7 discusses in detail the performance of the NURBS model.
2. Snakes definitions This section gives the important definitions pertaining to active contour models and the associated energy terms. In this paper, we assume two-dimensional (2D) images.
[0,1] £ [0,1) ! V that parameterises the image coordinates ðx; yÞ by the arc parameter s and a time parameter t [ ½0; 1 such that vðs; tÞ ¼ ðxðs; tÞ; yðs; tÞÞT In what follows, in a curve vðs; tÞ; we sometimes omit either the time parameter t or the arc parameter s; or both, to simplify the notations according to the context. Definition 3. Internal Energy. The internal energy Eint at a point vðsÞ on a curve v is defined as Eint ðvðsÞÞ ¼ alv0 ðsÞl2 þ blv00 ðsÞl2
ð1Þ
where a; b [ {0; R }; v ðsÞ ¼ ›vðsÞ=›s and v ðsÞ ¼ › vðsÞ= ›s2 : þ
0
00
2
Definition 4. External Energy. The external energy Eext at a point vðsÞ on a curve v is defined as Eext ðvðsÞÞ ¼ 2gl7IðvðsÞÞl2
ð2Þ
where g [ Rþ ; IðvðsÞÞ denotes the intensity at the point vðsÞ and l·l denotes the Euclidean norm. Definition 5. Balloon Force. The balloon force applied on a point vðsÞ is defined as [20] ! l7IðvðsÞÞl2 bðvðsÞÞ ¼ d 1 2 nðvðsÞÞ ð3Þ k7IðvðsÞÞk2 where d [ R; nðvðsÞÞ is the unit normal vector at vðsÞ and k·k denotes the maximum norm of a vector such that k7IðvðsÞÞk ¼ max l7IðvðsÞÞl: s
Definition 6. Total Energy. The total energy E of a curve point vðsÞ is defined as EðvðsÞÞ ¼ Eint ðvðsÞÞ þ Eext ðvðsÞÞ:
ð4Þ
Definition 7. Snake. A snake is a dynamic curve that evolves over time attempting to minimise the total energy of its contour. If vðs; tÞ represents a point on a dynamic curve at a discrete time instance t; the point at time t þ 1; vðs; t þ 1Þ; is given by vðs; t þ 1Þ ¼ arg min Eðvðs; tÞÞ þ bðvðs; tÞÞ vðs;tÞ[n
where n is the set of all possible curve points vðs; tÞ: Definition 1. Curve. If V denotes a 2D image domain, a curve is a one-parameter vector-valued mathematical function v: [0,1] ! V such that an arbitrary point on the curve v is parametrically expressed as vðsÞ ¼ ðxðsÞ; yðsÞÞT ; s [ ½0; 1: Definition 2. Dynamic Curve. A dynamic curve is a twoparameter vector-valued mathematical function v:
The energy of the snake is defined in such a way that starting from a user-specified configuration near an object boundary, it deforms itself to conform to the nearest salient boundary by minimising the combined internal and external energies. The internal energy of the snake, Eint ; models the elasticity of the contour and makes it behave as a membrane and/or thin plate in which the amount of stretching and
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
flexibility is characterised by the parameters a and b; respectively. The external energy, Eext ; is determined by the intensity statistics of the image and leads the snake toward the required edges of the image. The balloon force, b, is applied to an energy minimised point on the contour and it forces the snake to overcome noise and other spurious features on the image [21 – 23].
3. Basis function representation A parametric curve vðsÞ can be defined as a linear combination of a given set of points, referred to as control points, and a set of polynomial functions fi : ½0; 1 ! R; i ¼ 0; 1; …; n 2 1; such that vðsÞ ¼
nX 21
fi ðsÞpi
ð5Þ
i¼0
When the polynomial functions form a basis for the vector space of all the polynomials of the desired degree and continuity, the polynomials are referred to as basis functions. The basis function fi can be a piecewise polynomial function defined on an interval, say ½a; b; with specified continuity constraints such that when the interval is subdivided by a sequence of real values, si where a ¼ s0 , … , sn21 ¼ b; fi is a polynomial within each subinterval ½sj ; sjþ1 Þ: In such a case, fi is referred to as a spline basis function. Spline basis functions provide local control over the shape of a curve if each basis function, fi ; is made non-zero only within a limited number of subintervals of s: As each pi is multiplied by the corresponding basis function fi (Eq. 5), a movement of the control point pi affects only the segment of the curve lying in the subinterval where fi ðsÞ – 0: Therefore, such a spline basis function representation can be used to model snakes having local flexibility. Definition 8. B-spline Curve. A B-spline curve v : ½0; 1 ! V of degree k [ Zþ is defined by vðsÞ ¼
nX 21
Ni;k ðsÞpi
knot, si ; relates the parameter s to each control point pi such that a translation of pi affects only the curve segment vðsÞ; s [ ½si ; siþkþ1 Þ: This local control comes from the Bspline basis functions in which Ni;k ðsÞ ¼ 0 for s ½si ; siþkþ1 Þ: B-spline P basis functions hold the property of partition of unity; i Ni;k ðsÞ ¼ 1; making vðsÞ an affine combination of its control points. Moreover, since ;i ¼ 0; 1; …; n 2 1; Ni;k ðsÞ $ 0; the B-spline curve vðsÞ is a convex combination of its control points, which makes a snake formulated using a B-spline curve always lie inside the convex hull formed by the control points. The points on a uniform B-spline curve have equally spaced parameterisation. Conversely, the parameterisation of points on a ‘Non-Uniform’ B-spline curve may not necessarily be equally spaced. A non-rational B-spline function provides only an approximation to conics whereas a ‘Rational’ Bspline consists of an extra parameter, wi [ Rþ ; called a weight associated with each control point pi that gives not only an extra degree of freedom for shaping the curve but also defines precisely any of the conics. Incidently, NURBS is the acronym for Non-Uniform Rational Bsplines [19].
4. NURBS snake This section gives the definitions related to the proposed NURBS snakes model. Definition 9. Dynamic NURBS Curve. If pti and wti denote the ith control point and its weight at time t; respectively, a dynamic NURBS curve v of degree k is defined as Xn21 t wi Ni;k ðsÞpti vðs; tÞ ¼ Xi¼0 n21 t wi Ni;k ðsÞ i¼0 It can also be written as vðs; tÞ ¼
i¼0
where the ith B-spline basis function of degree k; Ni;k ; is recursively defined as ( 1; if s [ ½si ; siþ1 Þ Ni;0 ðsÞ ¼ 0; otherwise
553
n21 X
Rti ðsÞpti
ð6Þ
i¼0
where Rti ðsÞ ¼
wti Ni;k ðsÞ n21 X t wj Nj;k ðsÞ j¼0
and ðs 2 si ÞNi;k21 ðsÞ ðs 2 sÞNiþ1;k21 ðsÞ þ iþkþ1 Ni;k ðsÞ ¼ siþk 2 si siþkþ1 2 siþ1 where si ; i ¼ 0; 1; …; n þ k; are positive real numbers referred to as knots. The vector (s0 ; s1 ; …; snþk ÞT is called the knot vector of the B-spline curve v. The value of the ith
The NURBS curve in Eq. (6) can be expressed in matrix from as vðsÞ ¼ pT RðsÞ where p ¼ ðpt0 ; pt1 ; …; pTn21 ÞT and RðsÞ ¼ ðR0 ðsÞ; R1ðsÞ ; …; Rn21 ðsÞÞT :The internal energy, external energy and
554
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
the balloon force for a NURBS curve, computed using the definitions given in Section 2, are given, respectively, by Eint ðvðsÞÞ ¼ alpT R0 ðsÞl2 þ blpT R00 ðsÞl; 2 Eext ðvðsÞÞÞ ¼ 2g7IðpT RðsÞÞ and bðvðsÞÞ ¼ d 1 2
l7IðpT RðsÞÞl2 k7IðpT RðsÞÞk2
ð7Þ
! nðpT RðsÞÞ:
ð8Þ
If pi ¼ ðxi ; yi ÞT ; the gradient operation 7I is defined as s ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 › › T 7Iðp RðsÞÞ ¼ wpI þ wpI ›xT RðsÞ ›yT RðsÞ where x ¼ ðx0 ; x1 ; …; xn21 ÞT ; y ¼ ðy0 ; y1 ; …; yn21 ÞT ; w is a smoothing kernel and p denotes the convolution operator. Definition 10. NURBS Snake. A NURBS snake is a dynamic NURBS curve vðs; tÞ whose energy is minimised with respect to the points on the curve such that a curve point at time t þ 1; vðs; t þ 1Þ is given by vðs; t þ 1Þ ¼ arg min ðEint ðvðs; tÞÞ þ Eext ðvðs; tÞÞÞ þ bðvðs; tÞÞ vðs;tÞ
Fig. 1. Movement of a curve point vðsÞ towards the desired control point when the associated weight is increased.
of the control point ptþ1 is updated such that weight wtþ1 i i wtþ1 ¼ wti þ h i
Kt kkðs; tÞk
ð11Þ
where the parameter h [ R controls the amount of attraction of the curve towards the control point and kkðs; tÞk is the maximum curvature of the closed NURBS curve. Note that in Eq. (11), the curvature Kt is divided by the maximum curvature in order to normalise the weight modification.
4.1. Weight-based control 5. Implementation This section presents how the weights of a NURBS snake, without moving or duplicating the relevant control points, effectively control the local shape of the contour based on its curvature properties. Let vðs; wi Þ be the family of curves obtained by changing the weight wi at a control point pi : Also, let vðs; wi ¼ 0Þ be the set of curves when wi is zero. Then, v can be expressed as [19] vðs; wi Þ ¼ vðs; wi ¼ 0Þ þ tðpi 2 vðs; wi ¼ 0ÞÞ
ð9Þ
where
t¼
Ni;k ðsÞwi nX 21
ð10Þ
Nj;k ðsÞwj
j¼0
From Eqs. (9) and (10), it is clear that as the weight wi associated with the control point pi is increased, t is increased, and hence, the point vðs; wi Þ moves toward pi along the vector pi 2 vðs; wi ¼ 0Þ as shown in Fig. 1. Similarly, vðs; wi Þ moves away from pi if t is decreased. Let pti denote the position of a control point with weight t wi at time t and kðvðs; tÞÞ; the curvature at a point vðs; tÞ: If the maximum curvature within the spline segment vðsÞ; s [ ½si ; siþpþ1 Þ; at time t; Kt ; exceeds a pre-defined value (say, the average curvature along the NURBS curve), the new
5.1. Evolution of points In the beginning, the control points of the snake are selected near the object boundaries to ensure that the snake does not converge to a local minimum. The selection of these initial points is based on a priori information or an interactive manual placement. These points can lie inside (or outside) the object; the balloon force is able to push (or pull) the snake accordingly. At a given iteration, new control points are dynamically added to (or removed from) the control polygon based on both the Euclidean distance between two consecutive control points and the curvature properties within the relevant spline segment. The removal of control points is due to the fact that only a fewer points are needed for the snake to expand at low curvature regions (in B-spline snakes, however, new points must be added at high curvature regions). The positions of the control points are now moved according to the minimal energy configuration of the spline and a NURBS curve is constructed based on these new locations of the control points. The major components of our algorithm is illustrated in Fig. 2. We used dynamic programming, as highlighted in [24], to evolve the control points. Let d ¼ lpiþ1 2 pi l and kðsl Þ be the curvature at a sampled point vðsl Þ on the spline segment. Also, let Td and Tk be the average Euclidean distance
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
555
Algorithm 2. : Pseudo code of NURBS snakes. procedure NURBS_snakes initialize vðsÞ repeat Ð { minimise s EðvðsÞÞds apply balloon forces update locations of {pi }n21 do weight_adjustment construct NURBS curve based on {pi }n21 }until (snake converge)
Fig. 2. The block diagram of the NURBS snake algorithm.
between the control points and the average curvature of the NURBS contour, respectively. These values are used as the thresholding limits for the addition (or deletion) of the control points during the iteration process. The implementation of weight adjustments is given in Algorithm 1 while the pseudo code for the NURBS snake is in Algorithm 2. At each iteration, the energy function is calculated among the possible shifts of the control point within a neighbourhood of N pixels. Given a control point pi ; it can only be moved within a specific neighbourhood around it; in our algorithm, we consider a 3 £ 3 neighbourhood ðN ¼ 9Þ: Each control point is located in this window corresponding to the minimum energy configuration of the NURBS snake [24].
Algorithm 1. : Control point addition/deletion and weight adjustment. procedure weight_adjustment for i ¼ 0 to n þ k{ begin spline segment vðsÞ; s [ ½si ; siþpþ1 Þ { K ¼ maxl {vðsl Þ} if (d $ Td and K . Tk ) add new control point modify weight wi elseif (d $ Td and K . Tk ) piþ1 ! pi } end spline segment } endfor
After investigating the converging behaviour of the snake, we select a convergence criteria based on the varying characteristics of the total energy of the snake for a userspecified number of iteration steps. At the end of each iteration step, we check for lEðvðt þ 1ÞÞ 2 EðvðtÞÞl # e ; t # t0 where e is a given tolerance and t0 is the user-specified maximum number of iteration steps. If this condition is satisfied, the snake is considered to have reached equilibrium and hence, the iteration process is stopped.
6. Experiments Experiments with both synthetic images and real images obtained from MR, satellite, clinical and aerial sources are utilised to demonstrate the performance of NURBS snakes. The results are compared with previous snake models, namely, discrete snakes, T-snakes and B-spline snakes. The accuracy of the NURBS-based and other snake algorithm in extracting object boundaries is validated using the average distance error [25]. It is the average distance between each point on the deformed snake v and the closest point on the actual boundary vˆ. The average distance error e between v and vˆ is given by 1 ð eðv; v^ Þ ¼ min lvðsÞ 2 v^ ðs0 Þlds L s s0 [v^ where L is the length of v. The actual boundary was extracted by manually outlining the outer perimeter of the object of interest. The initial contour of the snake was manually placed close to the actual boundary and for all the snake models, the same initial contour was used for a given object in the image. 6.1. Synthetic images A synthetic image of size 256 £ 256 pixels was constructed using a brain mask as the foreground. The extraction process of the boundary of the brain was nontrivial as the mask had both convex and concave parts. Gaussian noise was added to have a signal-to-noise ratio (SNR) of 10 dB. The final boundaries were extracted using various snake models on the image and depicted in Fig. 3.
556
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
Fig. 3. Illustrations of snake models converging to the boundary of the foreground object in the synthetic image having an SNR ratio of 10 dB: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
As seen in this figure, the NURBS snake model gave the nearest actual boundary of the object compared with the other three. The model parameters used; g ¼ 0:02; d ¼ 0:42 and h ¼ 0:10: The evolutions of the four snake models are depicted in Fig. 4 which shows the change in number of control points in each snake model. Both B-spline and NURBS snakes added more control points to the curve. However, the NURBS snake had lesser number of control points at convergence than the former because B-spline snakes required more points to be added near high curvature regions to increase flexibility. The important feature in this situation is that the NURBS snake, even with lesser number of control points than B-spline snake, arrived at the energy minimising state (, 55 iterations) before the latter (, 65 iterations). The synthetic image is also used to analyse the effect of g on the average distance error e for the NURBS snake as depicted in Fig. 5. This weighting parameter, that governs the external energy term, made a significant impact on the accuracy of the segmentation due to large values of g disabling the snake being attracted to strong edges and making the snake wander away from the foreground object. Hence, a small positive value for g (< 0.02) is selected to increase the external energy near the strong edges.
As the direction of the balloon force can be changed by the sign of d (see Eq. 3), negative values of d inflates the snake so much that it moves away from the region of interest which results in a large distance error as seen in Fig. 6. Also,
Fig. 4. Behaviour of different snake models on the synthetic image with respect to the number of control points. The image had an SNR ratio of 10 dB.
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
557
Fig. 5. Effect of g on the performance of NURBS snakes for the synthetic image.
Fig. 7. Effect of parameter h on the performance of NURBS snakes for the synthetic image with SNR ¼ 10 dB.
a large value of þ d pulls the snake over the required boundary with a strong deflation force; it does not stop at the strong edge between the background and the brain mask. We explored the effect of h on the performance of the NURBS snake when applied to the synthetic image as shown in Fig. 7. It can be seen that e rapidly decreases (hence, accuracy increases) until , 40 iterations for the given h values. However, it is less sensitive to h when the iterations are continued for a long time. We also experimented the robustness of the snake models on the synthetic image with different SNR ratios where the values obtained for e are plotted against the SNR ratios for the four snake models as in Fig. 8. As observed from the graph, the performance indices exhibit approximately a linear
relationship with SNR values. The NURBS snake gave the lowest value for e in all the cases which confirms the superior performance of the proposed model even under low SNR ratios compared to other snake models. The parameter d is selected in proportion to the amount of noise in the image such that d ¼ 0:5e20:01SNR :
A coronal slice from a MR head scan was used to demonstrate the capabilities of different snake models in extracting a real object with complex boundaries: the brain contour. Fig. 9 illustrates the results of applying the snake algorithms to a representative MR head scan. Although
Fig. 6. Effect of parameters d on the performance of NURBS snakes for the synthetic image with SNR ¼ 10 dB.
Fig. 8. Performance of the four snake models at different SNR ratios for the synthetic image.
6.2. MR head scans
558
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
Fig. 9. Snake models converging to the brain surface on MR head scans: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
the B-spline snakes managed to converge to the required boundary satisfactorily, the proposed NURBS snake performed best by penetrating into deep sulci while avoiding inherent noise present in MR head scans. In this application, the model parameters: g ¼ 0:02; d ¼ 0:13 and h ¼ 0:10: It is evident from Table 1 that the NURBS snake model extracted the contour that is closest to the actual boundary with a greater accuracy than the other models. 6.3. Solar images The performance of the NURBS snakes was tested on representative solar images to segment sunspots as in Fig. 10. Such images are useful in testing the proposed algorithm because these spots have complex shapes coupled with noise. These satellite images are of vital significance in astronomical studies since the area of a sunspot gives an estimate of the changes in the magnetic field strength on the sun [26]. The values obtained for e in this application are listed in Table 1. In this case too, the NURBS snake managed to converge to the nearest actual boundary of the sunspot. The model parameters: g ¼ 0:02; d ¼ 0:21 and h ¼ 0:10:
6.4. Surgical video images The flexibility of the proposed NURBS snakes was further tested on surgical video images, as in Fig. 11, to detect the area of a melanoma; a skin lesion that arises in the cells producing pigments [27]. Measuring the physical dimensions of such a disease is important for a pathologist to track the spread of the infection and record the progress of healing. As seen in Table 1, the NURBS snake captured the boundary of the skin lesion with a higher precision than the other snake models which demonstrates the capability of the NURBS snake to extract a region of interest with smooth Table 1 Comparison of average distance error e for different snake models applied on real images Snake models
Discrete Topological adaptive B-spline NURBS
e values MR images
Satellite images
Video images
Aerial images
4.9 3.1
2.2 2.9
2.9 1.9
3.8 2.6
1.9 0.7
0.8 0.4
0.9 0.2
1.5 0.6
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
559
Fig. 10. Snake models converging to a sunspots in solar images: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
boundaries. The model parameters used were: g ¼ 0:02; d ¼ 0:22 and h ¼ 0:10: 6.5. Aerial images The performance of the snakes was tested on an aerial image of lake, as in Fig. 12 (Aerial Graphics, Michigan, US, http://www.aerialgraphics.com), because the region of interest consists of both simple and high curvature regions surrounded by random textures. In this application, the surface area measurements of lakes are considered vital for water resource management in agricultural land. It is apparent from Table 1 that the proposed NURBS-based active contour give the lowest average distance error in capturing the outer boundary of the lake. The model parameters: g ¼ 0:02; d ¼ 0:09 and h ¼ 0:10:
7. Discussion Due to the additional weighting parameter, the local shape of the NURBS-based snake can be altered without moving or duplicating the relevant control points. This could be considered as the major advantage of the proposed
NURBS snake over the earlier B-spline model; the latter needs duplication of control points to increase local flexibility. In other words, the B-spline formulation, although is a compact representation for curves, requires more control points, even a particular control point to be defined more than once, to attract the curve towards high curvature regions. Due to the spline basis functions representation, this weight-based flexibility of NURBS does not affect other parts of the curve while preserving the smoothness constraints globally. Furthermore, the modification of weights does not involve user-intervention because the curvature is automatically computed from the evolving curve itself. The novelty in our approach to increase local flexibility is attributed to the fact that the proposed model integrates automatic curvature-based shape control by weight modification; a feature unique for NURBS. For the synthetic images in our experiments, all the snakes are manually initialised with 50 control points. The number of control points is maintained at 50 for discrete snakes because it is incapable of multiplying the control points during evolution. One may argue that discrete snakes could have achieved the same performance as that of NURBS snakes if the number of control points of the former
560
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
Fig. 11. Snake models converging to a skin lesion in surgical video images: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
is increased. However, the major disadvantage of discrete snakes is that all the control points must manually be placed around the region of interest which is tedious, consumes time and subjective if the number of required points is very large. In contrast, NURBS snakes manage to increase the number of control points automatically although the initial number of control points, that are manually placed, is small. The Euclidean distance between two control points alone is not a good criteria to consider when dynamically adding (or removing) control points because the two points could be situated at a large distance apart at a low curvature region of the object of interest. Hence, we used both the Euclidean distance d and the curvature kðsÞ to increase (or decrease) the number of control points (notice in Algorithm 1 that a control point is added if both d and K exceed the relevant thresholding limits). In B-spline snakes, however, a new control point must be added even if d , Td to increase local flexibility which leads to increased computations. The estimation of model parameters a and b in Eq. (1) has been an open problem in active contour research due to difficulties in choosing the optimal configuration for a specific application. Even in most recent literature on snakes [28], these parameters are selected empirically to optimise
the convergence. Although some researchers have spelled several techniques to determine these values automatically [7,18], they are very much application specific. A cubic NURBS curve, on the other hand, allows us to overcome the problem of selecting optimal values for these smoothing parameters because cubic splines are optimal for cost functions involving smoothness constraints based on the first and second order derivatives [29]. Hence, the internal energy function of the traditional active contour model can be dropped by using a cubic NURBS curve. For the synthetic image, the parameter d is estimated based on its SNR ratio as given in Section 6.1 Experiments on the sensitivity of this parameter verify that the NURBS model produced satisfactory results when d is within , 0.1– , 0.4. Hence, we empirically selected a value within this range for the real images. A question may also arise as to why we have considered T-snakes for comparison purposes when all the objects in our experiments constitute the same topology. We must emphasize that our research deals more with extracting object boundaries accurately than topological adaptability. Therefore, we compared the accuracy of T-snakes with the proposed NURBS snakes within the context of extracting
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
561
Fig. 12. Snake models converging to the boundary of a lake in an aerial image: (a) discrete snake, (b) topological adaptive snake, (c) B-spline snake and (d) NURBS snake.
the required boundaries in objects of interest. Further research studies could be carried out to verify the adaptability of NURBS snakes on objects with varied topologies. The capability of NURBS snakes was demonstrated on a synthetic image and four types of real images. Also, its robustness to noise was demonstrated on the synthetic image with varying SNR ratios. As for the real images (MR head scans, satellite, video and aerial images), the snakes were applied to show the performance of NURBS snakes at different image features. The NURBS snake contains a balloon force that provides the curve a facility to overcome noise in the image. Notice in Eq. (3) that the balloon energy is dependent on the maximum gradient term along the curve; if k7IðvðsÞÞk ¼ 0; it is unlikely that noise pixel is encountered along the curve vðsÞ and hence, the balloon energy could be inhibited. Several aspects in our technique may be further explored in future. A method to select the optimal knot estimation during the addition (or deletion) of the control points needs to be investigated. When a control point is added, a new value has to be computed and inserted to the remaining knot vector. The difficulty in this scenario is that the shape of the curve depends on both the weight adjustment and spacing of
the knots which in turn effect the accuracy of the extracted boundary.
8. Conclusion We have revealed a new formulation for deformable contours based on NURBS curves. Though active contours are useful in many applications, the earlier models are unable to extract complex shapes accurately due to their inability to change the shape locally at high curvature regions. Compared to other active contour models, the NURBS snakes inherits a unique property that facilitate an accurate extraction of objects with complex boundaries; weight-based shape control. In contrast to the B-spline model, each control point of a NURBS curve is associated with a weighting parameter which makes the snake a combination of rational spline basis functions. The weights are made to change dynamically according to the curvature properties of the evolving curve so that the resulting flexibility influences the local shape of the contour to be attracted towards complex shapes without any user-intervention and duplication of control points.
562
R.G.N. Meegama, J.C. Rajapakse / Image and Vision Computing 21 (2003) 551–562
The rational approach presented in this paper is for the 2D implementation of snakes. However, extension of the proposed model to 3D is straightforward. Hence, the proposed algorithm provides the initial theoretical framework for flexible active contours which can be extended to 3D NURBS-based curvature dependent deformable surfaces.
References [1] M. Kaas, A. Witkin, D. Terzopoulos, Snakes: active contour models, International Journal of Computer Vision 1 (4) (1988) 321–331. [2] S. Atkins, B. Mackiewich, Fully automated segmentation of the brain in MRI, IEEE Transactions on Medical Imaging 17 (1) (1998) 99–107. [3] C. Davatzikos, J. Prince, An active contour model for mapping the cortex, IEEE Transactions on Medical Imaging 14 (1) (1995) 65– 80. [4] D. Xu, J. Hwana, A topology independent active contour tracking, Proceedings of IEEE Signal Processing Society Workshop on Neural Networks for Signal Processing, 1999, pp. 429–438. [5] T. Jones, P. Plassmann, An active contour model for measuring the leg ulcers, IEEE Transactions on Medical Imaging 19 (2) (2000) 1202–1210. [6] K. Lam, H. Yam, Locating head boundaries by snakes, International Symposium on Speech, Image Processing and Neural Networks (1994) 17–20. [7] S. Fenster, J. Kender, Sectored snakes: evaluating learned-energy segmentation, Proceedings of the International Conference on Computer Vision, 1998, pp. 420 –426. [8] A. Yezzi, S. Kichenessamy, P. Olver, A. Tannenbaum, A geometric snake model for segmentation of medical imagery, IEEE Transactions on Medical Imaging 16 (2) (1997) 199–209. [9] A.K.J. Gao, A. Kaka, A deformable model for human organ extraction, Proceedings of the International Conference on Image Processing, 3, 1998, pp. 323–327. [10] Y. Matsuzawa, T. Abe, Region extraction using competition of multiple active contour models, Proceedings of the International Conference on Image Processing, vol. 3, 1999, pp. 198–202, pp. 198–202. [11] D. Kang, Stable snakes algorithm for convex tracking of MRI sequences, IEEE Electronic Lettters 35 (1) (1999) 1070–1071. [12] K. Lai, R. Chin, Deformable contour modeling and extraction, IEEE Transactions on Pattern Analysis and Machine Intelligence 17 (11) (1995) 1084–1090.
[13] L. Ji, H. Yan, An attractable snake model for contour extraction in MRI images, International Conference on IEEE Engineering in Medicine and Biology Society 20 (2) (1998) 609– 612. [14] W. Neuenschwander, P. Fua, G. Szekely, O. Kubler, Making snakes coverage from minimal initialization, Proceedings of the International Conference on Pattern Recognition, vol. 1, 1994, pp. 613 –615. [15] T. McInerny, D. Terzopoulos, T-snakes: Topology adaptive snakes, Proceedings of the International Conference on Computer Vision, 1995, pp. 840 –845. [16] S. Menet, P. Marc, G. Medioni, Implementation and application to stereo, Image Understanding Workshop, 1990, pp. 720 –726. [17] H. Qin, D. Terzopoulos, D-NURBS: a physics-based framework for geometric design, IEEE Transactions on Visualization and Computer Graphics vol. 2 (1) (1996) 85–95. [18] F. Labroose, P. Willis, Towards continuous image representations, Proceedings of the International Conference in Central Europe on Computer Graphics, Visualization and Computer Vision (WSCG), 2001, pp. 206 –213. [19] L. Piegl, W. Tiller, The NURBS Book, Springer, Berlin, 1995. [20] L. Cohen, On active contour models and balloons, Computer Vision, Graphics and Image Processing: Image Understanding 53 (2) (1991) 211 –218. [21] T. McInerny, D. Terzopoulos, Deformable models in medical image analysis, Proceedings of the Workshop on Mathematical Methods in Biomedical Image Analysis, 1996, pp. 171– 180. [22] X. Li, J. Wang, Adaptive balloon models, International Conference on Computer Vision and Pattern Recognition 2 (1999) 434–439. [23] L. Cohen, I. Cohen, Finite element methods for active contour models and balloons for 2D and 3D images, IEEE Transactions on PAMI 15 (1) (1993) 1131–1147. [24] A. Klein, F. Lee, A. Amini, Quantitative coronary angiography with deformable spline models, IEEE Transactions on Medical Imaging 16 (5) (1997) 468–482. [25] L. Staib, J. Duncan, Model-based deformable surface finding for medical images, IEEE Transactions on Medical Imaging 15 (6) (1996) 859 –870. [26] E. Priest, Solar-Magneto Hydrodynamics, Dordrecht, Holland, 1984. [27] M. Brygel, Video Book of Surgery, North-East Valley Division of General Practice, 2000. [28] X. Pardo, M. Carreira, A. Mosquera, D. Cabello, A snake for CT image segmentation intergrating region and edge information, Image and Vision Computing 19 (1) (2001) 461–475. [29] P. Brigger, J. Hoeg, M. Unser, B-spline snakes: A flexible tool for parametric contour detection, IEEE Transactions on Image Processing 9 (9) (2000) 1484–1496.