**Automated Segmentation of Real-Time 3D Echocardiography Using an Incompressibility Constraint**

Yun Zhu, Xenophon Papademetris, Albert Sinusas and James Duncan *Yale University USA* 

#### **1. Introduction**

Coronary Heart Disease (CHD) is characterized by reduced blood flow and oxygen supply to the heart muscle, resulting from the occlusion of one or more major coronary arteries. CHD remains the most prevalent cause of mortality in developed countries and represents one of the major burdens on the healthcare systems today. Approximately every 25 seconds, an American will suffer from a coronary event, and about every minute, someone will die from one. According to the 2009 American Heart Association (AHS) report (Lloyd-Jones et al., 2009}, an estimated of 16.8 million American adults have CHD (extrapolated to US population in 2009 from National Health and Nutrition Examination Survey (NHANES) 2005-2006) and about 20% of total deaths in the United States are caused by CHD.

The advancements in noninvasive imaging techniques, such as real-time three-dimensional (RT3D) echocardiography, have enabled physicians to detect CHD in its earliest and most treatable stage. With cardiac imaging, physicians are able to evaluate essential global and local functional parameters, such as ejection fraction (EF), wall thickening, strain/strain rate, and etc. The rapid progress in cardiac imaging, however, has led to new challenges in handling of huge amount of image data involved in comprehensive functional patient studies. Manually analyzing these data sets becomes a formidable task for cardiologists, radiologists, and technicians in order to interpret the data and derive clinically useful information for diagnosis or decision support for surgical and pharmacological interventions. Also, manual analysis is subjective and therefore compromises the accuracy and reproducibility of quantitative measurements.

The abovementioned reasons have triggered a great demand for computerized techniques to automate the analysis of cardiac images. Various image-processing tasks need to be performed in order to recover diagnostically useful information, among which myocardial segmentation is one of the most important tasks. Myocardial segmentation aims to delineate the endocardial (ENDO) and epicardial (EPI) boundaries from cardiac images. Accurate segmentation of myocardial boundaries is essential for deriving cardiac global functional parameters such as ventricular mass/volume, ejection fraction, and wall thickening. It is also a prerequisite step for accurate myocardial deformation analysis.

Automated Segmentation of Real-Time

3D Echocardiography Using an Incompressibility Constraint 87

Fig. 1(b) shows a short-axis cross-section of the heart. The grey region represents the myocardium. The inner surface of the myocardium is called endocardium, and outer surface is called epicardium. As the LV is a powerful blood pump for the systemic circulation, the delineation of the LV endocardium and epicardium is often the object of interest for cardiologists, particularly because it is an important step for analyzing LV anatomical

Automatic LV segmentation algorithms can be broadly categorized into region-based and boundary-based methods. Region-based methods exploit the homogeneity of spatially dense information, e.g. pixel-wise grey level values, to produce segmented images. Examples of region-based segmentation methods include thresholding, region-growing, Markov random field-based approaches, and etc. In contrast to the region-based methods, boundary-based methods rely on the pixel-wise difference to guide the process of segmentation. They try to locate points of abrupt changes in the grey tone images. Examples of boundary-based methods include edge-detection, Hough transform, boundary-based snakes, and etc. However, in numerous medical imaging modalities, the LV boundaries cannot be reliably and accurately detected using the algorithms which rely only on boundary or edge information. Reasons for this include significant signal loss, noise, complicated geometric shapes. These problems are ever present for acquired in echocardiography, where the boundary detection problem is further complicated by the

In an effort to overcome these difficulties, various sources of image information have been used in the segmentation process. Most commonly used ones include grey level distribution, gradient, texture, and phase information. The choice of which information to use depends

Incorporation of imaging physics as prior knowledge has been proven to be useful for cardiac segmentation (Paragios, 2002; Pluempitiwiriyawe et al., 2005, Chen et al., 2008). The random intensity distribution in ultrasound images, also known as speckle, is caused by the interference of energy from randomly distributed scatters, too small to be resolved by the imaging system. When an acoustic pulse travels through tissue or any medium, backscattering from the scatters in the range cell contributes to the returned echo. This contribution to the echo from the scatters in the range cell has been treated as a random walk because the locations of the scatters are considered to be random. The backscattered echo, therefore, is complex valued with a real part *X* and imagery part *Y* . By virtue of

2 2 *I XY*

, where <sup>2</sup>

central limit theorem, *X* and *Y* have identical Gaussian distribution <sup>2</sup> *N* 0,

structure, quantifying cardiac function, and estimating LV motion.

presence of attenuation, speckle, shadows, confusing anatomical structures.

on imaging modality, image quality, and specific applications.

The envelop of the backscattered echo, *I* , is given by

It has been shown that the envelop *I* has a Rayleigh distribution

**2.1 Grey level distribution** 

is its variance.

However, myocardial segmentation is challenging, particularly for the EPI surface. This is due to the misleading, physically corrupted, and sometimes incomplete visual information of the EPI contours.


In this chapter, we present an automated coupled deformable model for the segmentation of LV myocardial borders from RT3D echocardiographic images. This approach was originally proposed in (Zhu et al., 2010). It incorporates the incompressibility property of myocardial, and therefore is able to handle fuzzy EPI borders. This model is formulated in a Bayesian framework, which maximizes the regional homogeneities of a cardiac image, while maintaining the myocardial volume during a cycle. By simultaneously evolving both ENDO and EPI surfaces, an automatic segmentation of the full myocardium is achieved from RT3D echocardiographic images.

The remaining sections are organized as follows. Section 2 presents background of echocardiographic image segmentation. Section 3 discusses the incompressibility property of the myocardium. Section 4 is the general framework and details of the method. Section 5 presents experimental results and finally section 6 is the conclusion.

#### **2. Background**

The heart is composed of four chambers: left atrium (LA), left ventricle, right atrium (RA), and right ventricle. Fig. 1(a) shows a long-axis cross-section of the heart. As shown in Fig. 1(a), the atria and ventricles are surrounded by muscle tissues called myocardium. Contraction and relaxation of the muscle fibers in the myocardium cause the pumping of the heart.

Fig. 1. Heart anatomy. (a) long-axis cross-section, (b) short-axis cross-section. The grey region represents the myocardium.

However, myocardial segmentation is challenging, particularly for the EPI surface. This is due to the misleading, physically corrupted, and sometimes incomplete visual information

 Myocardium and liver are anatomically close to each other and often share similar intensity values. Thus, there is usually no clear and visible edge between them; The myocardium/lung air contrast is lower than the blood pool/myocardium contrast. This implies that myocardium and lung air have similar intensity properties, and might

 The right ventricle (RV) myocardium, which has similar intensity values as the left ventricular (LV) myocardium, often merges into the LV myocardium at the LV/RV

In this chapter, we present an automated coupled deformable model for the segmentation of LV myocardial borders from RT3D echocardiographic images. This approach was originally proposed in (Zhu et al., 2010). It incorporates the incompressibility property of myocardial, and therefore is able to handle fuzzy EPI borders. This model is formulated in a Bayesian framework, which maximizes the regional homogeneities of a cardiac image, while maintaining the myocardial volume during a cycle. By simultaneously evolving both ENDO and EPI surfaces, an automatic segmentation of the full myocardium is achieved from RT3D

The remaining sections are organized as follows. Section 2 presents background of echocardiographic image segmentation. Section 3 discusses the incompressibility property of the myocardium. Section 4 is the general framework and details of the method. Section 5

The heart is composed of four chambers: left atrium (LA), left ventricle, right atrium (RA), and right ventricle. Fig. 1(a) shows a long-axis cross-section of the heart. As shown in Fig. 1(a), the atria and ventricles are surrounded by muscle tissues called myocardium. Contraction and

relaxation of the muscle fibers in the myocardium cause the pumping of the heart.

Fig. 1. Heart anatomy. (a) long-axis cross-section, (b) short-axis cross-section. The grey

presents experimental results and finally section 6 is the conclusion.

**(a) (b)** 

be misclassified as belonging in the same class;

of the EPI contours.

juncture.

echocardiographic images.

region represents the myocardium.

**2. Background** 

Fig. 1(b) shows a short-axis cross-section of the heart. The grey region represents the myocardium. The inner surface of the myocardium is called endocardium, and outer surface is called epicardium. As the LV is a powerful blood pump for the systemic circulation, the delineation of the LV endocardium and epicardium is often the object of interest for cardiologists, particularly because it is an important step for analyzing LV anatomical structure, quantifying cardiac function, and estimating LV motion.

Automatic LV segmentation algorithms can be broadly categorized into region-based and boundary-based methods. Region-based methods exploit the homogeneity of spatially dense information, e.g. pixel-wise grey level values, to produce segmented images. Examples of region-based segmentation methods include thresholding, region-growing, Markov random field-based approaches, and etc. In contrast to the region-based methods, boundary-based methods rely on the pixel-wise difference to guide the process of segmentation. They try to locate points of abrupt changes in the grey tone images. Examples of boundary-based methods include edge-detection, Hough transform, boundary-based snakes, and etc. However, in numerous medical imaging modalities, the LV boundaries cannot be reliably and accurately detected using the algorithms which rely only on boundary or edge information. Reasons for this include significant signal loss, noise, complicated geometric shapes. These problems are ever present for acquired in echocardiography, where the boundary detection problem is further complicated by the presence of attenuation, speckle, shadows, confusing anatomical structures.

In an effort to overcome these difficulties, various sources of image information have been used in the segmentation process. Most commonly used ones include grey level distribution, gradient, texture, and phase information. The choice of which information to use depends on imaging modality, image quality, and specific applications.

#### **2.1 Grey level distribution**

Incorporation of imaging physics as prior knowledge has been proven to be useful for cardiac segmentation (Paragios, 2002; Pluempitiwiriyawe et al., 2005, Chen et al., 2008). The random intensity distribution in ultrasound images, also known as speckle, is caused by the interference of energy from randomly distributed scatters, too small to be resolved by the imaging system. When an acoustic pulse travels through tissue or any medium, backscattering from the scatters in the range cell contributes to the returned echo. This contribution to the echo from the scatters in the range cell has been treated as a random walk because the locations of the scatters are considered to be random. The backscattered echo, therefore, is complex valued with a real part *X* and imagery part *Y* . By virtue of central limit theorem, *X* and *Y* have identical Gaussian distribution <sup>2</sup> *N* 0, , where <sup>2</sup> is its variance.

The envelop of the backscattered echo, *I* , is given by

$$I = \sqrt{X^2 + Y^2}$$

It has been shown that the envelop *I* has a Rayleigh distribution

Automated Segmentation of Real-Time

where *E* is the expectation.

images with large dropout (Qian et al., 2006).

and

2005).

**2.2 Gradient** 

**2.3 Texture** 

models.

3D Echocardiography Using an Incompressibility Constraint 89

<sup>2</sup> *E I*

Davignon *et al.* used a Nakagmi distribution to segment ultrasound images (Davignon et al.,

Apart from the theoretical models mentioned above, empirical models have also been used in ultrasound image segmentation. For example, Tao *et al.* modeled the ultrasound speckle with a Gamma distribution and incorporated it into a tunneling descent optimization framework to overcome local minima (Tao and Tagare, 2007). Xiao *et al.* used a log-normal distribution for modeling speckle in breast images (Xiao et al., 2002). Qian *et al.* incorporated a log-normal distribution into a level-set framework to segment rat echocardiographic

Intensity gradient looks for intensity discontinuities between subregions corresponding to different tissue types. Intensity gradient can be computed from an image with standard differential operators (e.g. Sobel operator). A voxel is considered to be a boundary voxel if its intensity gradient is above a threshold. Gradient-based segmentation, such as such as gradient-based snake (Kass et al., 1987), Active Shape Model (ASM) (Cootes et al., 1995), and level-set approaches (Malladi et al., 1995), has been extensively used in cardiac applications (Lynch et al., 2006; Lynch et al., 2008; Assen et al., 2008; Chen et al., 2002; Corsi et al., 2002). One limitation of intensity gradient is that it may suffer from spurious, missing, and discontinuous edges. This is especially evident in EPI segmentation because there are usually no clear and visible edges between the myocardium and background. In addition, intensity gradient is suboptimal for ultrasound images because a boundary response in ultrasound images is anisotropic and depends on the relative orientation of the transducer to the boundary. Therefore, intensity gradient information is often used in conjunction with

Although no formal definition of texture exists, it can be loosely defined as one or more local patterns that are repeated in a periodic manner. Texture analysis has been attempted in numerous ways in medical image analysis, particularly for ultrasound images. The three principal approaches to describe texture are image pyramids, random fields, and statistical

The idea behind image pyramids is to generate a number of homogeneous parameters that represent the response of a bank of filters at different scales and possibly different orientations. Different filters have been proposed including Gabor filters (Xie et al., 2005; Zhan and Shen, 2006), Gaussian derivatives (Chen et al., 1998) and wavelet transforms (Mojsilovic et al., 1998; Yoshida et al., 2003). The idea behind random fields is that the value

other image features (Pluempitiwiriyawej et al., 2005; Chen et al, 2003).

$$P\left(I; \sigma\right) = \frac{I}{\sigma^2} \exp\left(-\frac{I^2}{2\sigma^2}\right),$$

The Rayleigh distribution has been proven to be one of the most popular models in ultrasound image analysis. For example, a Rayleigh distribution was incorporated into a region-based level-set approach for the segmentation of ultrasound images (Sarti et al., 2005). Also Cohen *et al.* proposed to use Rayleigh distribution as an image formation model in a maximum likelihood motion estimation scheme for noisy ultrasound images (Cohen and Dinstein 2002). Steen *et al.* proposed to use Rayleigh distribution for anisotropic diffusion-based edge detection (Steen and Olstad 1994).

While Rayleigh distribution is extensively used in ultrasound image processing, it is not always not norm because it is valid only in the special case of a large number of randomly distributed scatters. In this condition, the speckle is called fully developed. It has been shown that the echo envelop has a Rayleigh distribution when signal noise ratio (SNR)is approximately 1.91. Deviations from such special scattering conditions are called pre-Rayleigh condition when SNR <1.91, and post-Rayleigh condition when SNR>1.91 (Tuthill et al., 1988). The K-distribution has been proposed to model the pre-Rayleigh condition by accounting for low effective scatter density (Shankar et al., 1993; Shankar, 1995), and Rice distribution has been proposed to model the post-Rayleigh condition by accounting for a coherence component due to the presence of a regular structure of scatters within the tissue (Insana et al., 1986; Tuthill et al., 1988). However, the Rician family fails to model the pre-Rayleigh condition, and K-distribution model does not take into account the post-Rayleigh condition. Generalized K-distribution (Jakeman, 1999), homodyned K-distribution (Dutt and Greenleaf, 1994), and Rician Inverse of Gaussian distribution (Eltoft, 2003) have been proposed as general models to encompass pre-Rayleigh, Rayleigh, and post-Rayleigh conditions. Unfortunately, the complex nature of these models limited their practical applications.

Recently, Nakagami distribution has been proposed as a simple generalized model to collectively represent pre-Rayleigh, Rayleigh, and post-Rayleigh distributions (Shankar, 2000). It is a two-parameter distribution with analytical simplicity, which makes it relatively easy to estimate it parameters. The probability density function (pdf) of the Nakagami distribution is given by

$$P\left(I; \mu, \Omega\right) = \frac{2\,\mu^{\mu}I^{2\mu-1}}{\Gamma\left(\mu\right)\Omega^{\mu}} \exp\left(-\frac{\mu}{\Omega}I^{2}\right),$$

where is the Gamma function, is the shape parameter, and is the scaling parameter. For 1 the pdf reduces to Rayleigh condition. For 1 , the pdf can be described as post-Rayleigh, which is similar to Rician distribution, while for 1 , the pdf can be described as pre-Rayleigh.

As shown in (Shankar, 2000), the shape parameter and scaling parameters can be obtained from the moments of the envelop as follows

$$\mu = \frac{\left[E\left(I^2\right)\right]^2}{E\left[I^2 - E\left(I^2\right)\right]^2}$$

and

88 Echocardiography – New Techniques

2 2 ; exp <sup>2</sup> *I I P I*

The Rayleigh distribution has been proven to be one of the most popular models in ultrasound image analysis. For example, a Rayleigh distribution was incorporated into a region-based level-set approach for the segmentation of ultrasound images (Sarti et al., 2005). Also Cohen *et al.* proposed to use Rayleigh distribution as an image formation model in a maximum likelihood motion estimation scheme for noisy ultrasound images (Cohen and Dinstein 2002). Steen *et al.* proposed to use Rayleigh distribution for anisotropic

While Rayleigh distribution is extensively used in ultrasound image processing, it is not always not norm because it is valid only in the special case of a large number of randomly distributed scatters. In this condition, the speckle is called fully developed. It has been shown that the echo envelop has a Rayleigh distribution when signal noise ratio (SNR)is approximately 1.91. Deviations from such special scattering conditions are called pre-Rayleigh condition when SNR <1.91, and post-Rayleigh condition when SNR>1.91 (Tuthill et al., 1988). The K-distribution has been proposed to model the pre-Rayleigh condition by accounting for low effective scatter density (Shankar et al., 1993; Shankar, 1995), and Rice distribution has been proposed to model the post-Rayleigh condition by accounting for a coherence component due to the presence of a regular structure of scatters within the tissue (Insana et al., 1986; Tuthill et al., 1988). However, the Rician family fails to model the pre-Rayleigh condition, and K-distribution model does not take into account the post-Rayleigh condition. Generalized K-distribution (Jakeman, 1999), homodyned K-distribution (Dutt and Greenleaf, 1994), and Rician Inverse of Gaussian distribution (Eltoft, 2003) have been proposed as general models to encompass pre-Rayleigh, Rayleigh, and post-Rayleigh conditions.

Unfortunately, the complex nature of these models limited their practical applications.

described as post-Rayleigh, which is similar to Rician distribution, while for

1 the pdf reduces to Rayleigh condition. For

Recently, Nakagami distribution has been proposed as a simple generalized model to collectively represent pre-Rayleigh, Rayleigh, and post-Rayleigh distributions (Shankar, 2000). It is a two-parameter distribution with analytical simplicity, which makes it relatively easy to estimate it parameters. The probability density function (pdf) of the Nakagami

> 2 1 <sup>2</sup> <sup>2</sup> ; , exp *<sup>I</sup> P I <sup>I</sup>*

As shown in (Shankar, 2000), the shape parameter and scaling parameters can be obtained

*E I*

*EI EI*

 

<sup>2</sup> <sup>2</sup>

<sup>2</sup> 2 2

 

is the shape parameter, and is the scaling

1 , the pdf can be

1 , the pdf

2

 

diffusion-based edge detection (Steen and Olstad 1994).

distribution is given by

parameter. For

where is the Gamma function,

from the moments of the envelop as follows

can be described as pre-Rayleigh.

<sup>2</sup> *E I*

where *E* is the expectation.

Davignon *et al.* used a Nakagmi distribution to segment ultrasound images (Davignon et al., 2005).

Apart from the theoretical models mentioned above, empirical models have also been used in ultrasound image segmentation. For example, Tao *et al.* modeled the ultrasound speckle with a Gamma distribution and incorporated it into a tunneling descent optimization framework to overcome local minima (Tao and Tagare, 2007). Xiao *et al.* used a log-normal distribution for modeling speckle in breast images (Xiao et al., 2002). Qian *et al.* incorporated a log-normal distribution into a level-set framework to segment rat echocardiographic images with large dropout (Qian et al., 2006).

#### **2.2 Gradient**

Intensity gradient looks for intensity discontinuities between subregions corresponding to different tissue types. Intensity gradient can be computed from an image with standard differential operators (e.g. Sobel operator). A voxel is considered to be a boundary voxel if its intensity gradient is above a threshold. Gradient-based segmentation, such as such as gradient-based snake (Kass et al., 1987), Active Shape Model (ASM) (Cootes et al., 1995), and level-set approaches (Malladi et al., 1995), has been extensively used in cardiac applications (Lynch et al., 2006; Lynch et al., 2008; Assen et al., 2008; Chen et al., 2002; Corsi et al., 2002).

One limitation of intensity gradient is that it may suffer from spurious, missing, and discontinuous edges. This is especially evident in EPI segmentation because there are usually no clear and visible edges between the myocardium and background. In addition, intensity gradient is suboptimal for ultrasound images because a boundary response in ultrasound images is anisotropic and depends on the relative orientation of the transducer to the boundary. Therefore, intensity gradient information is often used in conjunction with other image features (Pluempitiwiriyawej et al., 2005; Chen et al, 2003).

#### **2.3 Texture**

Although no formal definition of texture exists, it can be loosely defined as one or more local patterns that are repeated in a periodic manner. Texture analysis has been attempted in numerous ways in medical image analysis, particularly for ultrasound images. The three principal approaches to describe texture are image pyramids, random fields, and statistical models.

The idea behind image pyramids is to generate a number of homogeneous parameters that represent the response of a bank of filters at different scales and possibly different orientations. Different filters have been proposed including Gabor filters (Xie et al., 2005; Zhan and Shen, 2006), Gaussian derivatives (Chen et al., 1998) and wavelet transforms (Mojsilovic et al., 1998; Yoshida et al., 2003). The idea behind random fields is that the value

Automated Segmentation of Real-Time

as will be detailed in Section 4.3.

constraint can be expressed as follows

**4.1 General framework** 

deterministic constraint.

distribution from RT3D echocardiographic images.

**4.2 Data adherence** 

as shown in Fig. 2.

**4. Method** 

3D Echocardiography Using an Incompressibility Constraint 91

canine Dynamic Spatial Reconstructor data to study the changes in MV, obtaining a relatively constant volume that was consistent with Hamilton's findings (Hoffman and Ritman, 1987; Hoffman and Ritman, 1985). More recent work done by King *et al.*, who used freehand three dimensional (3-D) echocardiography to measure MV and mass, showed a 1.1% difference of volume between end-diastole(ED) and end-systole (ES) (King et al., 2002). Bowaman *et al.* found a variation of around 5% from ED to ES by using highresolution magnetic resonance (MR) imaging (Bowman and Kovacs, 2003). O'Donnell also analyzed the myocardial volume using MR imaging, and found the difference of ED and ES to be around 2.5% (O'Donnell and Funka-Lea, 1998). The common conclusion of these studies is that the myocardial is nearly incompressible and the variation is less than 5% during a cardiac cycle. This incompressibility property of the myocardium is used an as important cardiac structure constraint which is taken into consideration in our approach,

In this section, we formulate our segmentation problem in a maximum a posterior (MAP) framework combining image-derived information and the incompressibility property of the myocardium. Let *I* be a 3-D cardiac image, *S*in be the ENDO surface, and *S*out be the EPI surface. A MAP framework that combines image information and incompressibility

ˆ ˆ, arg max , | arg max | , ,

in out in out in out in out

Equation (1) is a probability function which adheres to image data, modulated by the prior knowledge of incompressibility property of the myocardium. In particular, *PIS S* | , in out is the probability of producing an image *I* given *S*in and *S*out by assuming the piecewise homogeneities of each region enclosed by *S*in and *S*out . *PS S* in out , is the incompressibility constraint which keeps the MV enclosed by *S*in and *S*out nearly constant during a cardiac cycle. Since the myocardium is only nearly incompressible, it is more reasonable to define it within a probabilistic framework rather than imposing a

Region-based deformable models have been successfully applied in the segmentation of images with weak boundaries, due to their robustness (Chan and Vese, 2001). In this work, we evolve a three-phase region-based deformable model based on the statistical intensity

To evolve a region-based model, we first need to determine the intensity distribution of each region within a cardiac image. It is obvious that the entire cardiac image is partitioned by *S*in and *S*out into three regions, namely, LV blood pool, LV myocardium, and background,

, , data adherence incompressibility constraint

*S S PS S I PI S S PS S* (1)

in out in out

*S S S S*

at each pixel/voxel is chosen by two-dimensional (2-D)/ three dimensional (3-D) stochastic process. Given a type of pdf of this stochastic process, one can estimate the value at a particular pixel/voxel given the values of other pixels/voxles in its neighborhood. The most commonly used random field is Markov Random Field (MRF) (Bouhlel et al., 2004). Statistical models analyze the spatial distribution of grey level values, by computing local features at each point in the image, and deriving a set of statistics from the distributions of the local features. Statistical approaches yield characterizations of textures as smooth, coarse, grainy, etc. Major statistical texture descriptors include co-occurrence matrix (Nicholas et al., 1986, Sahiner et al., 2004, Kuo et al., 2001), auto-correlation (Chen et al., 2000), edge frequency, run length, Law's texture energy, and fractal texture description. Several researchers have proposed to extract fractal texture features (Wu et al., 1992; Lee et al., 2003; Chen et al., 2005). In addition, several attempts have been made to combine multiple texture measures to improve discrimination abilities (Hao et al, 2001; Christodoulou et al., 2003; Stoitsis et al., 2006).

#### **2.4 Phase**

The local phase provides an alternative to intensity gradient to characterize structures in an image (Kovesi, 1996; Boukerroui et al., 2001; Mulet-Parada and Noble, 2000; Hacihaliloglu et al., 2008). Phase-based methods postulate that feature information is encoded at points where phase congruency is maximized, i.e. all the Fourier components are in phase. Generally, phase is estimated by quadrature filter bank. Thus, there is a link between phasebased methods and other wavelet methods. Phase-based methods are suggested to be more robust than intensity gradient for ultrasound images because intensity gradients in ultrasound images depend on the relative orientation of the transducer to the boundary. In addition, the presence of speckles and imaging artifacts might cause the variation of intensity gradient of equally significant features in the data. One limitation of phase-based methods, however, is that speckle also has its own phase signatures, and therefore an appropriate spatial scale has to be selected.

#### **3. Incompressibility of myocardium**

The heart is a remarkably efficient and during mechanical pump composed of complex biological materials. The main structural elements of the myocardium are inter-connected networks of muscle fibers and collagen fibers, as well as matrix that embeds them. The fibers are generally tangential to the ENDO and EPI surfaces, following a path of a righthanded helical geometry. The interstitial fluid carries only hydrostatic pressure, which in turn, is affected by the length and configuration changes of the fibers. These cause pressure gradients, which may result in the flow of the matrix. However, the fluid within the tissue is negligible for the duration of a cardiac cycle because the myocardium as low permeability. Consequently, the myocardium can be assumed to be nearly incompressible (Glass et al., 1990).

A few independent studies have been performed to quantitatively analyze the change of myocardial volume (MV) over an entire cardiac cycle. For example, Hamilton *et al.* performed experiments on frogs, turtles, and dogs, and discovered a relatively consistency of the MV during a cardiac cycle (Hamilton et al., 1932). Hoffman *et al.* used canine Dynamic Spatial Reconstructor data to study the changes in MV, obtaining a relatively constant volume that was consistent with Hamilton's findings (Hoffman and Ritman, 1987; Hoffman and Ritman, 1985). More recent work done by King *et al.*, who used freehand three dimensional (3-D) echocardiography to measure MV and mass, showed a 1.1% difference of volume between end-diastole(ED) and end-systole (ES) (King et al., 2002). Bowaman *et al.* found a variation of around 5% from ED to ES by using highresolution magnetic resonance (MR) imaging (Bowman and Kovacs, 2003). O'Donnell also analyzed the myocardial volume using MR imaging, and found the difference of ED and ES to be around 2.5% (O'Donnell and Funka-Lea, 1998). The common conclusion of these studies is that the myocardial is nearly incompressible and the variation is less than 5% during a cardiac cycle. This incompressibility property of the myocardium is used an as important cardiac structure constraint which is taken into consideration in our approach, as will be detailed in Section 4.3.

#### **4. Method**

90 Echocardiography – New Techniques

at each pixel/voxel is chosen by two-dimensional (2-D)/ three dimensional (3-D) stochastic process. Given a type of pdf of this stochastic process, one can estimate the value at a particular pixel/voxel given the values of other pixels/voxles in its neighborhood. The most commonly used random field is Markov Random Field (MRF) (Bouhlel et al., 2004). Statistical models analyze the spatial distribution of grey level values, by computing local features at each point in the image, and deriving a set of statistics from the distributions of the local features. Statistical approaches yield characterizations of textures as smooth, coarse, grainy, etc. Major statistical texture descriptors include co-occurrence matrix (Nicholas et al., 1986, Sahiner et al., 2004, Kuo et al., 2001), auto-correlation (Chen et al., 2000), edge frequency, run length, Law's texture energy, and fractal texture description. Several researchers have proposed to extract fractal texture features (Wu et al., 1992; Lee et al., 2003; Chen et al., 2005). In addition, several attempts have been made to combine multiple texture measures to improve discrimination abilities (Hao et al, 2001;

The local phase provides an alternative to intensity gradient to characterize structures in an image (Kovesi, 1996; Boukerroui et al., 2001; Mulet-Parada and Noble, 2000; Hacihaliloglu et al., 2008). Phase-based methods postulate that feature information is encoded at points where phase congruency is maximized, i.e. all the Fourier components are in phase. Generally, phase is estimated by quadrature filter bank. Thus, there is a link between phasebased methods and other wavelet methods. Phase-based methods are suggested to be more robust than intensity gradient for ultrasound images because intensity gradients in ultrasound images depend on the relative orientation of the transducer to the boundary. In addition, the presence of speckles and imaging artifacts might cause the variation of intensity gradient of equally significant features in the data. One limitation of phase-based methods, however, is that speckle also has its own phase signatures, and therefore an

The heart is a remarkably efficient and during mechanical pump composed of complex biological materials. The main structural elements of the myocardium are inter-connected networks of muscle fibers and collagen fibers, as well as matrix that embeds them. The fibers are generally tangential to the ENDO and EPI surfaces, following a path of a righthanded helical geometry. The interstitial fluid carries only hydrostatic pressure, which in turn, is affected by the length and configuration changes of the fibers. These cause pressure gradients, which may result in the flow of the matrix. However, the fluid within the tissue is negligible for the duration of a cardiac cycle because the myocardium as low permeability. Consequently, the myocardium can be assumed to be nearly incompressible (Glass et al.,

A few independent studies have been performed to quantitatively analyze the change of myocardial volume (MV) over an entire cardiac cycle. For example, Hamilton *et al.* performed experiments on frogs, turtles, and dogs, and discovered a relatively consistency of the MV during a cardiac cycle (Hamilton et al., 1932). Hoffman *et al.* used

Christodoulou et al., 2003; Stoitsis et al., 2006).

appropriate spatial scale has to be selected.

**3. Incompressibility of myocardium** 

**2.4 Phase** 

1990).

#### **4.1 General framework**

In this section, we formulate our segmentation problem in a maximum a posterior (MAP) framework combining image-derived information and the incompressibility property of the myocardium. Let *I* be a 3-D cardiac image, *S*in be the ENDO surface, and *S*out be the EPI surface. A MAP framework that combines image information and incompressibility constraint can be expressed as follows

$$P\left(\hat{S}\_{\text{in}}, \hat{S}\_{\text{out}}\right) = \underset{S\_{\text{in}}, S\_{\text{out}}}{\arg\max} P\left(S\_{\text{in}}, S\_{\text{out}} \mid I\right) = \underset{S\_{\text{in}}, S\_{\text{out}}}{\arg\max} \underbrace{P\left(I \mid S\_{\text{in}}, S\_{\text{out}}\right)}\_{\text{data\text{-}difference}} \underbrace{P\left(S\_{\text{in}}, S\_{\text{out}}\right)}\_{\text{incompressible}} \tag{1}$$

Equation (1) is a probability function which adheres to image data, modulated by the prior knowledge of incompressibility property of the myocardium. In particular, *PIS S* | , in out is the probability of producing an image *I* given *S*in and *S*out by assuming the piecewise homogeneities of each region enclosed by *S*in and *S*out . *PS S* in out , is the incompressibility constraint which keeps the MV enclosed by *S*in and *S*out nearly constant during a cardiac cycle. Since the myocardium is only nearly incompressible, it is more reasonable to define it within a probabilistic framework rather than imposing a deterministic constraint.

#### **4.2 Data adherence**

Region-based deformable models have been successfully applied in the segmentation of images with weak boundaries, due to their robustness (Chan and Vese, 2001). In this work, we evolve a three-phase region-based deformable model based on the statistical intensity distribution from RT3D echocardiographic images.

To evolve a region-based model, we first need to determine the intensity distribution of each region within a cardiac image. It is obvious that the entire cardiac image is partitioned by *S*in and *S*out into three regions, namely, LV blood pool, LV myocardium, and background, as shown in Fig. 2.

Automated Segmentation of Real-Time

3D Echocardiography Using an Incompressibility Constraint 93

Fig. 3. The histograms of echocardiographic images with their fitted distributions.

Let be a bounded open set of <sup>3</sup> *R* , and be partitioned by *S*in and *S*out into three regions, namely, LV blood pool, LV myocardium, and background, which are denoted as <sup>1</sup> , <sup>2</sup> , and <sup>3</sup> respectively. Thus, the data adherence term can be defined by a three-phase

> 3

*<sup>l</sup> l l*

<sup>0</sup> , *N V V*

**<sup>x</sup>** is the MV enclosed by the deforming ENDO and EPI surfaces. Parameter

 

0

*V V*

**<sup>x</sup>** (3)

(4)

1 log | , log ; ,

It is mentioned in Section 3 that the MV is nearly constant during cardiac cycle. Therefore,

<sup>2</sup>

in out 2 <sup>1</sup> , exp 2 2 *V V*

*V*<sup>0</sup> is the average volume, which can be calculated from the manual segmentation of the first frame. As mentioned in Section 3, the MV changes by less than 5% during a cardiac cycle.

*l PIS S PI d*

Fig. 3 shows the histogram of each region for echocardiographic images.

in out

we assume that the MV has a Gaussian distribution <sup>2</sup>

*PS S*

(a) LV blood pool, (b) LV myocardium, (c) background.

deformable model

where 2 *V d*

**4.3 Incompressibility constraint** 

Fig. 2. The short-axis view of a heart. Two dotted circles are ENDO- and EPI contours, which partitioned the entire image into three regions: LV blood pool, LV myocardium, and background. The background is inhomogeneous because it contains more than one tissues, such as RV blood pool, RV myocardium, lung air, and liver.

The LV blood pool and myocardium are homogeneous, and therefore can modeled with a single pdf. In this work, we use Nakagami distribution (Shankar, 1995) (see Section 2.1) to model the intensity distributions of LV blood pool and myocardium as follows

$$P\left(I; \mu\_{l'}, o\_l\right) = \frac{2\,\mu\_l^{\mu\_l}}{\Gamma\left(\mu\_l\right) o\_l^{\mu\_l}} I^{2\,\mu\_l - 1} \exp\left(-\frac{\mu\_l}{o\_l} I^2\right) \tag{2}$$

where *<sup>l</sup>* is the shape parameter of the Nakagami distribution, and *<sup>l</sup>* is the scaling parameter. Equation (2) describes the intensity distribution for LV blood pool when l=1, and intensity distribution for LV myocardium when l=2.

Unlike LV blood pool and myocardium, the background (see Fig. 2) because it contains more than one tissues (RV blood pool, RV myocardium, and liver), and therefore modeling it with a single distribution function would be insufficient because it contains a wide range of intensities. To handle this problem, we use a mixture model (McLachlan and Peel, 2000) to describe the intensity distribution of the background.

$$P\left(I; \mu\_{3\prime}, \alpha\_{3}\right) = \sum\_{k=1}^{M} \alpha\_{k} P\_{3,k}\left(I; \mu\_{3,k\prime}, \alpha\_{3,k}\right).$$

where <sup>3</sup> and <sup>3</sup> are the shape and scaling parameters of the component distributions. For ultrasound images, we choose 2 *M* because the background histogram has two peaks, as shown in Fig. 3 (c). The first peak corresponds to RV blood pool and lung air, while the second peak corresponds to RV myocardium and liver.

Fig. 2. The short-axis view of a heart. Two dotted circles are ENDO- and EPI contours, which partitioned the entire image into three regions: LV blood pool, LV myocardium, and background. The background is inhomogeneous because it contains more than one tissues,

model the intensity distributions of LV blood pool and myocardium as follows

*l l*

 

intensity distribution for LV myocardium when l=2.

to describe the intensity distribution of the background.

second peak corresponds to RV myocardium and liver.

*P I*

where

where

<sup>3</sup> and 

The LV blood pool and myocardium are homogeneous, and therefore can modeled with a single pdf. In this work, we use Nakagami distribution (Shankar, 1995) (see Section 2.1) to

> <sup>2</sup> 2 1 <sup>2</sup> ; , exp *<sup>l</sup> l l*

parameter. Equation (2) describes the intensity distribution for LV blood pool when l=1, and

Unlike LV blood pool and myocardium, the background (see Fig. 2) because it contains more than one tissues (RV blood pool, RV myocardium, and liver), and therefore modeling it with a single distribution function would be insufficient because it contains a wide range of intensities. To handle this problem, we use a mixture model (McLachlan and Peel, 2000)

> 3 3 3, 3, 3, 1 ; , ; , *M*

 

ultrasound images, we choose 2 *M* because the background histogram has two peaks, as shown in Fig. 3 (c). The first peak corresponds to RV blood pool and lung air, while the

*k*

 

*kk k k*

 

<sup>3</sup> are the shape and scaling parameters of the component distributions. For

*P I*

*P I I I* 

 

*<sup>l</sup>* is the shape parameter of the Nakagami distribution, and

*l l*

(2)

*<sup>l</sup>* is the scaling

*l l l*

such as RV blood pool, RV myocardium, lung air, and liver.

Fig. 3. The histograms of echocardiographic images with their fitted distributions. (a) LV blood pool, (b) LV myocardium, (c) background.

Fig. 3 shows the histogram of each region for echocardiographic images.

Let be a bounded open set of <sup>3</sup> *R* , and be partitioned by *S*in and *S*out into three regions, namely, LV blood pool, LV myocardium, and background, which are denoted as <sup>1</sup> , <sup>2</sup> , and <sup>3</sup> respectively. Thus, the data adherence term can be defined by a three-phase deformable model

$$\log P(I \mid S\_{\text{in}}, S\_{\text{out}}) = \sum\_{l=1}^{3} \int\_{\Omega\_{l}} \log P(I; \mu\_{l}, o\_{l}) \, d\mathbf{x} \tag{3}$$

#### **4.3 Incompressibility constraint**

It is mentioned in Section 3 that the MV is nearly constant during cardiac cycle. Therefore, we assume that the MV has a Gaussian distribution <sup>2</sup> <sup>0</sup> , *N V V*

$$P\left(S\_{\rm in}, S\_{\rm out}\right) = \frac{1}{\sqrt{2\pi}\sigma\_V} \exp\left(-\frac{\left(V - V\_0\right)^2}{2\sigma\_V^2}\right) \tag{4}$$

where 2 *V d* **<sup>x</sup>** is the MV enclosed by the deforming ENDO and EPI surfaces. Parameter *V*<sup>0</sup> is the average volume, which can be calculated from the manual segmentation of the first frame. As mentioned in Section 3, the MV changes by less than 5% during a cardiac cycle.

Automated Segmentation of Real-Time

**5. Experiments** 

volumetric data.

we define

**5.2 Quantitative measures** 

where *i i* , min *<sup>j</sup> <sup>j</sup>*

local similarities.

**5.1 Experimental setup** 

3D Echocardiography Using an Incompressibility Constraint 95

1 1

2 2

The ultrasound data were acquired using Philips echocardiographic system with a 4 MHz X4 xMatrix transducer. The transducer consists of 3000 miniaturized piezoelectric elements and offers steering in both azimuth and elevation of the beam, permitting real-time volumetric image acquisition and rendering (Philips, 2005). In this work, we acquired 11 sequences of canine images, and each sequence consists of 20-30 frames per cardiac cycle depending on the cardiac rate. Hence, we ran our program with a total of 286 sets of

To quantify the segmentation results, we used two distance error metrics and two area error metrics, namely, mean absolute distance (MAD), Hausdorff distance (HD), percentage of

Let *A* and *B*be two surfaces from automatic and manual segmentation, respectively. Suppose they are represented by point sets, i.e. *A* **aa a** 1 2 , ,..., *<sup>N</sup>* and *B* **bb b** 1 2 , ,..., *<sup>M</sup>* ,

1 1 <sup>1</sup> MAD , , , <sup>2</sup>

HD , max max , ,max , *i j i j A B dB d A* **a b**

Volume

Volume Volume

Volume

Volume

 1 1

*N M*

*i j AB d B d A N M* **a b**

*d B* **a ab** is the distance from point *<sup>i</sup>* **a** to the closest point on surface *B*. While MAD represents the global agreement between two contours, HD compares their

> 

 

*A A B B*

*A B A*

*i j*

; , *<sup>V</sup> P I V V*

; , *<sup>V</sup> P I V V*

2 in

2 out

(8)

(9)

 in 2 2 0

 out 3 3 <sup>0</sup>

 

 

 

 

 

; , log

true positives (PTP), and percentage of false positives (PFP).

Let *<sup>A</sup>* and *<sup>B</sup>* be the region enclosed by surface *A* and *B*, we define

PTP

PFP

*P I*

; , log

*P I*

Assuming that 0 0.025*V* is the maximum deviation and using the three rule, we have the following relationship:

$$3\sigma\_V = 0.025V\_0$$

$$\sigma\_V = \frac{1}{120}V\_0$$

Therefore, the incompressibility constraint is encoded into a probability function which favors small variance of the MV while penalizing large deviations.

Combining Equations (1), (3), and (4), we arrive at the following optimization problem

$$P\left(\hat{S}\_{\text{in}}, \hat{S}\_{\text{out}}\right) = \underset{S\_{\text{in}}, S\_{\text{out}}}{\arg\max} P\left(S\_{\text{in}}, S\_{\text{out}} \mid I\right) = \underset{S\_{\text{in}}, S\_{\text{out}}}{\arg\max} \left\{ \sum\_{l=1}^{3} \int\_{\Omega\_{l}} \log P\left(I; \mu\_{l}, a\_{l}\right) d\mathbf{x} - \frac{\left(V - V\_{0}\right)^{2}}{2\sigma\_{V}^{2}} \right\} \tag{5}$$

The maximization of Equation (5) can be identified by the coupled surface evolution equations

$$\frac{\partial S\_{\rm in}}{\partial \tau} = \left( \log \left( \frac{P\left(I; \mu\_2, \alpha\_2\right)}{P\left(I; \mu\_1, \alpha\_1\right)} \right) - \frac{V - V\_0}{\sigma\_V^2} \right) \mathbf{n}\_{\rm in} \tag{6}$$

$$\frac{\partial S\_{\text{out}}}{\partial \tau} = \left( \log \left( \frac{P\left(I; \mu\_3, \rho o\_3\right)}{P\left(I; \mu\_2, \rho o\_2\right)} \right) + \frac{V - V\_0}{\sigma\_V^2} \right) \mathbf{n}\_{\text{out}} \tag{7}$$

where **n**in and **n**out are the normals of *S*in and *S*out respectively, and is the propagation time step.

To implement the evolution of Equations (6) and (7), we embed the surface *S*in and *S*out in two higher dimensional functions in **x** and out **x** , which implicitly represents *S*in and *S*out as zero level sets, i.e. *S*in in **x x** : 0 and *S*out out **x x** : 0 .

To establish the relationship between surface evolution and level set function, we use the relation **n** . Given the curve/surface evolution form

$$\frac{\partial \mathbf{S}}{\partial \mathbf{r}} = F \mathbf{n}$$

the correspondent level set evolution equation reads (Osher and Paragios, 2003)

$$\frac{\partial \varphi}{\partial \tau} = F \left| \nabla \varphi \right|$$

where is the gradient operator. Hence, the level set evolution equation for Equations (8) and (9) are

$$\frac{\partial \rho\_{\rm in}}{\partial \tau} = \left( \log \left( \frac{P\left(I; \mu\_2, \alpha\_2\right)}{P\left(I; \mu\_1, \alpha\_1\right)} \right) - \frac{V - V\_0}{\sigma\_V^2} \right) \left| \nabla \rho\_{\rm in} \right| \tag{8}$$

$$\frac{\partial \rho\_{\rm out}}{\partial \tau} = \left( \log \left( \frac{P\{I; \mu\_3, \alpha\_3\}}{P\{I; \mu\_2, \alpha\_2\}} \right) + \frac{V - V\_0}{\sigma\_V^2} \right) \left| \nabla \, \rho\_{\rm out} \right| \tag{9}$$

#### **5. Experiments**

94 Echocardiography – New Techniques

<sup>0</sup> 3 0.025

Therefore, the incompressibility constraint is encoded into a probability function which

*<sup>V</sup> V*

0 1 120

*<sup>V</sup> V*

Combining Equations (1), (3), and (4), we arrive at the following optimization problem

in out in out <sup>2</sup> , , <sup>1</sup>

The maximization of Equation (5) can be identified by the coupled surface evolution

 in 2 2 0

 out 3 3 <sup>0</sup>

To implement the evolution of Equations (6) and (7), we embed the surface *S*in and *S*out in

To establish the relationship between surface evolution and level set function, we use the

*<sup>S</sup> <sup>F</sup>* 

*F* 

 

where is the gradient operator. Hence, the level set evolution equation for Equations (8)

 **n**

the correspondent level set evolution equation reads (Osher and Paragios, 2003)

and *S*out out **x x** : 0

 

 

*S P I V V P I* 

 

*S P I V V P I*

 

1 1

2 2

; , *<sup>V</sup>*

; , *<sup>V</sup>*

*S S S S l V*

**<sup>x</sup>** (5)

2 3

*<sup>l</sup>* 2 *l l*

 

2 in

2 out

rule, we have the

0

*V V*

**n** (6)

**n** (7)

is the propagation

out **x** , which implicitly represents *S*in and

.

Assuming that 0 0.025*V* is the maximum deviation and using the three

favors small variance of the MV while penalizing large deviations.

in out in out

where **n**in and **n**out are the normals of *S*in and *S*out respectively, and

**n** . Given the curve/surface evolution form

ˆ ˆ, arg max , | arg max log ; ,

*S S PS S I PI d*

; , log

; , log

in **x** and

following relationship:

equations

time step.

relation

and (9) are

two higher dimensional functions

*S*out as zero level sets, i.e. *S*in in **x x** : 0

#### **5.1 Experimental setup**

The ultrasound data were acquired using Philips echocardiographic system with a 4 MHz X4 xMatrix transducer. The transducer consists of 3000 miniaturized piezoelectric elements and offers steering in both azimuth and elevation of the beam, permitting real-time volumetric image acquisition and rendering (Philips, 2005). In this work, we acquired 11 sequences of canine images, and each sequence consists of 20-30 frames per cardiac cycle depending on the cardiac rate. Hence, we ran our program with a total of 286 sets of volumetric data.

#### **5.2 Quantitative measures**

To quantify the segmentation results, we used two distance error metrics and two area error metrics, namely, mean absolute distance (MAD), Hausdorff distance (HD), percentage of true positives (PTP), and percentage of false positives (PFP).

Let *A* and *B*be two surfaces from automatic and manual segmentation, respectively. Suppose they are represented by point sets, i.e. *A* **aa a** 1 2 , ,..., *<sup>N</sup>* and *B* **bb b** 1 2 , ,..., *<sup>M</sup>* , we define

$$\begin{aligned} \text{MAD}\left(A,B\right) &= \frac{1}{2} \left\{ \frac{1}{N} \sum\_{i=1}^{N} d\left(\mathbf{a}\_{i\cdot\cdot} B\right) + \frac{1}{M} \sum\_{j=1}^{M} d\left(\mathbf{b}\_{j\cdot\cdot} A\right) \right\} \\\\ \text{HD}\left(A,B\right) &= \max\left\{ \max\_{i} d\left(\mathbf{a}\_{i\cdot\cdot} B\right), \max\_{j} d\left(\mathbf{b}\_{j\cdot\cdot} A\right) \right\} \end{aligned}$$

where *i i* , min *<sup>j</sup> <sup>j</sup> d B* **a ab** is the distance from point *<sup>i</sup>* **a** to the closest point on surface *B*. While MAD represents the global agreement between two contours, HD compares their local similarities.

Let *<sup>A</sup>* and *<sup>B</sup>* be the region enclosed by surface *A* and *B*, we define

$$\text{PTP} = \frac{\text{Volume} \left(\Omega\_A \cap \Omega\_B\right)}{\text{Volume} \left(\Omega\_A\right)}$$

$$\text{PFP} = \frac{\text{Volume} \left(\Omega\_A\right) - \text{Volume} \left(\Omega\_A \cap \Omega\_B\right)}{\text{Volume} \left(\Omega\_B\right)}$$

Automated Segmentation of Real-Time

92.1 88.2 4.42%

constraint is 92.2 88.1 4.65% 88.1

incompressibility constraint.

95% confidence interval = [2.4%,35.2%]).

constraint, however, is 143.2 106.1 34.97%

106.1

88.2

3D Echocardiography Using an Incompressibility Constraint 97

As explained in Section 1, another challenge in the segmentation of the EPI boundary is the presence of the LV/RV myocardium junctures. The intensity similarity between the LV and RV myocardium makes the EPI boundary ambiguous at the juncture. When the incompressibility constraint was not applied, the EPI contour leaked out to segment RV myocardium. When the incompressibility constraint was applied, however, the LV

Fig. 6 compares the MV for an entire sequence from manual segmentation, automatic segmentation with incompressibility constraint, and automatic segmentation without incompressibility constraint. We took the mean value from three experts as the MV from manual segmentation. We observed that the variation of MV for manual segmentation is

, the variation of MV for automatic segmentation with incompressibility

incompressibility constraint is much larger than that from manual segmentation. This is because the EPI contour leaked into the background, leading to the over-estimation of the MV.

Fig. 6. Comparison of MV for an example sequence from manual segmentation, automatic segmentation with incompressibility constraint, and automatic segmentation without

In addition, we performed Bland-Altman analysis (Bland and Altman, 1986) to assess the agreement of the MV measurements from manual segmentation and automatic segmentation. The Bland-Altman analysis revealed small bias and good coherence between the MV measurements from manual segmentation and automatic segmentation with incompressibility constraint. As shown in Fig. 7, the bias is -0.2%, and 95% of the computer measurements of MV can be expected to differ from expert measurements by less than 6.4% below and 5.9% above the mean. In comparison, when the incompressibility constraint was not applied, the MV largely deviated from the MV from manual segmentation (bias = 18.8%,

. The variation of MV obtained without incompressibility

. Also, we noticed that the MV obtained without

myocardium was successfully separated from the RV myocardium at the juncture.

#### **5.3 Experimental results**

In Fig. 4, we show the long axis view of the automatically segmented ENDO- and EPI surfaces at frames 2, 4, 6, and 8 during ventricular systole.

Fig. 4. Long axis view of segmented ENDO- and EPI contours at frames 2, 4, 6, and 8 during cardiac systole.

In Fig. 5, we compare the segmentation results of the EPI surface with and without incompressibility constraint. For a fair comparison, we used the same data adherence term for both cases. We observed that while the ENDO border was correctly detected even without the incompressibility constraint, the EPI contour leaked into other tissues, such as liver, that look like myocardium. This is because the LV myocardium/background contrast is lower than the LV blood pool/myocardium contrast, especially for the lateral and inferior sectors of myocardium. This low contrast obscures the exact location of the EPI boundary, making EPI segmentation more challenging than ENDO segmentation. When the incompressibility constraint was applied, however, the coupling of ENDO and EPI contours forced the evolution of the EPI contour to be consistent with that of the ENDO contour, thus preventing the leakage of the EPI contour.

Fig. 5. Comparison of EPI segmentation with and without incompressibility constraint. (a) with incompressibility constraint; (b) without incompressibility constraint. Arrow: myocardium juncture.

In Fig. 4, we show the long axis view of the automatically segmented ENDO- and EPI

Fig. 4. Long axis view of segmented ENDO- and EPI contours at frames 2, 4, 6, and 8 during

In Fig. 5, we compare the segmentation results of the EPI surface with and without incompressibility constraint. For a fair comparison, we used the same data adherence term for both cases. We observed that while the ENDO border was correctly detected even without the incompressibility constraint, the EPI contour leaked into other tissues, such as liver, that look like myocardium. This is because the LV myocardium/background contrast is lower than the LV blood pool/myocardium contrast, especially for the lateral and inferior sectors of myocardium. This low contrast obscures the exact location of the EPI boundary, making EPI segmentation more challenging than ENDO segmentation. When the incompressibility constraint was applied, however, the coupling of ENDO and EPI contours forced the evolution of the EPI contour to be consistent with that of the ENDO contour, thus

Fig. 5. Comparison of EPI segmentation with and without incompressibility constraint. (a) with incompressibility constraint; (b) without incompressibility constraint. Arrow:

**5.3 Experimental results** 

cardiac systole.

preventing the leakage of the EPI contour.

myocardium juncture.

surfaces at frames 2, 4, 6, and 8 during ventricular systole.

As explained in Section 1, another challenge in the segmentation of the EPI boundary is the presence of the LV/RV myocardium junctures. The intensity similarity between the LV and RV myocardium makes the EPI boundary ambiguous at the juncture. When the incompressibility constraint was not applied, the EPI contour leaked out to segment RV myocardium. When the incompressibility constraint was applied, however, the LV myocardium was successfully separated from the RV myocardium at the juncture.

Fig. 6 compares the MV for an entire sequence from manual segmentation, automatic segmentation with incompressibility constraint, and automatic segmentation without incompressibility constraint. We took the mean value from three experts as the MV from manual segmentation. We observed that the variation of MV for manual segmentation is

92.1 88.2 4.42% 88.2 , the variation of MV for automatic segmentation with incompressibility

constraint is 92.2 88.1 4.65% 88.1 . The variation of MV obtained without incompressibility

constraint, however, is 143.2 106.1 34.97% 106.1 . Also, we noticed that the MV obtained without

incompressibility constraint is much larger than that from manual segmentation. This is because the EPI contour leaked into the background, leading to the over-estimation of the MV.

Fig. 6. Comparison of MV for an example sequence from manual segmentation, automatic segmentation with incompressibility constraint, and automatic segmentation without incompressibility constraint.

In addition, we performed Bland-Altman analysis (Bland and Altman, 1986) to assess the agreement of the MV measurements from manual segmentation and automatic segmentation. The Bland-Altman analysis revealed small bias and good coherence between the MV measurements from manual segmentation and automatic segmentation with incompressibility constraint. As shown in Fig. 7, the bias is -0.2%, and 95% of the computer measurements of MV can be expected to differ from expert measurements by less than 6.4% below and 5.9% above the mean. In comparison, when the incompressibility constraint was not applied, the MV largely deviated from the MV from manual segmentation (bias = 18.8%, 95% confidence interval = [2.4%,35.2%]).

Automated Segmentation of Real-Time

**6. Conclusion** 

significantly improved.

**7. References** 

3D Echocardiography Using an Incompressibility Constraint 99

Table 2. Comparison of automatic outline to three observers' outline of EPI boundaries.

constrains the variation of myocardial volume within 5%.

In this chapter, we have presented a novel approach to segmenting the full myocardial volume from cardiac MR and ultrasound images. Motivated by the incompressibility property of myocardium during a cardiac cycle, we coupled the propagation of the ENDO and EPI surfaces according to image-derived information that maximized the piecewise homogeneities of a cardiac image, as well as the incompressibility constraint which

To validate our approach, we computed the MAD, HD, PTP, and PFP between the contours from automatic algorithm and manual segmentation. We observed that when the incompressibility constraint was not applied, the EPI boundaries leaked into other tissues. When incompressibility constraint was applied, however, the MAD, HD, PTP, and PFP were

Assen, H; Danilouchkine, M.; Dirksen, M.; Reiber, J. & Lelieveldt, B. (2008). A 3-D Active

Bouhlel, N.; Sevestre-Ghalila, S. & Rajhi, H. & Hamza, R. (2004). New Markov Random Field

Boukerroui, D.; Noble, A.; Robini, M. & Brady, M. (2001). Enhancement of Contrast Region

Bowman, A. & Kovacs, S. (2003). Assessment and Consequences of the Constant-Volume

Chan, T. & Vese, L. (2001). Active Contours Without Edges. *IEEE Transactions on Image* 

Chen, C.; Pau, L. & Wang, P. (1998). Handbook of Pattern Recognition and Computer Vision

Chen, D.; Chang, R. & Huang, Y. (2000). Breast cancer diagnosis using self-organizing map for sonography. *Ultrasound in Medicine and Biology*, vol.26, no.3 (2000), pp.405-411

*Ultrasound in Medicine & Biology*, vol.27, no.12 (2001), pp.1583-1594

Methods of Clinical Measurement. *Lancet*, vol.1 (1986)

*Circulation Physiology*, vol.285 (2003), pp.H2027-2033

*Processing*, vol.10, no.2 (2001), pp.266-277

(2nd edition). World Scientific

S.Y. Emelianov, pp.363-372

Shape Model Driven by Fuzzy Inference: Application to Cardiac CT and MR. *IEEE Trans. on Information Technology in Biomedicine*, vol.12, no.5 (2008), pp.595-605 Bland, J. & Altman, D. (1986). Statistical Methods for Assessing Agreement Betweeen Two

Model Based on K-Distribution for Textured Ultrasound Image. In *Proceedings of Medical Imaging: Ultrasonic Imaging and Signal Processing, ser. SPIE*, W.F. Walker and

in Suboptimal Ultrasound Images with Application to Echocardiography.

Attribute of the Four-Chamber Heart. *American Journal of Physiology - Heart* 

Fig. 7. Bland-Altman analysis showing the agreement between the MV measurements from manual segmentation and automatic segmentation.

Table 1 shows that the ENDO boundaries were detected with sufficient accuracy even if the incompressibility constraint was not applied. This is because that the ENDO boundaries are relatively clear compared to the EPI boundaries. However, we observed from Table 2 that for EPI boundaries, when the incompressibility constraint was not applied, the automaticmanual MAD was 2.78mm larger than the manual-manual MAD, the automatic-manual HD was 3.82mm larger than the manual-manual HD, and the automatic-manual PTP&PFP were 18%-20% worse than the manual-manual PTP&PFP. This is because of the leakage problem of EPI boundaries without incompressibility constraint. When the incompressibility constraint was applied, however, the automatic-manual MAD decreased by 3.38mm, the automatic-manual HD decreased by 4.11mm, and the automatic-manual PTP&PFP decreased by 16%-19%. We can see from Tables 1 and 2 that the automatic algorithm produced results with comparable accuracy to a manual segmentation. Furthermore, we observed from Tables 1 and 2 that the variability of manual-manual segmentation was smaller for the ENDO surface than for the EPI surface. This is probably because the EPI boundaries are more ambiguous for observers to detect, which was also the reason why we have multiple observers instead of a single one.


Table 1. Comparison of automatic outline to three observers' outline of ENDO boundaries.


Table 2. Comparison of automatic outline to three observers' outline of EPI boundaries.

#### **6. Conclusion**

98 Echocardiography – New Techniques

Fig. 7. Bland-Altman analysis showing the agreement between the MV measurements from

Table 1 shows that the ENDO boundaries were detected with sufficient accuracy even if the incompressibility constraint was not applied. This is because that the ENDO boundaries are relatively clear compared to the EPI boundaries. However, we observed from Table 2 that for EPI boundaries, when the incompressibility constraint was not applied, the automaticmanual MAD was 2.78mm larger than the manual-manual MAD, the automatic-manual HD was 3.82mm larger than the manual-manual HD, and the automatic-manual PTP&PFP were 18%-20% worse than the manual-manual PTP&PFP. This is because of the leakage problem of EPI boundaries without incompressibility constraint. When the incompressibility constraint was applied, however, the automatic-manual MAD decreased by 3.38mm, the automatic-manual HD decreased by 4.11mm, and the automatic-manual PTP&PFP decreased by 16%-19%. We can see from Tables 1 and 2 that the automatic algorithm produced results with comparable accuracy to a manual segmentation. Furthermore, we observed from Tables 1 and 2 that the variability of manual-manual segmentation was smaller for the ENDO surface than for the EPI surface. This is probably because the EPI boundaries are more ambiguous for observers to detect, which was also the reason why we

Table 1. Comparison of automatic outline to three observers' outline of ENDO boundaries.

manual segmentation and automatic segmentation.

have multiple observers instead of a single one.

In this chapter, we have presented a novel approach to segmenting the full myocardial volume from cardiac MR and ultrasound images. Motivated by the incompressibility property of myocardium during a cardiac cycle, we coupled the propagation of the ENDO and EPI surfaces according to image-derived information that maximized the piecewise homogeneities of a cardiac image, as well as the incompressibility constraint which constrains the variation of myocardial volume within 5%.

To validate our approach, we computed the MAD, HD, PTP, and PFP between the contours from automatic algorithm and manual segmentation. We observed that when the incompressibility constraint was not applied, the EPI boundaries leaked into other tissues. When incompressibility constraint was applied, however, the MAD, HD, PTP, and PFP were significantly improved.

#### **7. References**


Automated Segmentation of Real-Time

(1999), pp.784-789

3D Echocardiography Using an Incompressibility Constraint 101

Jakeman, E. (1999). K-distributed Noise. *Journal of Optics A: Pure and Applied Optics*, vol.1

Kass, M.; Witkin, A. & Terzopoulos, D. (1987). Snakes: Active Contour Models. *International* 

King, D.; Coffin, L. & Maurer, M. (2002). Noncompressibility of Myocardium During Systole

Kovesi, P. (1996). Invariant Measures of Image Features from Phase Information. PhD

Kuo, W.; Chang, R.; Chen, D. & Lee, C. (2001). Data mining with decision trees for diagnosis

Lee, W.; Chen, Y. & Hsieh, K. (2003). Ultrasonic liver tissues classification by fractal feature

Lloyd-Jones, D.; et al. (2009). Heart Disease and Stroke Statistics 2009 Update: A Report

Lynch, M.; Ghita, O. & Whelan, P. (2006). Left-Ventricle Myocardium Segmentation Using a

Lynch, M; Ghita, H & Whelan, P (2008). Segmentation of the Left Ventricle of the Heart in 3-

Malladi, R.; Sethian, J. & Vemuri, B. (1995). Shape Modeling with Front Propagation: A

Mojsilovic, A.; Popovic, M.; Markovic, S. & Krstic, M. (1998). Characterization of visually

Nicholas, D.; Nassiri, D.; Garbutt, P.; Hill, C. (1986). Tissue characterization from ultrasound B-scan data. *Ultrasound in Medicine & Biology*, vol.12, no.2 (1986), pp.135-145 O'Donnell, T. & Funka-Lea, G. (1998). 3-D Cardiac Volume Analysis Using Magnetic

Osher, S. & Paragios, N. (2003). Geometric Level Set Methods in Imaging, Vision, and

Paragios, N (2002). A VariationalApproach for the Segmentation of the Left Ventricle in Cardiac Image Analysis. *International Journal of Computer Vision*, vol.50, no.3 (2002), pp.345-362 Philips (2005). http://www.medical.philips.com/main/products/ultrasound/cardiology/.

Pluempitiwiriyawej, C.; Moura, J.; Wu, Y. & Ho, C. (2005). STACS: New Active Contour

Scheme for Cardiac MR Image Segmentation. *IEEE Transactions on Medical Imaging*,

Echocardiography. *Medical Image Analysis*, vol.4 (2000), pp.21-30

With Freehand Three-Dimensional Echocardiography. *Journal of the American* 

of breast tumor in medical ultrasonic images. *Breast Cancer Research and Treatment*,

vector based on M-band wavelet transform. *IEEE Transactions on Medical Imaging*,

From the American Heart Association Statistics Committee and Stroke Statistics

Coupled Level-Set With a Prior Knowledge. *Computerized Medical Imaging and* 

D+t MRI Data Using an Optimized Nonrigid Temporal Model. *IEEE Transactions on* 

Level Set Approach. *IEEE Transactions on Pattern Analysis and Machine Intelligence*,

similar diffuse diseases from B-scan liver images using nonseparable wavelet transform. *IEEE Transactions on Medical Imaging*, vol.17, no.4 (1998), pp.541-549 Mulet-Parada, M. & Noble, A. (2000). 2D+T Acoustic Boundary Detection in

Resonance Imaging. In *Proceedings of the Fourth IEEE Workshop on Applications of* 

*Journal of Computer Vision*, vol.1 (1987), pp.321-331

Thesis, University of Western Australia

*Graphics*, vol.30, no.4 (2006), pp.255-262

*Medical Imaging*, vol.27, no.2 (2008), pp.195-203

*Computer Vision (WACV)* (1998), pp.240-241

Philips Medical Systems Corporate Website

vol.24, no.5 (2005), pp.593-603

vol.66, no.1 (2001), pp.51-57

vol.22, no.3 (2003), pp.382-392

vol.17, no.2 (1995), pp.158-175

Graphics. Springer

*Society of Echocardiogrpahy*, vol.15, no.12 (2002), pp.1503-1506

Subcommittee. *Circulation*, vol.119, (2009), pp.e21-e181


Chen, D.; Chang, R.; Chen, C.; Ho, M.; Kuo, S.; Chen, S.; Hung, S. & Moon, W. (2005).

Chen, T.; Babb, J.; Kellman, P.; Axel, L. & Kim, D. (2008). Semiautomated Segmentation of

Chen, Y.; Huang, F.; Tagare, H.; Rao, M.; Wilson, D. & Geiser, E. (2003). Using Prior Shape and

*IEEE Transactions on Medical Imaging*, vol.27, no.8 (2008), pp.1084-1094 Chen, Y.; Tagare, H.; Thiruvenkadam, S; Huang, F; Wilson, D; Gopinath, K; Briggs, R &

vol.29 (2005), pp.235-245

*Imaging* vol.22, no.7 (2003), pp.902-912

*Imaging*, vol.21, no.9 (2002), pp.1202-1208

(2003)

(2008), pp.287-295

pp.559-565

pp.241-257

(1986), pp.743-748

Classification of breast ultrasound images using fractal feature. *Clinical Imaging*,

Myocardial Contours for Fast Strain Analysis in Cine Displacement-Encoded MRI.

Geiser, E (2002). Using Prior Shapes in Geometric Active Contours in a Variational Framework. *International Journal of Computer Vision*, vol.50, no.3 (2002), pp.315-328 Christodoulou, C.; Pattichis, C.; Pantziaris, M. & Nicolaides, A. (2003). Texture-Based

Classification of Atherosclerotic Carotid Plaques. *IEEE Transactions on Medical* 

Intensity Prior Profile in Medical Image Segmentation. In *Proceedings of International Conference on Computer Vision and Pattern Recognition*, vol.2 (2003), pp.1117-1125 Cohen, B. & Dinstein, I. (2002). New Maximum Likelihood Motion Estimation Schemes for Noisy Ultraound Images. *Pattern Recognition*, vol.48, no.2 (2002), pp.455-463 Cootes, T. & Taylor, C.; Cooper, D. & Graham, J. (1995). Active Shape Models - Their Training and Application. *Computer Vision and Image Understanding*, vol.61, no.1 (1995), pp.38-59 Corsi, C; Saracino, G.; Sarti, A. & Lamberti, C. (2002). Left Ventricular Volume Estimation

for Real-Time Three-Dimensional Echocardiography. *IEEE Transactions on Medical* 

Davignon, F.; Deprez, J. & Basset, O. (2005). A Parametric Imaging Approach for the Segmentation of Ultrasound Data. *Ultrasonics*, vol.43 (2005), pp.789-801 Dutt, V. & Greenleaf, J. (1994). Ultrasound Echo Envelop Analysis Using A Homodyned Kdistribution Signal Model. *Ultrasonic Imaging*, vol.16, no.4 (1994), pp.265-287 Eltoft, T. (2003). Speckle: Modeling and Filtering. In *Norwegian Signal Processing Symoposium*

Glass, L.; Hunter, P. & McCulloch, A. (1990). Theory of Heart: Biomechanics, Biophysics,

Hacihaliloglu, I.; Abugharbieh, R.; Hodgson, A. & Rohling, R. (2008). Bone Segmentation

Hamilton, W. & Rompf, J. (1932). Movement of the Base of the Ventricle and the Relative

Hao, X.; Bruce, C.; Pislaru, C. & Greenleaf, J. (2001). Segmenting High-Frequency Intracardiac

Insana, M.; Wagner, R.; Garra, B.; Brown D. & Shawker, T. (1986). Analysis of Ultrasound

*IEEE Transactions on Medical Imaging*, vol.20, no.12 (2001), pp.1373-1383 Hoffman, E. & Ritman, E. (1985). Invariant Total Heart Volume in the Intact Thorax. *American Journal of Physiology - Heart Circulation Physiology*, vol.249 (1985), pp.H883-890 Hoffman, E. & Ritman, E. (1987). Heart-Lung Interaction: Effect of Regional Lung Air

and Fracture Detection in Ultrasound Using 3D Local Phase Features. In *Proceedings of International Conference on Medical Imaging and Computer Assisted Intervention*

Constancy of the Cardiac Volume. *American Journal of Physiology*, vol.132 (1932),

Ultrasound Images of Myocardium Into Infarcted, Ischemic, and Normal Regions.

Content and Total Heart Volume. *Annual Biomedical Engineering*, vol.15 (1987),

Image Texture via Generalized Rician Statistics. *Optical Engineering*, vol.25, no.6

and Nonlinear Dynamics of Cardiac Function. Springer-Verlag


**6** 

*Japan* 

**Evaluation of Subvalvular Apparatus** 

Yasushige Shingu, Suguru Kubota and Yoshiro Matsui

*Hokkaido University Graduate School of Medicine,* 

*Department of Cardiovascular Surgery* 

**by 2-D Transthoracic Echocardiography** 

Functional mitral regurgitation (MR) is one of the major contributing factors for heart failure and hospitalization in both nonischemic and ischemic dilated cardiomyopathy. The limitation of undersized mitral annuloplasty has been recognized. Recurrent late MR is associated with continued left ventricular (LV) remodeling and enhanced papillary muscle displacement outside posterior ring after mitral annuloplasty. Understanding of the tethering mechanism provided both annular and subvalvular targets for therapy and several

Kron et al. reported a novel approach to reduce tethering length and improve coaptation by approximating displaced posterior papillary muscle toward the annulus. Hvass et al. evolved a papillary muscle sling procedure using an ePTFE tube graft around the muscles and reported improved mitral tethering. Messas et al. proposed a chordal cutting procedure. Langer et al. reported transaortic repositioning of posterior papillary muscle to the midseptal fibrous annulus in a procedure named "RING plus STRING". Our group has reported papillary muscle approximation and its suspension to the mitral artificial ring.

However, there exists no standard method to determine the indication and effects of those various procedures and it makes the comparison between the procedures difficult. Although it would be ideal to determine 3-dimensional anatomy of submitral apparatus, we discuss here the usefulness and importance of 2-D transthoracic echocardiography in this area.

We developed our own submitral procedure, i.e., papillary muscle approximation and its suspension to the mitral ring. According to the method reported by Matsunaga et al., we measure the pre and postoperative values of those parameters relating to the submitral apparatus using 2-D transthoracic echocardiography in patients who need this submitral

Fig. 1 shows our own surgical submitral procedure. Our recent surgical strategy for functional MR is to reconstruct the annulus and subvalvular apparatus of the mitral valve.

**1. Introduction**

**2. Methods** 

procedure for functional MR.

**2.1 Submitral procedure** 

procedures have been reported.

