DOI 10.1515/bmt-2013-0111      Biomed Tech 2014; 59(3): 219–229

Mirela (Visan) Punga, Rahul Gaurav and Luminita Moraru*

Level set method coupled with Energy Image features for brain MR image segmentation Abstract: Up until now, the noise and intensity inhomogeneity are considered one of the major drawbacks in the field of brain magnetic resonance (MR) image segmentation. This paper introduces the energy image feature approach for intensity inhomogeneity correction. Our approach of segmentation takes the advantage of image features and preserves the advantages of the level set methods in region-based active contours framework. The energy image feature represents a new image obtained from the original image when the pixels’ values are replaced by local energy values computed in the 3 × 3 mask size. The performance and utility of the energy image features were tested and compared through two different variants of level set methods: one as the encompassed local and global intensity fitting method and the other as the selective binary and Gaussian filtering regularized level set method. The reported results demonstrate the flexibility of the energy image feature to adapt to level set segmentation framework and to perform the challenging task of brain lesion segmentation in a rather robust way. Keywords: active contour model; intensity inhomogeneity; lesion segmentation; MRI. *Corresponding author: Luminita Moraru, Faculty of Sciences and Environment, Department of Chemistry, Physics and Environment, Dunarea de Jos University of Galati, 47 Domneasca St., 800008 Galati, Romania, Phone: +40745649014, Fax: +40236461353, E-mail: [email protected] Mirela (Visan) Punga: Aurel Vlaicu High School, 1 Decembrie 1918 St., 800511 Galati, Romania; and Faculty of Sciences and Environment, Department of Chemistry, Physics and Environment, Dunarea de Jos University of Galati, 47 Domneasca St., 800008 Galati, Romania Rahul Gaurav: Advanced Research and Techniques for Multidimensional Imaging Systems, Télécom Sud Paris, 9 rue Charles Fourier, 91011 Evry Cedex, Paris, France

Introduction Over the last few years, medical image segmentation has seen drastic changes and growing vigorously in an unimaginable way. But the underlying difficulty in this task due to the presence of noise, poor resolution, low contrast,

and intensity inhomogeneity makes it a rather challenging topic. The intensity inhomogeneity is the main challenge in our attempt to segment an image into homogeneous objects. There are various difficulties related to the boundary leakage between the gray matter (GM) and white matter (WM), due to intensity variation across image or low contrast. Another striking drawback is the anatomical variability of the brain, starting with the age parameter and ending with the tissues altered by pathologies placed in various locations in the brain [16]. Undoubtedly, it is one of the most interesting areas of research today, but still the necessity for development of effective automatic image analysis techniques to segment and quantify brain lesions in a robust way is imperative. The level set formulation in the region-based active contour framework offers the advantages of nonconstrained initialization and the ability to automatically capture the boundaries based on the topological changes. However, it is a large consensus that the level set method uses a complex set of controlling parameters, and this impairs the method. Thus, it becomes crucial for segmentation coupled to level set framework to properly define the conditions of the final location of the contour, namely, the level set in terms of evolution speed, mean curvature, contour length, and area. Following the mathematical foundation, the level set methods belong to the deformable segmentation models. In the last decade, level set functions are widely used [2, 6–8, 10, 12, 14, 17, 19–23]. In [21], the authors proposed a region-based active contour in a variational level set formulation for a more flexible initialization of the contours. They implement the algorithm for both two-phase segmentation and for multiphase segmentation specific to the brain MR images. However, they looked only to segment the WM and GM without taking various brain lesions into account. Wang et al. [20] used the region-based active contour model and the level set framework to handle the local intensity of images by using the Gaussian distributions with various mean and variance parameters. The local intensity means and variances drive the local Gaussian distribution and are the variables of the energy functional. They regularize the level set function by penalizing its deviation from a signed distance function or its length to derive a smooth contour. This method is able to segment images showing intensity

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

220      M. (Visan) Punga et al.: Level set method coupled with Energy Image features inhomogeneity and noise. A level set method based on Selective Binary and Gaussian Filtering Regularized Level Set is proposed in [23]. This algorithm allows the segmentation of objects showing weak edges encompassing both the property of local or global segmentation. The Gaussian filter regularizes the level set, and the standard deviation of the Gaussian filter is a critical parameter. Li et al. [7] proposed a fuzzy level set algorithm for automated medical image segmentation, which incorporates spatial information and eliminates the intermediate morphological operations. They have shown that their algorithm is less noise sensitive, and it is suitable for medical image segmentation. We focused here on the method that can enclose local statistics in a global framework analysis. In this paper, we proposed to run the level set method dealing with the local and global intensity fitting and with selective binary and Gaussian filtering segmentation algorithms in energy image feature framework in order to perform the lesion’s boundaries extraction. Our goal is to provide an efficacy segmentation of the various brain lesions. The challenge is that such lesions are not necessarily complete because the intensities of the lesion vary between values characteristic to undamaged tissue and values characteristic to the cerebrospinal fluid (revealing completely damaged areas). Furthermore, the lesion area shows inhomogeneity as having completely damaged core parts and minor damage in peripheral portions. As the inhomogeneity often affect the segmentation accuracy and severely impede the quantitative analysis, improvements are possible if energy image feature is included. While the particular algorithms used in our work are not new, their combination with energy image feature is unique and proves to be very efficient for inhomogenous noisy images. We evaluated our method in terms of its ability to accurately capture the brain lesion’s boundaries when it was tested on various scans of the same subject. Compared with the initial segmentation algorithms [20, 23], our method is improved mainly due to the incorporation of the spatial information during the mask scanning operation. The segmentation algorithms were developed in academic environment and are freely available for download.

The data sets used are freely available for download from “The whole brain atlas” database of the Harvard Medical School’s Whole Brain Atlas (http://www.med.harvard. edu/AANLIB/home.html). The variability of anatomic shapes was taken into consideration, and the following types of lesions were analyzed: bleed (code sample 018), infarction (code 012), tumor (code 010), and multiple sclerosis plaque (code 013). Supplementary, the codes T1 and T2 were added to specify the samples in T1w and T2w view, respectively. As shown in Figure 1, the proposed framework consists of two steps. First is the manual segmentation process. Second is the automatic process together with performance analysis and optimization criterion. We tested the segmentation on different MR images using the same sets of parameters and the same initialization setup. The MRI and CT data sets were of axial size (256 × 256) with an inplane resolution of 0.86 mm and a 3-mm slice thickness. As a matter of fact, there are no gold standards for brain lesions so, in order to evaluate the performance of the proposed segmentation method, we utilized the expertise of a radiologist to manually segment the lesions. These data sets have been previously labeled and segmented via a labor-intensive manual technique by one independent radiologist in order to locally hand draw the regions of anatomical interest. The most widely used software codes in neuroimaging analysis encompass a set of skull stripping, denoising operation, and automated segmentation routines. Among them, the mathematics-oriented models such as active contours with selective local or global segmentation Selective Binary and Gaussian Filtering Regularized Level Set (SBGFRLS) [23] written and developed at the Department of Computing, The Hong Kong Polytechnic University, and Local Gaussian Distribution local intensity of images by using the Gaussian distributions (LGIF) fitting energy [20] written by Li Wang at the University of North Carolina at Chapel Hill, USA, have been examined in this paper. This study will provide a quantitative and qualitative assessment for segmentation algorithms in the codes.

Overview of contributions

Materials and methods Data sets and the framework of segmentation algorithm Three sets of axial T1-weighted (T1w), T2-weighted (T2w), and CT images of various lesions were used in our study.

Mumford-Shah (MS) segmentation functional model considers a gray image I:Ω→R, where Ω is a 2D image space. An optimal contour Γ that segments the image I into disjoined regions is found via minimization of the functional energy FMS(f, Γ) = ∫Ω|f-I|2dx+μ∫Ω/Γ|∇f|2dx+v|Γ| where μ, v > 0 are constant parameters to balance the smoothing term and to ensure that the contour has minimal length [12].

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

M. (Visan) Punga et al.: Level set method coupled with Energy Image features      221

Start the algorithm: MRI or CT image acquisition I(x, y)

Manual segmentation

Image feature approach applying E=

2 Σ G–1 i=0 [p(i)] and 2D

convolution with a 3×3 mask Segementation process

SBGFRLS

LGIF

Get the contour of the segmentation results

Post process: Performance analysis (OR, USR, OSR)

N

Optimization criterion satisified?

Y Algorithm termination

Get the final contour of the segmentation results and statistical information along region

Figure 1 Flowchart explaining the steps involved in the segmentation framework.

|Γ| is the contour length, and f approximates the original image. The main drawback of this model consists of the existence of multiple local minima of the nonconvex function FMS. In 2001, Chan and Vese [3] proposed the “active contour without edges” model based on simplified MS segmentation functional [12]. To surpass the weakness of the MS model, energy functional having two smoothing functions defined in two exclusive subregions (inside C and outside C) is proposed:

FCV ( C , c1 , c2 ) = µ( length( C )) + v( area( insideC ))

+ λ1 ∫

inside ( C )

|c1 -I | 2 dx + λ2 ∫

inside ( C )

|c2 -I | 2 dx



(1)

where μ, v  ≥  0 and λ1, λ2 > 0 are constant parameters. c1 and c2 denote the average intensities inside and outside the contour, respectively. In order to minimize the functional in Eq. (1), the level set method was proposed by Osher and Sethian [13] as follows. The level set function or scalar Lipschitz function

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

222      M. (Visan) Punga et al.: Level set method coupled with Energy Image features φ:Ω→R; embeds the moving surface, and its zero level set φ(x, y, z, 0) = 0 is defined as [23]:  C = { x ∈Ω: φ( x ) = 0 }  inside ( C ) = { x ∈Ω: φ( x ) >0 }  outside ( C ) = { x ∈Ω: φ( x ) < 0 } 



(2) 

The components in Eq. (1) are reformulated in the level set framework as [1]:

length ( C ) = length ( φ = 0 ) =∫ Ω | ∇H ( φ )| dx area ( C ) = area( φ < 0 ) =∫ Ω H ( φ ) dx Average of I inside the curve C:∫Ω|c1-I|2H(φ)dx Average of I outside the curve C: ∫Ω|c2-I|2(1-H(φ))dx

Eq. (1) is transformed in the level set framework and for the entire domain Ω [1, 21]: F ( φ, c1 , c2 ) = λ1 ∫ | c1 -I | 2 H ( φ ) dx + Ω

λ2 ∫ | c2 -I | 2 ( 1-H ( φ )) dx +

(3)



µ∫ | ∇H ( φ )| + ∫ H ( φ ) dx









By minimizing Eq. (3) for a fixed level set function φ, the mean values c1 and c2 are: c1 ( φ ) = c2 ( φ ) =









I ( x ) H ( φ( x )) dx





H ( φ( x )) dx

and (4)

I ( x )( 1-H ( φ( x )) dx





( 1-H ( φ( x )) dx



The variational level set formulation is obtained minimizing the energy functional in Eq. (3), by keeping c1 and c2 formulated in Eq. (4) fixed. It is the moving equation to find the minimal surface:



   ∇φ  ∂φ = δ( φ )  µ div  - v + λ 1 | c 1 -I | 2 - λ 2 | c 2 -I | 2   ∂t  | ∇φ |   



(5)

Wang et al. [20] regularized the level set function by penalizing its deviation from a signed distance function 1 by ℘ ( φ ) = ∫ (| ∇φ( x )| −1) 2 dx and also by penalizing its Ω2 length to derive a smooth contour during evolution as L( φ ) = ∫ | ∇φ( x )| dx . Ω

Usually, the models taken into analysis, regularized the energy functional by using the approximated Heavi x  1 2 side, H ε ( x ) =  1 + arctan   and Dirac functions, 2  π  ε   1 ε δε ( x ) = H ε′( x ) = , x ∈ℜ. π ε2 + x 2

In [20], the experimental parameters are set to the following values: in a binary step function approach, the level set function φ is initialized as a constant c0 = 2, the spread of the Gaussian smoothing function σ = 3, time step evolution Δt = 0.1, the weighting coefficient of the penalty term ℘(φ) is μ = 1, and the weighting coefficient of the penalty term L( φ ) or the balloon force is v = 0.0008 × 255 × 255 for brain tumors. If the regulator parameter for Dirac function ε is too small, the energy functional will fall into a local minimum. A similar situation is encountered if ε is large. In [23], the level set formulation has a simplified approach as the level set function satisfied |∇φ| = 1 and the regularized term in Eq. (5) becomes Δφ. In this case, the authors used the Gaussian filtering process to regularize the level set function. Moreover, the standard deviation of the Gaussian filter is the same as the weighting coefficient of regularization strength μ. As a consequence, the term that incorporates the curvature information div(∇φ/|∇φ|)‧|∇φ| can be discarded. The following experimental parameters were used: ρ = 1, ε = 1.5, σ = 1, the Gaussian kernel size K = 5, and the time step Δt = 1. The velocity term α is the evolution speed parameter or balloon force. It controls the contour shrinking or expanding and also the propagation speed of the contour, and it is set for each image.

Image feature proposed algorithm The image feature approach is proposed as an image processing technique allowing a proper fit of a contour written in terms of the level set function. Both level set methods discussed above are characterized by the local properties of the Gaussian kernel functions, and such localization properties will generate many local minimums of the energy functional. In the brain MR images, the intensity inhomogeneity and noise are present, and multiple, adjacent regions exist, so the computational task increases accordingly. Starting at this point, we discuss how to couple the boundary smoothness and the increasing segmentation efficacy in our method with a lower computational cost. An intuitive approach to alleviate these limitations is to modify the inhomogeneity level by addressing to the second-order statistical feature, which intercorrelates the intensity value concentrations occurring at specific locations. In this formulation, the local texture homogeneity properties are taken into account. Energy image feature is estimated using an area of small size to reduce the possibility of mixing statistical information along region

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

M. (Visan) Punga et al.: Level set method coupled with Energy Image features      223

borders. The texture in the new image is smoothed by convolution but the region border steeped is maintained. Because in a previous study [11], the 3 × 3 mask size has been proven to provide salient image features and strong edges, we perform 2D convolution of an image with a 3 × 3 mask that represent the chosen feature pattern. We propose energy (or angular second moment in Gray Level Co-Occurrence Matrix) as relevant second-order statistical feature: E = ∑ Gi=−01 [ p( i )] 2 , where the histogram of the intensity levels is p(i) = h(i)/MN. It corresponds to an image represented by the function I(x, y), where x = 0, 1, …, N-1 and y = 0, 1, …, M-1. The function I(x, y) takes the value i = 0, 1, …G-1, G being the total number of the intensity level in the image. The choice of the energy feature kernel had the following main motivation. This kernel tends to be less sensitive to nearby edges as it largely weights the pixels near the sample position and less weights the distant pixels. A Matlab code was developed for energy image feature generating. This problem has not yet been addressed in this segmentation framework and is the novelty of our study.

Performance metric of proposed segmentation methods The segmented areas provide the basis for subsequent performance analysis. The performance analysis (or precision) of a segmentation algorithm allows setting the best segmentation parameters for efficacy segmentation, and also, a comparison between different segmentation approaches can be done. The region-based measures are used to evaluate the proposed segmentation algorithm [9]. The metrics used in our analysis were as follows. Overlap rate (or Jaccard similarity measure) [4, 5] is defined as the ratio of the intersection of segmented lesion area A and ground truth lesion area (or manually segmented areas) B to the union of segmented lesion area A and ground truth area B: | A∩ B | | A∪ B | Higher values of OR measure indicate the good probability of segmentation. Under segmentation rate defines the proportion of the unsegmented lesion area. It is defined as: OR =

|U | | B| where U = |B-(A∩B)| is the unsegmented lesion area. USR is expected to be low for superior performance of segmentation. USR =

Oversegmentation rate is defined as the ratio of the segmented nonlesion area V and the ground truth area B:

OSR =

|V | | A|



where V = |A-(A∩B)| is the segmented nonlesion area. For lower values of OSR, the segmentation performance is superior.

Experimental results The proposed energy image feature coupled with the level set-based segmentation method was checked for efficacy by comparing the results with manual segmentation. Thus, our method was tested in images with different types of lesions and injuries. The overlap or difference between our proposed method and manual segmentation method was computed. As suggested by [23], the setting parameters used in SBGFRLS method are σ ranges from 0.3 to 1.5 according to the noise level in the image and the type of the brain lesion, ρ = 1, ε = 1.5, and time step Δt = 1 assures the stability of the curve evolution. The σ parameter was set according to the different experiments. The velocity term was α = 10 and 15. The number of iterations was varied between 80 and 200. The zero level, which is used to represent the object contour, is in squared shape, and its size and location can be changed according to each object of interest in the image. The Matlab source code of the SBGFRLS method can be found in http://www4.comp.polyu.edu.hk/~cslzhang/. For the LGIF method [20], the following setting parameters were used: σ is 3 or 5, c0 = 2, ε is 0.7 and 0.8, λ1 = λ2 = 1, v ranges from 0.0003 × 255 × 255 to 0.0005 × 255 × 255, the regularization term μ = 0.1. The initial level set function is a circle. The Matlab source code of the LGIF method can be found in http://www.unc.edu/~liwa/. All these default setting parameters are kept unaltered when the segmentation methods are used in energy image features. The framework of the proposed approach is illustrated in Figure 1. Figures 2 and 3 show some segmentation results by LGIF and SBGFRLS methods and, also, the segmentation results by energy image feature coupled with LGIF and SBGFRLS methods, for various brain injuries and various types of images. Some test images yield satisfying segmentation results. For other images, considerable failure in contour delineation is visually provided because the boundary of this image is very weak and blurred. Note that the samples 010CT (SBGFRLS

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

224      M. (Visan) Punga et al.: Level set method coupled with Energy Image features

A

B

A

010T1

013T1

010T2

013T2

010CT

B

018C2

Figure 2 Examples of (A) original and (B) manual segmented images in axial views. The images were freely downloaded from “The whole brain atlas” database (http://www.med.harvard.edu/AANLIB/home.html).

coupled with energy image feature), 013T1 (LGIF coupled with energy image feature), 012T1, T2 and CT (LGIF) and 018T2 (SBGFRLS and SBGFRLS coupled with energy image feature) have not been segmented correctly. The segmentation performance results of the analyzed methods for representative sample of brain MR lesions and injuries are shown in Figure 4. A comparison between “manual vs. automated” delineation of boundaries is provided. The average segmentation performance results of 320 extracted lesions of brain for both analyzed and proposed segmentation algorithm are calculated and tabulated in Tables 1–3. The following remarks were deduced after careful visual analysis: –– In the case of 013T1 and 013T2 samples showing multiple sclerosis plaques, the small area manually delimitated by the physician has been encompassed in the large second area when the segmentation has been automatically provided. The local homogeneity criteria that describe the similarity of adjacent image objects led to the merging action. Small, neighbor plaques are merged when the homogeneity criterion is fulfilled. So, the size of the resulting image objects is higher. –– Focusing on the blurred boundaries, for example, 012T1, 012T2, and 012CT samples showing brain infarction, visual inspection of the brain images reveals the difficulties encompassed by both techniques (manual and automated) to delineate the area of the lesion.

–– The quality in CT images is affected by artifacts, low contrast and spatial resolution, and high noise. Overall, the efficiency or practical viability of the segmentation method in CT images is very low.

Discussions Only real brain MR images were used in these experiments to verify that our method is applicable in clinical research settings. For the comparisons, the performance metrics (Jaccard similarity measure and under- and over segmentation rate) are computed between the segmentations. Bleed, infarction, tumor, and multiple sclerosis plaque were segmented using our method against corresponding wellknown level set methods in region-based active contours framework and manual segmentation of these structures performed by a trained expert. Overall, our method was able to closely reproduce the accuracy of the manual segmentations because it constrains the boundaries to be smooth by employing a local statistics in a global framework analysis. The first challenge posed by this study was that there are no universally correct criteria to choose the parameters of segmentation because choices depend to some extent on the image type. Choosing the parameters to ensure an efficient segmentation was a time-consuming operation. Moreover, the presence of weak and diffused edges (due to edema around the tumor) leads to over segmentation. Despite the continuous efforts to enhance

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

M. (Visan) Punga et al.: Level set method coupled with Energy Image features      225

Figure 3 Examples of segmented images in axial views. Rows: 1, brain tumor T1w (sample 010T1); 2, brain tumor T2w (sample 010T2); 3, brain tumor CT (sample 010CT; the sample in column 3 is not segmented); 4, multiple sclerosis plaque T1w (sample 013T1; the sample in column 4 is not segmented); 5, multiple sclerosis plaque T2w (sample 013T2); 6, bleed T2w (sample 018T2; the samples in column 1 and 3 are not segmented). Columns: 1, SBGFRLS; 2, LGIF; 3, SBGFRLS coupled with energy image feature; 4, LGIF coupled with energy image feature.

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

226      M. (Visan) Punga et al.: Level set method coupled with Energy Image features Table 1 Mean and standard deviation of the overlap rate provided by the proposed method (SBGFRLS and LGIF coupled with energy image feature) and level set segmentation (SBGFRLS and LGIF) methods for analyzed dataset. Samples  

SBGFRLS

010T1 010T2 010CT 012T1 012T2 012CT 013T1 013T2 018T2 018CT

0.522 ± 0.109  0.319 ± 0.024  0.402 ± 0.088  0.423 ± 0.115  0.320 ± 0.178  0.470 ± 0.227  0.219 ± 0.062  0.330 ± 0.103  NS   NS  

                   



LGIF



0.770 ± 0.058  0.448 ± 0.071  0.626 ± 0.102  NS   NS   NS   0.585 ± 0.130  0.517 ± 0.368  0.426 ± 0.068  0.370 ± 0.117 

SBGFRLS coupled with   energy image feature

LGIF coupled with energy image feature

0.583 ± 0.045 0.333 ± 0.073 NS NS 0.332 ± 0.110 NS 0.537 ± 0.124 0.619 ± 0.016 0.472 ± 0.088 NS

0.791 ± 0.022 0.548 ± 0.064 0.533 ± 0.038 0.646 ± 0.120 0.506 ± 0.099 NS NS 0.781 ± 0.047 0.524 ± 0.078 0.523 ± 0.173

                   

NS, not segmented. Table 2 Mean and standard deviation of the undersegmentation rate provided by proposed method (SBGFRLS and LGIF coupled with energy image feature) and level set segmentation (SBGFRLS and LGIF) methods for analyzed dataset. Samples



SBGFRLS



010T1 010T2 010CT 012T1 012T2 012CT 013T1 013T2 018T2 018CT

                   

0.306 ± 0.065  0.260 ± 0.070  0.477 ± 0.110  0.502 ± 0.086  0.289 ± 0.081  0.188 ± 0.094  0.589 ± 0.022  0.314 ± 0.078  NS   NS  

LGIF



0.155 ± 0.045  0.368 ± 0.062  0.185 ± 0.042  NS   NS   NS   0.329 ± 0.039  0.268 ± 0.098  0.255 ± 0.046  NS  

SBGFRLS coupled with   energy image feature

LGIF coupled with energy image feature

0.268 ± 0.043 0.328 ± 0.034 NS NS 0.204 ± 0.056 NS 0.260 ± 0.076 0.193 ± 0.060 0.036 ± 0.009 NS

0.107 ± 0.014 0.255 ± 0.042 0.247 ± 0.027 0.297 ± 0.052 0.541 ± 0.121 NS NS 0.054 ± 0.020 0.055 ± 0.010 0.332 ± 0.066

                   

NS, not segmented. Table 3 Mean and standard deviation of the oversegmentation rate provided by proposed method (SBGFRLS and LGIF coupled with energy image feature) and level set segmentation (SBGFRLS and LGIF) methods for analyzed dataset. Samples



SBGFRLS



010T1 010T2 010CT 012T1 012T2 012CT 013T1 013T2 018T2 018CT

                   

0.364 ± 0.083  0.672 ± 0.098  0.456 ± 0.058  0.325 ± 0.062  0.648 ± 0.105  0.623 ± 0.122  0.447 ± 0.072  0.440 ± 0.084  NS   NS  

LGIF



0.352 ± 0.056  0.509 ± 0.075  0.143 ± 0.031  NS   NS   NS   0.313 ± 0.056  0.243 ± 0.040  0.206 ± 0.019  0.523 ± 0.053 

SBGFRLS coupled with   energy image feature

LGIF coupled with energy image feature

0.493 ± 0.094 0.715 ± 0.116 NS NS 0.606 ± 0.110 NS 0.055 ± 0.020 0.240 ± 0.044 0.335 ± 0.046 NS

0.263 ± 0.033 0.489 ± 0.081 0.415 ± 0.045 0.307 ± 0.065 0.201 ± 0.028 NS NS 0.121 ± 0.026 0.251 ± 0.018 0.217 ± 0.031

                   

NS, not segmented.

the performance of the level set methods, still there are inaccurate contour detections especially for boundaries presenting small gaps. Our goal was to keep the setting

parameters used in the analyzed methods unchanged over the segmented images. The balloon force parameters (v in [20] and α in [23]) and the spread of Gaussian

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

M. (Visan) Punga et al.: Level set method coupled with Energy Image features      227

A 0.8

Overlap rate

SBGFRLS

SBGFRLS

SBGFRLS_Energy

1.0

LGIF

LGIF SBGFRLS_Energy

0.8

LGIF_Energy

0.6

Overlap rate

LGIF_Energy

0.6 0.4

0.4 0.2

0.2

0.0 0.0

B

1-1 1-2 1-3 1-1 1-2 1-3 1-1 1_2 1_3 013T1

010T1 010T2 010CT 012T1 012T2 012CT 018T2

Under segmentation rate

0.8

SBGFRLS_Energy

0.6 0.4

SBGFRLS LGIF

0.6

LGIF LGIF_Energy

018CT

Under segmentation rate

0.8

SBGFRLS

013T2

SBGFRLS_Energy LGIF_Energy

0.4 0.2

0.2 0.0 1-1 1-2 1-3 1-1 1-2 1-3 1-1 1_2 1_3

0.0

013T1

010T1 010T2 010CT 012T1 012T2 012CT 018T2

C

Over segmentation rate

0.8

SBGFRLS_Energy LGIF LGIF_Energy

0.4

018CT

Over segmentation rate SBGFRLS

0.6

013T2

0.8

SBGFRLS LGIF

0.6

SBGFRLS_Energy LGIF_Energy

0.4 0.2

0.2 0.0 1-1 1-2 1-3 1-1 1-2 1-3 1-1 1_2 1_3

0.

013T1

010T1 010T2 010CT 012T1 012T2 012CT 018T2

013T2

018CT

Figure 4 Performance results of the segmentation algorithms for representative brain lesion areas. (A) Overlap rate. (B) Undersegmentation rate. (C) Oversegmentation rate.

function σ were tuned for different brain images in order to enhance the segmentation ability. The optimization problem was employed within each frame. Segmentation setting parameters were optimized as a result of a supervised segmentation. For this purpose, the overlap rate between manual and segmented frames was used. Also, the visual assessment is considered. However, these operations increase the computational burden. For this point of view, the image feature method improved the computational time in average with 30%. In order to avoid the reinitialization and the irregularities of level set functions during evolution, Zhang et al. [23] used a Gaussian filter to regularize the level set function, and Wang et al. [20] regularize the level set function by penalizing its deviation from a signed distance function and its length to derive a smooth contour during evolution. These penalty terms may cause the numerical

instability when the level set function becomes very flat or steep near the zero level set. Instead of imposing the smoothing step in algorithm, we proposed that the algorithms run in images with smoothing boundaries. In this case, the normal derivative (provided by Neumann boundary condition) can exist and allows the minimization of the energy functional by solving the descent gradient flow equation. Thus, the main advantage of the energy image feature consists of its sensitivity to the image scale due to smoothing operation. Segmentation is most successful when there is little overlap in distributions of pixel values from the different objects in an image. Noise is one cause of overlap, and it is reduced by using our method. Wang et al. [21] used a multiphase segmentation based on two initial contours in 2D MR images. In order to reduce the number of iterations, they accomplished a preliminary segmentation of the background and cerebral spinal fluid

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

228      M. (Visan) Punga et al.: Level set method coupled with Energy Image features (CSF) and another preliminary segmentation of the GM and WM followed by dilatation and erosion operations. These additional operations led to an increased computational cost. Supplementary, each class (i.e., GM, WM, and CSF) has an associated probability density function (pdf). Generally, the segmentation methods assume that pdf depends on the region but is invariant in the region. This will raise difficulties when the intensity inhomogeneity manifests. As a consequence, the Gaussian distribution with local intensity mean and standard deviation is used to define the pdf. We note an important overlap between Gaussian pdf of WM and GM. As a result, many misclassifications between WM and GM exist. To overcome this drawback, we propose taking the spatial information into account. Looking at brain images with various pathologies, the following issues are highlighted. Multiple sclerosis lesion is a disease of WM, and the segmentation methods produce reliable result. The performance metrics for our methods are slightly better. On the contrary, the analysis of the brain infarction is a very difficult task. This problem is a nonformal issue because the lesion may be connected to other biological objects with similar intensities. Cerebral infarction may include the cortical GM and also CSF. The morphology and tissue-characterization parameters may vary strongly in ischemic regions, indicating various pathological processes. In the case of 012 samples, the LGIF algorithm fail to segment the lesion in this series, but it can be observed that the LGIF coupled with energy image feature produces better results for these difficult cases. We also tested our method for CT images. It is evidently that the segmentation results are undesirable, mainly due to the weak boundaries and very low contrast in CT images. When the zero level set reaches the object boundary, the driving force will be close to zero, and this should make the zero level set stopping at the object boundary. However, the weak edge object can impede the regularization force by zero, making the zero level set continue to evolve, and finally causing the boundary leakage problem and failing the segmentation. Clearly, these pathologies with global extension like cerebral infarction raise problems to segmentation. The drawback of the studied methods consist of very unstable balance between the time step, the epsilon values of the Heaviside and Dirac functions, and the evolution speed parameter in order to prevent the boundary leakage. Owing to pronounced intensity inhomogeneity of brain images, these methods are also sensitive to the initialization of the contour. In our experiment, the mask size is an important parameter. In real images, the noise is usually in high

frequencies and can create false edge pixels. Also, high gradient pixels are usually not neighbors and not connected. This impedes the accurate segmentation. The segmentation works by selecting the local maxima, which are proportional to the degree of smoothing. If the degree of smoothing is small, the edge detection is very accurate but very sensitive to noise. On the contrary, if the degree of smoothing is high, the localization is less accurate but less sensitive to noise. We performed the image transformation using a 3 × 3 mask, and as result, the texture is smoothed by convolution, but the region border steeped is maintained. In this manner, the image is smoothed without destroying the edges, and our method demonstrates consistency if the right segmentation parameters are available. A possible limitation of our method (if we compare with other reported methods) is the lack of preprocessing intervention as skull stripping or non-brain tissue removal. If we compared our results with those reported by Shattuck et al. [18], the Jaccard metric is similar with the results reported for Brain Surface Extractor (BSEv08a) method for default setting parameters. Thus, our average Jaccard coefficient of 0.597 is similar to the best score by BSVv08a (namely, 0.595). Saad et al. [15] presented brain lesion segmentation based on the thresholding technique and gray level co-occurrence matrix. The intensity was enhanced using two different algorithms, namely, the γ-law transformation and contrast stretching. They reported the following values for area overlap parameter: 0.489 for brain tumor, 0.701 for infarction, and 0.692 for hemorrhage when the γ-law transformation is used and 0.264 for brain tumor, 0.649 for infarction, and 0.676 for hemorrhage when contrast stretching is done. Obviously, the comparison presented here is not exhaustive, and it would be interesting to include more image features. However, the comparison seems to indicate that energy statistical feature perform well. An evaluation of the feature descriptors in the context of texture smoothness and evaluation of the segmentation accuracy will be a useful and valuable addition to our work.

Conclusions The proposed method provides encouraging results for segmentation of brain lesions even if it is in its initial stages of development. We have shown that the brain MR images can be segmented using the concept of energy image features. The primary intent of this study was to segment brain lesions, which are typically intensity

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

M. (Visan) Punga et al.: Level set method coupled with Energy Image features      229

inhomogeneous images. Signal-intense lesions (e.g., certain brain tumors) might be segmented as well. Problems are encountered for lesions, which mixed low- and high-intensity subregions and in CT images. Experiment results show that our method is feasible to medical images and deserves further research. In particular, the improvements achieved with our method in terms of accuracy compared to local and global intensity fitting and selective binary and Gaussian filtering regularized

level set methods are demonstrated by brain segmentation experiments. Acknowledgments: The authors would like to thank Dr. Adina-Geanina Nămoianu, St. Andrew Emergency Hospital, Galati, Romania, for useful discussions. Received October 17, 2013; accepted February 10, 2014; online first March 4, 2014

References [1] Angelini ED, Song T, Mensh BD, Laine A. Segmentation an quantitative evaluation of the brain MRI data with a multi-phase three-dimensional implicit deformable model, In: Fitzpatrick JM, Sonka M, editors. Medical imaging: image processing. Proceedings of SPIE, vol 5370, Bellinham, WA 2004. [2] Balafar MA, Ramli AR, Saripan MI, Mashohor S. Review of brain MRI image segmentation methods. Artif Intell Rev 2010; 33: 261–274. [3] Chan TF, Vese LA. Active contours without edges. IEEE Trans Image Process 2001; 10: 266–277. [4] Ciecholewski M, Chocholowicz J. Gallbladder shape extraction from ultrasound images using active contour models. Comput Biol Med 2013; 43: 2238–2255. [5] Elter M, Held C, Wittenberg T. Contour tracing for segmentation of mammographic masses. Phys Med Biol 2010; 55: 5299–5315. [6] Eskildsen SF, Coupé P, Fonov V, et al. BEaST: brain extraction based on nonlocal segmentation technique. NeuroImage 2012; 59: 2362–2373. [7] Li BN, Chui CK, Chang S, Ong SH. Integrating spatial fuzzy clustering with level set methods for automated medical image segmentation. Comput Biol Med 2011; 41: 1–10. [8] Ma B, Wu Y, Li P. Level set segmentation using image second order statistics. In: Zhang T, editor. Automatic target recognition and image analysis. Proceedings of SPIE, vol 8003, Nong Sang, Guilin, China 2011. [9] Mallikarjuna PB, Guru DS. Performance Evaluation of Segmentation and Classification of Tobacco Seedling Diseases. Int J Mach Intell 2011; 3: 204–211. [10] Mitchell IM. The flexible, extensible and efficient toolbox of level set methods. J Sci Comput 2008; 35: 300–329. [11] Moraru L, Moldovanu S. Comparative study on performance of textural image features for active contour segmentation. Sci China Life Sci 2012; 55: 637–644. [12] Mumford D, Shah J. Optimal approximation by piecewise smooth functions and associated variational problems. Commun Pure Appl Math 1989; 42: 577–685.

[13] Osher S, Sethian JA. Fronts propagating with curvaturedependent speed: algorithms based on Hamilton-Jacobi formulation. J Comput Phys 1988; 79: 12–49. [14] Paragios N. A level set approach for shape-driven segmentation and tracking of left ventricle. IEEE Trans Med Imaging 2003; 22: 773–776. [15] Saad NM, Abu-Bakar SAR, Abdullah AR, Salahuddin L, Muda S, Mokji M. Brain Lesion segmentation from diffusion weighted MRI based on adaptive thresholding and gray level co-occurrence matrix. Journal of Telecommunication, Electronic and Computer Engineering 2011; 3: 1–14. [16] Schlaggar BL, Brown TT, Lugar HM, Visscher KM, Miezin FM, Petersen SE. Functional neuroanatomical differences between adults and school-age children in processing of single words. Science 2002; 296: 1476–1479. [17] Schmidt P, Gaser C, Arsic M, et al. An automated tool for detection of FLAIR-hyperintense white-matter lesions in multiple sclerosis. NeuroImage 2012; 59: 3774–3783. [18] Shattuck DW, Prasad G, Mirza M, Narr KL, Toga AW. On line resource for validation of brain segmentation methods. NeuroImage 2009; 45: 431–439. [19] Suri JS, Liu K, Singh S, Laxminarayan SN, Zeng X, Reden L. Shape recovery algorithms using level sets in 2-D/3-Dmedical imagery: a state-of-the-art review. IEEE Trans Inf Technol Biomed 2002; 6: 8–28. [20] Wang L, He L, Mishra A, Li C. Active contours driven by local Gaussian distribution fitting energy. Signal Process 2009; 89: 2435–2447. [21] Wang L, Li C, Sun Q, Xia D, Kao C-Y. Active contours driven by local and global intensity fitting energy with application to brain MR image segmentation. Comput Med Imaging Graph 2009; 33: 520–531. [22] Wang L, Shi F, Lin W, Gilmore JH, Shen D. Automatic segmentation of neonatal images using convex optimization and coupled level sets. NeuroImage 2011; 58: 805–817. [23] Zhang K, Zhang L, Song H, Zhou W. Active contours with selective local or global segmentation: a new formulation and level set method. Image Vision Comput 2010; 28: 668–676.

Brought to you by | University of Manitoba Authenticated Download Date | 6/10/15 9:43 PM

Level set method coupled with Energy Image features for brain MR image segmentation.

Up until now, the noise and intensity inhomogeneity are considered one of the major drawbacks in the field of brain magnetic resonance (MR) image segm...
3MB Sizes 2 Downloads 3 Views