Segmentation of Left and Right Ventricles in Cardiac MRI Using Active Contours

Segmentation of left and right ventricles plays a crucial role in quantitatively analyzing the global and regional information in the cardiac magnetic resonance imaging (MRI). In MRI, the intensity inhomogeneity and weak or blurred object boundaries are the problems, which makes it difficult for the intensity-based segmentation methods to properly delineate the regions of interests (ROI). In this paper, a hybrid signed pressure force function (SPF) is proposed, which yields both local and global image fitted differences in an additive fashion. A characteristic term is also introduced in the SPF function to restrict the contour within the ROI. The overlapping dice index and Hausdorff-Distance metrics have been used over cardiac datasets for quantitative validation. Using 2009 LV MICCAI validation dataset, the proposed method yields DSC values of 0.95 and 0.97 for endocardial and epicardial contours, respectively. Using 2012 RV MICCAI dataset, for the endocardial region, the proposed method yields DSC values of 0.97 and 0.90 and HD values of 8.51 and 7.67 for ED and ES, respectively. For the epicardial region, it yields DSC values of 0.92 and 0.91 and HD values of 6.47 and 9.34 for ED and ES, respectively. Results show its robustness in the segmentation application of the cardiac MRI.


Introduction
Cardiac MRI is a noninvasive imaging methodology, which is used to obtain the anatomical data of a heart for clinical diagnosis of cardiovascular analysis [1]. In cardiac MRI, volumetric analysis of left ventricle (LV) and right ventricle (RV) is an initially quantified method to diagnose cardiac contractile function. A detailed understanding of the cardiac contractility is essential in the quest to prevent, diagnose, and treat heart-related disorders [2,3]. Therefore, segmentation of LV and RV plays a crucial role in the detection and prevention of the heart attacks, which is a common cause of mortality in this century. Manual segmentation is a hectic and time-consuming job for both radiologists and cardiologists. Moreover, it can lead to a high number of false positives due to human errors, such as fatigue and distractions. Therefore, a computer-aided diagnosis (CAD) system is needed, which can assist both radiologists and cardiologists to accurately segment LV and RV boundaries in cardiac MRI.
To date, numerous segmentation methods [4][5][6][7][8][9][10][11][12][13][14][15][16] have been proposed to segment either one or both LV and RV. In MRI, weak or blurred edges and intensity inhomogeneities are the problems, which make it difficult for the intensitybased segmentation methods to properly extract the regions of interests (ROI). Therefore, accurate segmentation of LV and RV is still an open challenge to the researchers in the area of medical image segmentation.
Active contours are one of those methods used to segment LV and RV in cardiac MRI. Active contour model also known as snakes was introduced by Kass et al. in late 80s in which a curve is evolved towards object boundaries to delineate a region of interest [17]. Since then, active contours have been adapted abundantly in numerous image segmentation techniques. The basic idea of the active contours is to control the deformable curve and restrict it to the object boundary by minimizing a force known as a balloon force.
There are many variants of active contours, which are classified as edge-based [15][16][17][18][19][20] and region-based methods 2 Computational and Mathematical Methods in Medicine [21][22][23][24][25][26][27][28][29][30]. These methods have been contemplated for image segmentation, where their points of interest are beneficial. Edge-based methods use image gradient information to deform the level set curve towards the object boundary. However, these types of methods are not effective on images with weak edges, noise, and blurred boundaries. On the other hand, region-based models use a region-based descriptor that exploits image statistical information to evolve the contour. Among all region-based methods, Chan-Vese [23] is a broadly utilized method that is able to properly segment images with homogeneous regions. However, this method is unable to generate acceptable segmentation results in the presence of intensity inhomogeneity, which is a well-known problem in medical imaging.
A region-based active contour method proposed by Min et al. adopts an intensity-based global division algorithm [31]. In this method, two intensity means are computed for inner and two for the outer region. Therefore, it showed better performance compared to Chan-Vese method [23]. However, this method is also unable to produce acceptable segmentation results in the presence of intensity inhomogeneity. Numerous methods have proposed a viable solution to segment inhomogeneous regions by introducing image local information in their models [22][23][24][25][26].
A local binary fitting (LBF) method for image segmentation is proposed by Yezzi Jr. et al. in the context of intensity inhomogeneity [32]. In this method, a Gaussian kernel is introduced in the energy formulation to exploit image local information. A localized active contour method (LAC) is devised by Lankton and Tannenbaum in which global region-based methods are reformulated by replacing global means with image local information [28]. These methods can segment intensity inhomogeneous regions, unlike their global counterparts. However, the techniques explained in these methods are sensitive to the position of initial contour. Moreover, they also have high computational cost due to the complicated local information in their formulation.
Recently, hybrid methods [27][28][29][30][31] gained popularity among region-based methods. These methods either combine both region (local or global) and edge information or both local and global region information in their energy formulations. In [33], Zhang et al. proposed a method which combines the advantages of edge-based and regionbased active contours. In this paper, a region-based signed pressure force (SPF) function is also introduced which utilizes image global intensity means from Chan-Vese method [23]. This method adapts similar approach from geodesic active contour (GAC) model. However, in their model the edge-indicator function is replaced with a region-based SPF function; moreover, the traditional regularization function is also replaced with a Gaussian smoothing. This method only uses global image intensity information; therefore, it is unable to segment intensity inhomogeneous images.
In [34], Wang et al. introduced a new energy formulation in which image local and global information from LBF and Chan-Vese methods are incorporated in an additive manner. This method is capable of handling intensity inhomogeneities and yields better segmentation results compared to the stateof-the-art methods. However, this method is sensitive to the position of initial contour. Moreover, for different types of images, the scaling constant for the additive local and global intensity distributor varies; therefore, it is very difficult to choose the correct value. Recently, Liu e al. [13] and Yang et al. [14] proposed active contour segmentation methods for cardiac MRI data. These methods segment endocardial and epicardial boundaries by using level set formulation.
In this paper, a hybrid region-based active contour method is proposed to segment left and right ventricles in a cardiac MRI. The main contribution of this work is the formulation of new SPF function, which incorporates both local and global image intensity information. Local and global intensity information is obtained from Lankton and Tannenbaum [28] and Min et al. [31] methods, respectively. Global intensity term in the proposed SPF function helps to segment objects with big intensity difference and the local intensity term plays its role to segment objects with intensity inhomogeneity. The integration of local and global intensities has overcome the limitations of Lankton and Tannenbaum and Min et al. methods, that is, sensitivity towards the initial position of contour and inability to segment intensity inhomogeneous regions. In the proposed method, a Gaussian kernel is also used to regularize the level set and eliminate the need of expensive reinitialization.
In order to understand the qualitative results discussed in the result section, it is necessary to understand the physical structure of both ventricles inside cardiac MRI. In short axis view of heart MRI, the myocardium is the dark area between two concentric circles surrounded by LV and RV. Figure 1 shows LV and RV along with other important regions in cardiac MRI. Epicardium is a wall between myocardium, surrounded organs, and tissues. Endocardium wall covers the LV and RV cavity. In the presence of grey level blood flow and intensity inhomogeneity, segmentation of endocardium is a quite challenging task for both cardiologists and radiologists to quantify the cardiac contractile function. In this paper, a segmentation method is proposed which is able to segment intensity inhomogeneous regions and delineate weak object boundaries. Consequently, it helps to segment endocardium and epicardium walls for both left and right ventricles. This paper is organized as follows. Background and previous methods are briefly discussed in Section 2. The main idea and formulation of the proposed method are presented Computational and Mathematical Methods in Medicine 3 in Section 3. Experimental results and comparisons are shown in Section 4. Quantitative analysis with the state-ofthe-art methods is discussed in Section 5. Finally, conclusion and future work are given in Section 6.

Related Work and Background
2.1. Chan-Vese Method. Chan and Vese [23] proposed a region-based active contour method based on Mumford-Shah [22]. Let : Ω ⊂ 2 be an input image, : Ω ⊂ 2 a level set function, and a closed curve corresponding to the zero level set: = { ∈ Ω | ( ) = 0}. The Chan-Vese formulation is defined as follows: where ( ) represents the length of the curve, which is used to regularize the contour . (in( )) represents the area inside the curve . 1 , 2 are image intensities inside and outside of the curve . ≥ 0, V ≥ 0, and ( 1 , 2 ) > 0 are scaling constants and ( ) is the regularized Heaviside function which is defined as where constant controls the smoothness of the Heaviside function. By minimizing (1) with respect to 1 , 2 , and using the steepest gradient descent [35] the following definitions and solution equations are acquired: In (5), ( ) is a smooth version of the Dirac delta function, which is defined as where constant controls the width of the Dirac delta function. This method is widely used to segment the images with uniform intensity distribution. However, it cannot properly segment images with intensity inhomogeneity.

Min et al. Method.
In [31], a region-based active contour method with modified global region term is proposed to solve complicated intensity difference problem in Chan-Vese method. The proposed energy function computes two intensity means (big and small) for inside and two for outside of curve . This extra information from new means provides better segmentation experience compared to Chan-Vese method. However, this method is sensitive to noise. Their energy functional based on a novel region-based term is defined as follows: where ( ( ) − ) is a division function based on the intensity means from both inside and outside of object at = 1, 2. By minimizing (7) with respect to 11 , 12 , 21 , and 22 , using the steepest gradient descent [35], the following definitions are acquired: where 11 and 12 are big and small intensity means inside the contour. Similarly, 21 and 22 represent the big and small intensity means outside the contour. By using four values of intensity means, unlike Chan-Vese method, this method fairly improves segmentation accuracy. Figure 2 shows segmentation of a complicated (texture like) region using both Chan-Vese [23] and Min et al. [31] methods. Figure 3(a) shows that Chan-Vese method produces an unacceptable segmentation result. This method takes small black dots as separate regions that lead to an unacceptable segmentation of the big rectangular region, which is the actual region of interest. In turn, Min et al. method is able to properly segment the rectangular region as shown in Figure 2(b).

Zhang et al. Method.
Primarily an edge-indicator function was proposed in GAC method [19] to segment an object by evolving the level set curve towards object boundaries. However, this method was not able to segment global structure of the object. In [26], Li et al. proposed a region-based segmentation method, which combines the advantages of GAC and Chan-Vese methods. In this method, a region-based signed pressure force (SPF) function is introduced to replace the edge-indicator function in GAC method. SPF function helps to segment global structure of the given image using intensity means from inside and outside of the curve. A Gaussian kernel is used to regularize the level set, which also removes the need of its reinitialization. Let : Ω ⊂ 2 be the given image and : Ω ⊂ 2 a level set curve; then the solution to their energy functional is defined as where > 0 is a scaling parameter and spf( ) is the SPF function, which is defined as In (13), the value of the SPF function is in the range [−1, 1]. It shrinks the contour when it is defined outside and expands when defined inside of the object. 1 and 2 are the image intensity means defined in (3) and (4), respectively. The energy functional is defined as follows:

Localized Chan-Vese
where ( ) and ( ) are smooth versions of Heaviside and Dirac delta functions, defined in (2) and (6). ( , ) is a mask function which helps to detect the local regions in terms of radius defined as follows: ( , ) will be 1 when point is within the mask of radius centered at ; otherwise, it will be 0. In Figure 3, the local neighborhood is shown in green and the contour is shown in red. The process to compute the local intensity information for both regions inside and outside of the contour is shown in Figures 3(a) and 3(b), respectively.
By substituting ( ( ), ( )) in (14) following localized energy function is formulated: where 1 and 2 are the local intensity means inside and outside of the contour, respectively. These intensities localized by mask function ( , ) at point are defined as This method is capable of segmenting images with intensity inhomogeneities. However, it is sensitive to the initial position of the contour.

The Proposed Method
In this paper, a novel region-based active contour method is formulated, which uses both local and global image statistical information. An energy functional is devised using a new region-based SPF function based on Lankton and Tannenbaum [28] and Min et al. [31] intensity terms, which drives the zero level set curve towards the object boundary. During the curve evolution, its movement depends on the sign of the SPF function, which moves inwards if SPF is positive and outwards if it is negative. Let : Ω ⊂ 2 be the given image and : Ω ⊂ 2 a level set curve; then the proposed energy functional is defined as where V ≥ 0 and > 0 are the scaling constants.
( ) is a length and ( ) is an area term, defined as In the above equations, ( ) and ( ) are smooth versions of the Heaviside and Dirac delta functions, which are defined in (2) and (6), respectively. The energy ( ) controls the inner force, which computes the region information inside and outside of the curve and evolves the curve towards the object boundary. In turn, ( ) term deals with the curvature of the object boundary and regularizes the curve. spf ( ) is the proposed SPF function, which incorporates local and global characteristics from Kass et al. [17] and Lankton and Tannenbaum [28] methods, respectively. Like a traditional SPF function, the value of the proposed SPF function is in the range [1, −1]. The only difference is that it uses both local and global intensity means in its formulation. The proposed SPF function is defined as where spf ( ) and spf ( ) are global and local SPF functions, respectively. is a scaling constant whose value lies in 0 ≤ ≤ 1. It decides which SPF term will play a key role during the contour evolution. When is close to 0 then spf ( ) with the local terms will be dominant. In turn, when is close to 1 then spf ( ) with the global terms will be dominant. Selection of the parameter depends on the type of image used for the segmentation. If the given image has a uniform intensity distribution then should be close to 1. In turn, if the given image has intensity inhomogeneity then should be close to 0. The global SPF function spf ( ) used in (22) is defined as follows: where GFI is a global fitted image defined as GFI = ( 11 + 12 ) ( ( )) where 11 , 12 , 21 , and 22 are the global intensity means, which are defined in (8), (9), (10), and (11), respectively. Similarly, the local SPF spf ( ) used in (22) is defined as where 1 and 2 are the local intensity means, which are defined in (18) and (19), respectively. In (23) and (25), a characteristic term is used to restrict the evolution of the contour inwards (towards the region of interest), which is defined as where is a nonnegative integer. By minimizing (20) with respect to using the steepest gradient descent [30], the following solution is obtained: where spf ( ) is the proposed SPF function which incorporates both local and global intensity information in an additive manner as defined in (22). In the traditional level set methods [12,14,15], the level set function is reinitialized after each time step for smooth transitions during the curve evolution process. In [13], an energy penalization term is introduced by Liu et al., which maintains the level set function as a signed distance function (SDF) to remove the need of reinitialization. In this paper, a Gaussian kernel is used, which not only regularizes the level set curve but also removes the computationally expensive reinitialization step. For outmoded level set methods, it is essential to initialize level set function as a signed distance function (SDF). The initial level set function is defined as where > 0 is a constant (in this paper = 2). Finally, the iterative steps of the proposed method are summarized as follows: (1) Initialize the level set function using ( , = 0) from (28). (2) Initialize the characteristic term using 0 from (25).

Experimental Results
The proposed method is validated on 15 training and 15 validation datasets from the MICCAI 2009 [36] for left ventricle segmentation. Using Test 1 dataset from MICCAI 2012 [37] proposed method is also validated for right ventricle segmentation, respectively. These datasets contain cardiac MRI short axis volumetric data along with their respective ground truths. The proposed method is implemented in MATLAB 8.5 version in Windows 8 environment using 2.97 GHz Intel Core-i7 processor with 4 GB RAM.

Parameters Selection.
For the proposed method the following parameters are used for all the experiments: = 1, V = 0.002 × 255 × 255, Δ = 1.0, = 10, = 0.03, and = 1.5. It is critical to properly tune the parameter , which controls the level segmentation of homogeneous and inhomogeneous regions. ranges between 0 and 1; it is chosen small for intensity inhomogeneous objects while for homogenous intensity regions is chosen near to 1.
Initially, the proposed method has been tested on both synthetic and real images. Figure 4 shows the qualitative segmentation comparison with the state-of-the-art using the synthetic images. Column (a) shows the original images with the initial contour, columns (b), (c), and (d) show the segmentation results using the Lankton and Tannenbaum [28], LBF [25], and proposed method, respectively. Similarly, the segmentation results using the real images are shown in Figure 5. Column (a) shows the original images with the initial contour, columns (b), (c), and (d) show the segmentation results using the Lankton and Tannenbaum [28], LBF [32], and proposed method, respectively. Results show that the proposed method yields high accuracy compared to the stateof-the-art methods, which are unable to properly segment the regions of interest.

Endocardium Segmentation.
The proposed method is used to delineate the endocardial boundaries of both LV and RV from cardiac MRI. Figure 6 demonstrates the segmentation results of both left and right ventricles from the cardiac MRI using two different databases [33,38], where columns (a) and (d) show the images with the initial contours from the MICCAI 2009 [36] and MICCAI 2012 [37] databases, respectively. Columns (b) and (e) show intermediate segmentation results using the proposed method. In turn, columns (c) and (f) show final segmentation results using the proposed method. Figure 7 shows the left ventricle segmentation using the cardiac MRI from the MICCAI 2009 [36] database. Similarly, Figure 8 shows the right ventricle segmentation using the cardiac MRI from the MICCAI 2012 [37] database. The qualitative comparison of the segmentation results of the proposed method with the LBF [25] and Lankton and Tannenbaum [28] methods is also shown in Figures 7 and  8. In Figures 7 and 8, column (a) shows the original images with the initial contour, column (b) shows the ground truths, columns (c), (d), and (e) show segmentation results of the LBF [25], Lankton and Tannenbaum [28], and proposed method, respectively. It is evident from the results that the proposed method yields better segmentation than both Lankton and Tannenbaum and LBF methods.

Epicardium Segmentation.
The proposed method is also used to segment epicardial contours of both ventricles (LV and RV). In [13], Liu et al. adapted the anatomical knowledge of the heart in his model and proposed a two-step method, where both epicardial and endocardial contours are represented by a single level set function. Inspired by the work of Liu et al. [13], the proposed method can also segment epicardial and endocardial contours at the same time but, instead of two steps, it can properly delineate both contours (ED and EP) in single step level set evolution. Figure 9 explains the complete segmentation process. The contour is initialized with two circular level sets manually as shown in Figure 9(a), one inside ventricles and one near the endocardial boundary. After the evolution process, the final contours (epicardial and endocardial) will capture both RV and LV as shown in Figure 9(b). Furthermore, more experiments are also performed on single ventricle in order to extract epicardial and endocardial contours. Figure 10 illustrates the process of single ventricle segmentation, where Figure 10(a) shows the initial contours, Figure 10(b) shows the ground truth, and Figure 10(c) shows the final segmentation result using proposed method.

Quantitative Analysis. Dice coefficient (DSC) and
Hausdorff-Distance (HD) metrics are used for the quantitative analysis and comparison with the state-of-the-art methods. DSC measures how well segmentation overlaps the ground truth . Segmentation results have high accuracy when DSC value is close to 1. DSC is defined as where Ω is the ground truth region and Ω is the segmented region. HD is the second statistical measure used for the quantitative evaluation in this paper. It is the distance between the segmented region and the ground truth contour. Segmentation results have high accuracy when the HD value is close to 0. HD is defined as where the ground truth and the segmented region contain a group of points = { 1 , 2 , 3 , . . . , } and = { 1 , 2 , 3 , . . . , }, respectively, and is the distance from to the nearest point on the region . Figure 11 shows DSC and HD value of LV cardiac MRI segmentation based on Figure 7. Similarly, Figure 12 shows the DSC and HD value of RV cardiac MRI segmentation based on Figure 8. It shows that proposed method gets minimum HD values and highest DSC value compared to other intensity-based methods.
The segmentation accuracy of the proposed method is being validated with cardiac application based methods using DSC and HD metrics as shown in Table 1 using the training dataset from MICCAI 2009 [36]. Table 2 shows a segmentation accuracy comparison between the proposed method and the state of the art using DSC metrics on validation dataset from MICCAI 2009 [36]. It shows that the proposed method yields highest DSC value among all methods, which means the segmentation results of proposed method are closest to their respective ground truths. In Table 3, cardiac MRI MICCAI 2012 [37] test dataset is used to compute DSC and HD metrics to compare the accuracy of the proposed method with the state of the art. It shows that the proposed method yields better segmentation accuracy compared to previous methods.

Selection of Parameters and .
Parameter plays an essential role in the proposed segmentation method. This parameter deals with the measure of local and global intensity force during the segmentation process. The local SPF force is scaled with a (1 − ) parameter and the global SPF force is scaled with parameter, where 0 ≤ ≤ 1. When the input image is affected by the intensity inhomogeneity, the value decided for should be near 0 to diminish the interference of global SPF. Similarly, if the image has uniform intensity distribution and it is affected by the noise, then the selected value of should be near 1 to make global SPF force dominant.
The parameter > 0 plays an important role in segmenting local structure of an image. If the image has one or more objects with complex structures, then a small value of is chosen. In turn, for images with large objects and fairly simple structure, big value of is chosen. In this paper, for the images which have the objects with complex structures, the chosen value of is between 4 and 8. In turn, for large objects, a value between 12 and 20 is chosen.

Conclusion
In this paper, a new region-based active contour method is presented to segment both left and right ventricles in  Finally, a Gaussian kernel is used to regularize the level set curve after each time step, which also eliminates the need of expensive reinitialization. Experimental results demonstrate that the proposed method can efficiently segment LV and RV separately or together in cardiac MRI. Furthermore, a quantitative comparison of the proposed method with the state-of-the-art methods demonstrates that the proposed method yields better segmentation accuracy.

Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.