Localized Patch-Based Fuzzy Active Contours for Image Segmentation

This paper presents a novel fuzzy region-based active contour model for image segmentation. By incorporating local patch-energy functional along each pixel of the evolving curve into the fuzziness of the energy, we construct a patch-based energy function without the regurgitation term. Its purpose is not only to make the active contour evolve very stably without the periodical initialization during the evolution but also to reduce the effect of noise. In particular, in order to reject local minimal of the energy functional, we utilize a direct method to calculate the energy alterations instead of solving the Euler-Lagrange equation of the underlying problem. Compared with other fuzzy active contour models, experimental results on synthetic and real images show the advantages of the proposed method in terms of computational efficiency and accuracy.


Introduction
Image segmentation has been one of the most fundamental and important tasks in image analysis and computer vision. Its purpose is to divide the given image into several regions where the inner pixels are homogeneous with respect to some characteristic. Over the past decades, many approaches have been developed to improve the performance of the image segmentation algorithms. Variational formulation [1][2][3] has become one of the most effective algorithms for image segmentation because it can minimize an objective function containing terms that embedded description of its regions and boundaries. However, the segmentation procedure is still considered as an essential and difficult process, especially due to the variety and complexity of the images.
Active contour models (ACMs) [2][3][4][5][6][7][8] have been proved to be an efficient framework for image segmentation. The existing ACMs can mainly be categorized into two classes: edge-based models and region-based models. The edge-based ACMs [2,3] utilize image gradient to construct an edge stopping function to stop the contour evolution on the object boundaries. These models are highly dependent on the initial contours and easily suffer from serious boundary leakage problems in the position of weak boundaries. The regionbased ACMs [4][5][6][7][8] exploit the statistical image intensity information (e.g., intensity, color, and texture) to ensure that the energy functional achieves minima when the contours reach the object boundary. Chan-Vese (CV) model [5], which is based on Mumford-Shah (M-S) model [4], has become one of the well-known region-based ACMs. This model uses the global intensity difference to guide the contour and has succeeded in detecting the objects of which the boundaries are not necessarily defined by gradient. However, in these models, the periodical reinitialization of the level set function (LSF) causes a lot of computations and some numerical errors.
In order to improve the segmentation performance of the region-based ACMs, many multiphase level set methods [7][8][9][10] are proposed to segment images with many objects, which leads to a complicated expression of energy functional and greatly increases the computational complexity. Local region-based ACMs [11][12][13][14], which combine the region-based techniques with the benefits of local information, have been proposed to segment images with intensity inhomogeneity.
Using the local intensity information, Li et al. [11,12] proposed an efficient region-based level set method driven by a local binary fitting (LBF) energy and has achieved promising results. But the LBF model needs to perform four convolution operations at each iteration, which greatly increases the computational complexity. Brox and Cremers [13] derive a statistical interpretation of the Mumford-Shah functional on local region statistics. Thanks to the analytical expression of the smooth approximation via Gaussian convolution, the coordinate descent can be replaced by a true gradient descent. A characteristic function defined by the radius parameter [14] is used to extract the local information. Three specific energies in the model are introduced that can be inserted as internal energy measures: uniform modeling, means separation, and histogram separation energy.
Later on, the fuzzy energy-based active contour (FEAC) model is introduced by Krinidis and Chatzis [15] to reject local minima. In the model, a fast optimization algorithm is applied to minimize the fuzzy energy function instead of traditional methods solving Euler-Lagrange equation. And Pereira et al. [16] proposed global and local fuzzy energybased active contours (GL-FEAC) to deal with intensity images with inhomogeneity. Besides, the updating of average prototypes could be easily influenced by noise and outliers. Wu et al. [17] propose a novel region-based fuzzy active contour model with kernel metric by the minimization of a predefined energy function. A fuzzy global and local energy [18] is proposed and the local energy is constructed by considering both local spatial and gray level/color information.
In this paper, inspired by the FEAC model [15], we proposed a novel localized patch-energy active contour. In the present work, we make two main contributions.
(1) By incorporating local patch-energy functional along each pixel of the evolving curve into the fuzziness of the energy, we construct a patch-based energy functional without regularization term. Its significant improvement is that objects which have heterogeneous statistics can be successfully segmented with localized fuzzy patch-based energies while corresponding global fuzzy region-based energies fail. In addition, the model not only avoids the periodical initialization during the evolution but also reduces the effect of noise. (2) To reject the minimal of the energy functional, we utilize a direct method to calculate the energy alterations instead of solving the Euler-Lagrange equation of the underlying problem.

Previous Work
Let ( ) : Ω → be a given gray level image and be a close contour. The piecewise constant energy functional in the Chan-Vese model [5] is defined by where 1 , 2 , and are three positive constants and 1 and 2 are two constants that approximate the image intensity ( ) inside and outside . This energy function CV ( , 1 , 2 ) can be represented by a level set formulation, and the energy minimization problem can be converted to solving a level set evolution equation.
To solve this minimization of the energy functional, the level set method [5] is represented by the unknown curve as the zero level set of a Lipschitz function ( ), such that Thus, the energy functional CV ( 1 , 2 , ) can be reformulated in terms of the level set function ( ) as follows: where denotes the Heaviside function and (⋅) denotes the Dirac delta function defined as follows: The Euler-Lagrange equations are used to solve this minimization problem in (3) and update the level set function by the gradient descent method. It is clear that the Chan-Vese model can deal with the detection of objects whose boundaries are either smooth or not necessarily defined by gradient. They do not require image smoothing and thus cannot efficiently process the images with noise. That is to say, the model is more sensitive to noise and cannot handle objects with ill-defined boundaries.
The FEAC model [15] combining the fuzzy sets with the active contour methodology aims to segment regions of interest into a two-phase image. It is provided to the algorithm with an initial partition of the image such as an initial contour. This partition defines the curve that will iteratively evolve by a minimization process. The evolving contour is implicitly represented as the pseudo zero level set of the LSF such that we have the following definitions expressed by  ( , 1 , 2 , ) = ⋅ Length ( ) where the parameters ≥ 0, 1 , 2 > 0 are fixed parameters. The fuzzy membership functions ( ) ∈ [0, 1] and (1 − ( )) represent the membership value that pixel belongs to 1 and 2 , respectively. is the weighting exponent on each fuzzy membership and usually set to 1 or 2. 1 and 2 are average prototypes of the original image inside and outside .
Keeping ( ) fixed and minimizing the energy function ( , 1 , 2 , ) with respect to 1 and 2 , it is easy to get the equations by updating the following values 1 and 2 : Furthermore, keeping 1 and 2 and fixed and minimizing the energy ( , 1 , 2 , ) with respect to , it is easy to express variable in the following way: Specifically, for a certain pixel , we compute the fuzzy membership in this pixel using (8). Then, according to the change of Δ caused by the single change of the fuzzy membership in pixel , we decide whether this new fuzzy membership replaces the old one in this pixel or not. If the change of Δ becomes negative, then the new fuzzy membership is adopted. If not, the old one is kept. However, by updating , pixels in the background could be easily labeled as pixels belonging to the object region if their intensities are very close to the average prototype of the object region. Also, these approaches still need localized information to achieve reliable performance.
The remainder of this paper is organized as follows. Section 3 describes the proposed model, including the model description, numerical approximation, and the description steps of the proposed model. Experiments and results, including experimental results and the scale parameter of the localization radius, are discussed in Section 4. Finally concluding remarks are given in Section 5.

The Proposed Model
3.1. The Model Description. Let the local patch (a circle region) Ω be centered at location ; the spatial variable ∈ Ω is a single point and independent of spatial variable .
The local patch Ω with the radius is represented as Ω = {| − | ≤ , ∈ Ω}. We define the mask function ( , ) to describe the local patch Ω as described in Figure 1: It is clear that the mask function ( , ) will be 1 when the point is within the local patch Ω and 0 otherwise. It is noticed that the value of radius should be selected properly so as to capture enough local intensity information. In our work, we assume that the energies can be constructed of a family of localized patch-based energies at each point along the curves instead of global energies in the whole image. In order to formulate the local patch-based functional, we treat each point separately and split local neighborhoods into local interior and local exterior by the evolving curve.
In the model, the evolving contour is implicitly represented as the pseudo LSF based on the membership function in (5) similar to the FEAC model [15]. For a given pixel of image domain along the curves, the local patch-based functional by incorporating the fuzzy set is given as follows: where two constants 1 and 2 represent the intensity averages of ( ) > 0.5 and ( ) < 0.5. The local patch-based energy function is an internal measure functional to express the local adherence in the proposed model at each point along the evolving contour. The pseudo LSF ( ) represent the degree of the pixel ( ) belonging to the interior region. The overall localized patch-based functional can be formulated as follows: To compute the energy functional , we ignore the contributions from the points far away from the current computed point since computing a big number of the points in the whole image requires more computation. To decrease the impact caused by inhomogeneity that may arise far away, we widen a range of objects because narrow range only captures limited objects. Note that we do not consider the regularization term of the energy functional. For any point in the local patch, we use the mask function ( , ) to make sure that only on local patch information centered at can be computed. Thus, the contribution of the total energy functional is the sum of all neighborhood points along the evolving contour. In the following, we will use two steps to solve functional (11) in the following.
Step 1. Keeping ( ) fixed and computing the minimization of functional (11) with regard to 1 and 2 , we can easily get these constant functions by Step 2. Keeping the variables 1 and 2 fixed and computing the minimization of energy (11) with regard to ( ), we can get the following Euler-Lagrange equation: From the above equation, the variable ( ) can be expressed in the following way: ) .
The variable ( ) is then updated based on the change of the energy .

Numerical Approximation.
Since the energy functional (11) is nonconvex, it is difficult to solve the minimization in practice. Generally, the gradient descent schema driven by the Euler-Lagrange equation is usually applied to explicit time marching and causes local minimal. To solve this problem, inspired by the scheme developed by Krinidis and Chatzis [15] and Pereira et al. [16], we apply a fast numerical scheme to make its time step unconstrained in the explicit time marching. The algorithm can calculate the energy directly and judge if the degree of membership for any point is changed instead of solving the partial differential equation.

Lemma 1.
Let ∈ be a given point in the local patch Ω , the intensity value of point be 0 , and the corresponding degree of membership for this point be 0 . Correspondingly, let the new value of the degree of membership at the same point be ; the values of 1 and 2 will be changed to new two ones:̂1 and̂2. The new value of̂1 and̂2 could be calculated as follows: Proof. The two old constants 1 and 2 , which approximate the image intensity in local patch corresponding to ( ) > 0.5 and ( ) < 0.5 in (12), respectively, are written in the following forms: Assuming that we change the degree of membership for only one point when we compute the new degree Computational and Mathematical Methods in Medicine 5 of membership for the point , the constant̂1 will be obtained bŷ . In a similar way, we obtain the new valuê2 aŝ . Thus, the changed values Δ 1 =̂1 − 1 and Δ 2 =̂2 − 2 for the point can be very easily computed using formulations (17) and (18), respectively. This completes the proof.

Lemma 2.
Let the old total energy functional be and the new total energy functional bêwhen we change the degree of membership for the point into . Correspondingly, the changed energy Δ between the new and old total energy functional is given as follows: Proof. For the fixed point in the model, the change of the degree membership will lead to the change of the new energy functional. To compute the alteration Δ =̂− between the new and old total energy functional, we firstly calculate the new energy functional̂, which is written in the following form: From above, we can see that we should separately com-puter̃and̃in order to calculate the alteration. Sõ We will first compute the following equations ( ( ) −̃1) 2 and ( 0 −̃1) 2 : Then we will insert (22) into (20), we have the following equatioñ1: Computational and Mathematical Methods in Medicine For the image domain Ω , we havẽ In the same way, we can get=

Combining (20), (24), and (25), the new total energy functional is rewritten aŝ
Thus, it is very easy to compute the changes Δ =̂− by updating 1 and 2 , when the change at the degree of membership of a point is occurred.

Remark 3.
It is required that local region statistics must be calculated for all the points along the evolving curve. So it increases the complexity of the algorithm and the computation time. To improve the computational efficiency, we only update the membership ( ) in a narrow band region around the pseudo level set function ( ) = 0.5. In this paper, the computation of local statistics is separated into two parts. In the first part, the local region-based method begins by initializing every pixel in the narrow band with the local interior and exterior statistics. In the other part, the statistical models of all pixels within the narrow band neighborhoods are updated when any initialized pixel is crossed by the contour moving it from the interior to the exterior or vice versa.

Description of Algorithm
Steps. Here, the segmentation procedure of the proposed model is summarized as follows.
Step 1. Give an initial partition of the image, set 0 > 0.5 for one part and 0 < 0.5 for the other.
Step 3. Assume that the intensity value of point be 0 and the corresponding degree of membership for the point be 0 . Calculate the new degree of membership using (14) for current pixel 0 and the difference between the new and old energy Δ using (26).
If Δ ≥ 0, then change 0 with value; otherwise, keep the old value 0 .
Step 4. Repeat Step 3 to compute the total energy within the narrow band neighborhood using Jacobi iterations. When all pixels within the neighborhood image have been swept one time, one iteration is finished. The updated values of the current iteration are used for the next iteration.
Step 5. Repeat Steps 2-4 till the total energy remains unchanged.

Experiments and Results
In the section, we will validate that localized patch-based model can improve the performance of a given global energy instead of specially comparing with the fuzzy region-based global energies. To demonstrate the performance of the proposed model, we also test the segmentation results using the proposed model, the FEAC model [15], and the GL-FEAC model [16] on different synthetic images and real images. The experiment results will demonstrate that only the localized model can obtain a correct segmentation in these cases. All experiments are performed on a 1.86-GHz Intel dual-core notebook computer with Memory 3 GB using the MATLAB programming language. In these experiments, the parameter is 2. The code will be uploaded to the website: http://fangjx2005.com/. Figure 2 shows that the proposed model segments a synthetic image with three different objects. In the experiment, the parameter of the localization radius is set to 20. For each image region, the image includes different intensities and holes in the interior of the objects. In Figure 4, we applied the global FEAC model to segment the synthetic image with the same initial curve. The intermediate curve image, the segmentation result, and the evolution of the contours are shown in Figure 4 The next experiment is to test the robustness to noise of the proposed model and the results are shown in Figure 5. It is carried out on synthetic images mixed with different Gaussian noises using the FEAC model and the proposed model. The objects in the image include two kinds of shapes (circle and rectangle). From left to right, Figure 5 shows initial contour, the segmentation results using the FEAC model, the segmentation results using the proposed model, and the final object regions using the proposed model, respectively. From Figures 5(b)-5(d), it shows the segmentation results on conditions that the images are mixed with Gaussian noise with variances 0.01, 0.10, and 0.20, respectively. The results show that more noise causes too more small segmentation regions using the FEAC model. From these segmentation results, the proposed model can accurately extract the boundary while the FEAC model cannot get satisfactory result even if the image is populated by severe noise. Thus, we can see that our method shows the better robustness to noise.

Experiment Results.
In Figure 6, we further compare the proposed model with these global fuzzy active contour models, such as the FEAC model and the GL-FEAC model, on BIRD, MONKEY, and medical images. In the GL-FEAC model, the parameters are set as follows: = 0.5, 1 = 2 = 1. The images in Figure 6   show objects and backgrounds which are multimodal but that have intensities that change smoothly and quickly. Figure 6(a) shows the initial curve image. In Figures 6(b), 6(c), and 6(d), the boundary is obtained using the FEAC model, the GL-FEAC model, and the proposed model, respectively. From the segmentation results, we can see that the global energy finds only the brightest parts of the image while the localization stops on object boundaries.
In this paper, we use the popular error metric, the Dice coefficient [19], to quantitatively evaluate the performance of the competing methods. The Dice coefficient between two regions and is calculated as ( , ) = 2 × | ∩ |/(| | + | |), where | ∩ |, | |, and | | denote the pixel number of their union areas and and and , respectively. Obviously, the closer the Dice coefficient value is to 1, the better the segmentation results we will get. Table 1 depicts the Dice coefficient which gives the quantitative comparison scores in Figure 6.
We next illustrate the advantage of the proposed contour model on the PELVIS image in order to extract two object regions (the left obturator foramen and right obturator foramen) and the results are shown in Figure 7. The localization radius is 20. The PELVIS image with intensity inhomogeneity has blur boundary in some parts. The segmentation  Figure 7, the image has more seriously blurred boundary and is much more complicated. Figure 8 shows the segmentation results  using the FEAC model, the GL-FEAC model, and the proposed model with the same local radii. Figure 8(a) shows the initial contour image with two contours (two rectangles). In  Figure 9. Table 2 reports the Dice measure for segmentation results by the FEAC model, GL-FEAC model, and the proposed model in Figures 7 and 8. In the following experiment, Figure 10 shows different segmentation results for cervical spine images with noise. The typical schemes include the fuzzy active contour model with kernel metric (FAC-Ker) [17] and the global and local fuzzy active contour with the information (GLFAC) [18]. These models have the same initial positions. The initial contour image, the segmentation result, and the object image using the proposed model are shown in Figure 10(a). The results indicate that only the proposed model can extract two objects.
In Figures 10(b) and 10(c), it shows the segmentation results using the FAC-Ker and GLFAC models, respectively. From the segmentation results, we can see that the FAC-Ker model can lead to poor separation of different regions because the model is based on global image information. The GLFAC model can reduce the noise, but it cannot extract the objects.

The Scale
Parameter: Localization Radius. The important parameter in our model is localization radius which plays a key role in how local object region(s) the proposed model will extract. As such, it should be selected based on the scale of the extracted object of interest and approximation of the surrounding region of clutter. A small localization radius is selected when we attempt to extract small objects with nearby clutter and vice versa.
The example of lumbar spine image segmentation in Figure 11 illustrates the effect of different localization radii. In all these experiments, the number of the iterations is 200. The segmentation results are shown in Figures 11(b)-11(f) using the proposed model with localization radius 10, 15, 20, 25, and 30, respectively. With the same initialization, Figure 11 shows different results using five different local radii. From the segmentation results in Figures 11(b) and 11(c), the smallest radius size results in an incorrect segmentation that is too local, whereas the largest radius in Figures 11(e) and 11(f) leads to an incorrect energy value that is too global for the task at hand. Thus, the localization radius should be correctly chosen in terms of the nature of the objects so that the neighborhood is large enough to detect the desired boundary from the initialization. Table 3 describes the Dice measure for segmentation results with different localization radius.     Figure 7 216 × 151 0.164674 0.196280 0.916548 Lumbar spine image in Figure 8 199 × 171 0.128035 0.184621 0.854213

Conclusion
In this paper, we propose a localized patch-energy fuzzy active contour model by incorporating local image statistics for each pixel into the fuzziness of the energy, which can avoid local minima of the energy functional. In addition, we use a fast numerical method directly computing the changed value of the energy functional to update the curve evolution. Experiments show that the evolution of contour in our model is more stable which results in not only faster computation efficiency but also better performance of segmentation. Moreover, the alterations of the energy in the energy functional are calculated directly instead of solving the Euler-Lagrange equation. Thus, the active contour converges quickly to the object boundaries. Experimental results on synthetic and real images have validated the effectiveness of the proposed model.