**Speckle Detection in Echocardiographic Images**

Arezou Akbarian Azar, Hasan Rivaz and Emad M. Boctor

*The Russell H. Morgan Department of Radiology and Radiological Science The Johns Hopkins University School of Medicine, Baltimore USA* 

#### **1. Introduction**

28 Will-be-set-by-IN-TECH

68 Echocardiography – New Techniques

Weickert, J. (1999). Coherence-enhancing diffusion filtering, *Int. J. Comput. Vision*

Weickert, J., Romeny, B. & Viergever, M. (1998). Efficient and reliable schemes for nonlinear

Yu, Y. & Acton, S. (2004). Edge detection in ultrasound imagery using the instantaneous

Yue, Y., Croitoru, M., Bidani, A., Zwischenberger, J. & Clark, J. J. (2006). Nonlinear multiscale

Zhang, F., Yoo, Y. M., Koh, L. M. & Kim, Y. (2007). Nonlinear diffusion in Laplacian pyramid

Zong, X., Laine, A. & Geiser, E. (1998). Speckle reduction and contrast enhancement of echocardiograms via multiscale nonlinear processing, 17(4): 532–540.

wavelet diffusion for speckle suppression and edge enhancement in ultrasound

Yu, Y. & Acton, S. (2002). Speckle reducing anisotropic diffusion, 11(11): 1260–1270.

domain for ultrasonic speckle reduction, 26(2): 200–211.

31(2-3): 111–127.

images, 25(3): 297–311.

diffusion filtering, 7(3): 398–410.

coefficient of variation, 13(12): 1640–1655.

B-scan Echocardiographic images unlikely carry unlikely images by scattering of ultrasound beams backed from structures within the body organ that is being scanned. Two major types of scatterings are diffused and coherent scatterings. Diffuse scatterings are caused when there are a large number of scatterers is with random phase within the resolution cell of the ultrasound beam and causes speckles in the reconstructed image; whereas the e coherent scattering arises when the scatterers in the resolution cell are in phase and causes light or dark spots in the image. Rayleigh distribution is the most common statistical model for the envelope signal and assumes that a large number of scatterers per resolution cell exist. However in some ultrasonic imaging fields such as echocardiography, Rayleigh distribution fits to reflect properties of reflections from blood but fails with complex structures such as myocardial tissue. The K distribution , on the other hand, was initially designed for the envelope signal and have been proposed to model different kinds of tissue in ultrasound envelope imaging. This distribution also has the advantage to model both fully and partially developed speckle. The first-order envelope statistics have been thought to follow a Rayleigh distribution, but recent work has shown that more general models, such as the Nakagami, K, and homodyned K distributions better describe envelope statistics.

With the current digital ultrasound imaging, the radio-frequency (RF) signal has gained more interest as it may contain more information than the envelope echo. When there are a large number of scatterers per range cell it yields Gaussian statistics for the RF signal, but the statistics of the RF signal in the case of partially developed speckle don't follow the Gaussian distribution. Therefore, in this study to model statistical behaviour of the RF data, we used K distribution framework, described in and for such statistics applied them to the RF data. By splitting the ultrasound image to image patches, statistical features for image patches can be extracted using the statistical modelling of the RF signal. These features could overlap for some tissues and the pattern classification approaches should be utilized to classify tissues based on the extracted statistical features. Over past decades, several supervised and unsupervised classification and segmentation algorithms have been proposed to process the medical images. Some of these techniques are listed in. Because of above mentioned problems (overlap between statistical features of tissues) and the fact that in our application (speckle classification), we cannot have enough training material and the data size (=number of image patches) is finite and small, we only focused on the unsupervised clustering techniques in this study.

Speckle Detection in Echocardiographic Images 71

*y py e*

The underlying parameter of the Rayleigh distribution, Σ={σi,j}, associated with each pixel intensity of the RF image, yi,j , is related to the acoustic properties at the corresponding location (i, j), in particular, the so-called echogenicity. Let Z={zi,j} be a N ×M B-mode image corrupted by speckle where each pixel is generated according to the following Log-

where (a, b) are unknown parameters which account for the contrast and brightness

,

*z z i j*

. As pointed out in , (3) defines a double exponential distribution with

*i j*

(, ) <sup>ˆ</sup> <sup>24</sup> *<sup>z</sup> i j i j <sup>a</sup>* 

where σz(i, j) is the SD of the observations inside the window w, centred at the (i, j)th pixel. To estimate the parameter b, we first consider the minimum of the observed pixels zi,j given

> min{ , } log , (min{ } 1) log( 1) *i j s ij a y b at b*

> > *b sZ a t L* ( ) lo g( ( , ) 1)

*L p b sZ tt e t a*

, ,, ,

*b b kpb k s*

*i j i j i j i j*

ˆ ( ) ( ( )| , )

. An estimator of b is found by computing the expected value of bi,j using a

<sup>2</sup> ( | ( ), ) ( 1)

1

*k*

*L*

, 2

known standard deviation (SD) analytical expression, yielding an estimate for a:

,

( 1) ( ) ()

*dy y y pz py <sup>e</sup> dx <sup>a</sup>*

 

*i j*

*i j*

, 2 ,

( )

respectively. Given (1), the distribution of the observed pixels z,

*i j*

with Z = {zi,j}. The distribution of b, derived in , is:

Compression law:

where 1 *z a <sup>b</sup> y e* 

which means that:

where 1

*s b <sup>a</sup> t e* 

numerical approach, such that:

by:

*i j*

2 , 2 , , 2

*i j i j y*

(1)

zi,j = a log(yi,j + 1) + b (2)

2 , 2 ,

*y*

2 <sup>2</sup> 2

(7)

*L*

(8)

*i j i j*

(4)

(3)

(5)

(6)

, , 2

,

*ij ij*

#### **2. Methods**

#### **2.1 Speckles artifact**

One of inherent characteristics of coherent imaging techniques including ultrasound imaging is the presence of speckle-type- noise. Speckle is a random and deterministic pattern in the image formed by the use of the coherent radiation of a medium containing many scatterers. Although the texture of the speckle pattern does not correspond to underscanning structure, the local brightness of the speckle pattern reflects the local echogenicity of the under-scanning scatterers. As can be seen in figure 1, each pixel in an ultrasound image is formed by the back scattered echoes from an approximately ellipsoid called the resolution cell. If each resolution cell in an image patch has many scatterers, the corresponding patch is called fully developed speckle (FDS). Speckle has a negative impact on ultrasound imaging. It has been shown that the detectability of lesion reduced approximately by a factor of eight due to the presence of speckle in the image. To track speckles as well as to estimate probe movement and improve performance of algorithms for adaptive speckle suppression and the elevational separation of B-scans by speckle decorrelation, we needed to model both fully speckle (blood pool) and partially developed speckle (tissue area). To this end, in the next subsections we investigated the ability of various statistical modelling of the RF signal and different unsupervised clustering techniques.

Fig. 1. Ultrasound beam and resolution cell.

#### **2.2 Statistical features for speckle classification**

#### **2.2.1 Rayleigh distribution**

Given the assumption of fully developed speckle, the envelope RF image patch, Y = {yi,j}, is modeled by Rayleigh statistics, where the probability density function (pdf) is given by:

One of inherent characteristics of coherent imaging techniques including ultrasound imaging is the presence of speckle-type- noise. Speckle is a random and deterministic pattern in the image formed by the use of the coherent radiation of a medium containing many scatterers. Although the texture of the speckle pattern does not correspond to underscanning structure, the local brightness of the speckle pattern reflects the local echogenicity of the under-scanning scatterers. As can be seen in figure 1, each pixel in an ultrasound image is formed by the back scattered echoes from an approximately ellipsoid called the resolution cell. If each resolution cell in an image patch has many scatterers, the corresponding patch is called fully developed speckle (FDS). Speckle has a negative impact on ultrasound imaging. It has been shown that the detectability of lesion reduced approximately by a factor of eight due to the presence of speckle in the image. To track speckles as well as to estimate probe movement and improve performance of algorithms for adaptive speckle suppression and the elevational separation of B-scans by speckle decorrelation, we needed to model both fully speckle (blood pool) and partially developed speckle (tissue area). To this end, in the next subsections we investigated the ability of various statistical modelling of the RF signal and

**2. Methods** 

**2.1 Speckles artifact** 

different unsupervised clustering techniques.

Fig. 1. Ultrasound beam and resolution cell.

**2.2.1 Rayleigh distribution** 

**2.2 Statistical features for speckle classification** 

Given the assumption of fully developed speckle, the envelope RF image patch, Y = {yi,j}, is modeled by Rayleigh statistics, where the probability density function (pdf) is given by:

$$p\left(y\_{i\_{i,j}}\right) = \frac{y\_{i\_{i,j}}}{\sigma \frac{2}{i\_{i,j}}} e^{-\frac{\left\|\frac{y\_{i,j}^2}{2\sigma\sigma\_{i,j}^2}\right\|}{2\left\|\sigma\right\|^2}}\tag{1}$$

The underlying parameter of the Rayleigh distribution, Σ={σi,j}, associated with each pixel intensity of the RF image, yi,j , is related to the acoustic properties at the corresponding location (i, j), in particular, the so-called echogenicity. Let Z={zi,j} be a N ×M B-mode image corrupted by speckle where each pixel is generated according to the following Log-Compression law:

$$\text{zi}\_{\text{j}} = \text{a } \log(\text{yi}\_{\text{j}} + 1) + \text{b} \tag{2}$$

where (a, b) are unknown parameters which account for the contrast and brightness respectively. Given (1), the distribution of the observed pixels z,

$$p(z\_{i,j}) = \left\{ \left| \frac{dy}{dx} \right| p(y) \right\}\_{z = z\_{i,j}} = \frac{y\_{i,j} (y\_{i,j} + 1)}{a \sigma\_{i,j}^2} e^{-\frac{y\_{i,j}^2}{2 \sigma\_{i,j}^2}} \tag{3}$$

where 1 *z a <sup>b</sup> y e* . As pointed out in , (3) defines a double exponential distribution with known standard deviation (SD) analytical expression, yielding an estimate for a:

$$
\hat{a}\_{i,j} = \sqrt{24} \frac{\sigma\_z(i,j)}{\pi} \tag{4}
$$

where σz(i, j) is the SD of the observations inside the window w, centred at the (i, j)th pixel.

To estimate the parameter b, we first consider the minimum of the observed pixels zi,j given by:

$$\begin{aligned} s &= \min\{i, j\} = a \log(\min\{y\_{i,j}\} + 1) + b \\ &= a \log(t + 1) + b \end{aligned} \tag{5}$$

which means that:

$$b = \text{s}(Z) - a \log(\text{t}(\sigma, L) + 1) \tag{6}$$

with Z = {zi,j}. The distribution of b, derived in , is:

$$p(b \mid s(Z), \sigma) = \frac{L}{a\sigma^2} t(t+1) e^{-\frac{L}{2\sigma^2}} t^2 \tag{7}$$

where 1 *s b <sup>a</sup> t e* . An estimator of b is found by computing the expected value of bi,j using a numerical approach, such that:

$$\hat{b}\_{i,j} = \sum\_{k=1}^{L} b\_{i,j}(k) p(b\_{i,j}(k) \mid s\_{\prime} \beta\_{i,j}) \tag{8}$$

Speckle Detection in Echocardiographic Images 73

The corresponding distribution is called KRF distribution in the following. This pdf may thus provide the basis for segmentation of ultrasonic images in the case of partially developed speckle. However, KRF distribution has the following disadvantages based on the literature: Numerical simulation show that estimation bias grows rapidly as parameter n increases,

Equation (2) implies repeated evaluation of a Bessel function, increasing the computational

The coherent signal to the diffuse signal energy ratio (r) for a patch in an ultrasound image and the effective number of scatterers per resolution cell (µ) could be used as statistical features to identify speckles and characterize tissues. Speckle detection is a useful input for adaptive speckle suppression algorithms and for use in decorrelation algorithms to estimate

To find r and µ, we need to model envelope signal behaviour. Having found r and µ we can then label as speckle patches with high µ and low k. Based on we calculate these statistics on

*mean A*

*Std A A*

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

 

( ) ( ) ( ) ( )

*A A*

*A A A A*

*A A*

Where Std means standard deviation, <. > =mean. Based on the results in , values of v more than one is suggested to perform well. As explained in the previous subsections, each statistical models for the RF signal has some advantages and some disadvantages. This suggests that we must consider both statistical models (Rayleigh and KRF) to better characterize statistical behaviour of the RF signal. Therefore for each image patch A, we

After extracting features for each image patch A, we can use the classification scheme shown

propose to compute statistical features in (10) and Maximum Likelihood (ML)

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

 

 

 

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

  3

(12)

*i j* , in (8)

4

yielding unreliable estimates in blood regions (i.e. υ>> 1);

**2.2.3 Statistical features for speckle classification** 

the elevational distance between neighbouring B-scans.

*R*

*Skewness*

*Kurtosis*

estimation of the data A following the Rayleigh distribution.

in figure 2 to classify each image patch to FDS and non-FDS.

Fig. 2. Pipeline for Speckle detection.

arbitrary powers v of the image patch A.

cost of the algorithm.

where bi,j(k) = k s/(L − 1) and k = 0, 1, ...,L − 1 are L uniformly distributed values in the interval [0, s], since b ≥ 0 and from (6), b ≤ s. In (8),

$$\mathcal{B}\_{i,j} = \sqrt{\frac{1}{2\text{nm}} \Omega\_{k,l} q\_{k,l}^2} \tag{9}$$

is Maximum Likelihood (ML) estimation of σi,j from the pixels inside the window [*w,qk,l*].

The estimators of a and b, considered constant across the image, are obtained by averaging the estimates , ˆ*<sup>i</sup> <sup>j</sup> a* and , ˆ *<sup>i</sup> <sup>j</sup> <sup>b</sup>* , such that: , <sup>1</sup> ˆ ˆ *N M i j i j a a NM* and , <sup>1</sup> ˆ ˆ *N M i j i j b b NM* .

These parameters ( ˆ *a b* ˆ, ) are used to retrieve the envelope RF image according to:

$$y\_{i,j} = e \frac{z\_{i,j} - b}{a} - 1 \tag{10}$$

In some ultrasonic imaging fields such as echocardiography Rayleigh distribution fits to reflect properties of t model reflections from blood but fails with complex structures such as myocardial tissue. Another statistical model named K distribution have been proposed to model different kinds of tissue in ultrasound envelope imaging.

#### **2.2.2 K distribution**

The K distribution has been developed for the envelope signal. The interest of such distribution in ultrasonic images relies on its ability to model both fully speckle (blood pool) and partially developed speckle (tissue area) situations. In this section, we briefly describe the K distribution and the corresponding pdf for the RF signal. The backscattered ultrasonic signal results from the individual energy contributions of each scatterer embedded in the resolution cell. This situation can mathematically be described as a random walk in the complex plane. From this random flight model, the analytic signal can be expressed as a random process depending on the number of scatterers present inside the resolution cell, their relative position and contribution. Thus, a joint density function of the envelope and phase can be obtained by expressing both statistical properties of the phase and amplitude of each scatterer. This results in a K distribution when the scatterers phase is assumed to be uniformly distributed and when their amplitude is modelled as a K distribution itself.

The RF signal corresponds to the real part of the analytic signal. The power density function (pdf) of the RF signal thus corresponds to the marginal distribution obtained by integrating the pdf corresponding to the analytic signal with respect to its imaginary part, which yields the following expression (see for details):

$$pdf^{RF}(\mathbf{x}) = \frac{b}{\sqrt{\pi}\,\Gamma(\nu)} (\frac{b|\mathbf{x}|}{2})^{\nu - 0.5} K\_{\nu - 0.5}(b|\mathbf{x}|) \tag{11}$$

where is the Gamma function and *Kv*0.5 is the modified Bessel function of the second kind of order *v* 0.5 .

This expression is completely specified by its two parameters υ & b, such that n controls the shape and b the scale of the pdf.

where bi,j(k) = k s/(L − 1) and k = 0, 1, ...,L − 1 are L uniformly distributed values in the

, , , 1 <sup>2</sup> *<sup>i</sup> <sup>j</sup> kl kl <sup>q</sup> nm*

> *i j i j*

, , <sup>1</sup> *i j*

In some ultrasonic imaging fields such as echocardiography Rayleigh distribution fits to reflect properties of t model reflections from blood but fails with complex structures such as myocardial tissue. Another statistical model named K distribution have been proposed to

The K distribution has been developed for the envelope signal. The interest of such distribution in ultrasonic images relies on its ability to model both fully speckle (blood pool) and partially developed speckle (tissue area) situations. In this section, we briefly describe the K distribution and the corresponding pdf for the RF signal. The backscattered ultrasonic signal results from the individual energy contributions of each scatterer embedded in the resolution cell. This situation can mathematically be described as a random walk in the complex plane. From this random flight model, the analytic signal can be expressed as a random process depending on the number of scatterers present inside the resolution cell, their relative position and contribution. Thus, a joint density function of the envelope and phase can be obtained by expressing both statistical properties of the phase and amplitude of each scatterer. This results in a K distribution when the scatterers phase is assumed to be uniformly distributed and when their amplitude is modelled as a K distribution itself.

The RF signal corresponds to the real part of the analytic signal. The power density function (pdf) of the RF signal thus corresponds to the marginal distribution obtained by integrating the pdf corresponding to the analytic signal with respect to its imaginary part, which yields

> *b b x pdf x K b x*

where is the Gamma function and *Kv*0.5 is the modified Bessel function of the second

This expression is completely specified by its two parameters υ & b, such that n controls the

  0.5

*v*

(11)

0.5 () ( ) ( ) .() <sup>2</sup>

*z b*

*a*

*a a NM* and ,

is Maximum Likelihood (ML) estimation of σi,j from the pixels inside the window [*w,qk,l*]. The estimators of a and b, considered constant across the image, are obtained by averaging

<sup>1</sup> ˆ ˆ *N M*

These parameters ( ˆ *a b* ˆ, ) are used to retrieve the envelope RF image according to:

*i j*

*y e*

*<sup>i</sup> <sup>j</sup> <sup>b</sup>* , such that: ,

model different kinds of tissue in ultrasound envelope imaging.

2

(9)

*i j i j b b NM* .

(10)

<sup>1</sup> ˆ ˆ *N M*

interval [0, s], since b ≥ 0 and from (6), b ≤ s. In (8),

ˆ

the following expression (see for details):

kind of order *v* 0.5 .

shape and b the scale of the pdf.

*RF*

the estimates , ˆ*<sup>i</sup> <sup>j</sup> a* and ,

**2.2.2 K distribution** 

The corresponding distribution is called KRF distribution in the following. This pdf may thus provide the basis for segmentation of ultrasonic images in the case of partially developed speckle. However, KRF distribution has the following disadvantages based on the literature:

Numerical simulation show that estimation bias grows rapidly as parameter n increases, yielding unreliable estimates in blood regions (i.e. υ>> 1);

Equation (2) implies repeated evaluation of a Bessel function, increasing the computational cost of the algorithm.

#### **2.2.3 Statistical features for speckle classification**

The coherent signal to the diffuse signal energy ratio (r) for a patch in an ultrasound image and the effective number of scatterers per resolution cell (µ) could be used as statistical features to identify speckles and characterize tissues. Speckle detection is a useful input for adaptive speckle suppression algorithms and for use in decorrelation algorithms to estimate the elevational distance between neighbouring B-scans.

To find r and µ, we need to model envelope signal behaviour. Having found r and µ we can then label as speckle patches with high µ and low k. Based on we calculate these statistics on arbitrary powers v of the image patch A.

$$\begin{aligned} R &= \frac{mean}{Std} = \frac{\left\{A^{\nu}\right\}}{\sqrt{\left\{A^{2\nu}\right\} - \left\{A^{\nu}\right\}^2}}\\ Skewness &= \frac{\left\{A^{\nu} - \left\{A^{\nu}\right\}\right\}^3}{\left(\left\{A^{2\nu}\right\} - \left\{A^{\nu}\right\}^2\right)^{3/2}}\\ Kurtsis &= \frac{\left\{A^{\nu} - \left\{A^{\nu}\right\}\right\}^4}{\left(\left\{A^{2\nu}\right\} - \left\{A^{\nu}\right\}^2\right)^{4/2}} \end{aligned} \tag{12}$$

Where Std means standard deviation, <. > =mean. Based on the results in , values of v more than one is suggested to perform well. As explained in the previous subsections, each statistical models for the RF signal has some advantages and some disadvantages. This suggests that we must consider both statistical models (Rayleigh and KRF) to better characterize statistical behaviour of the RF signal. Therefore for each image patch A, we propose to compute statistical features in (10) and Maximum Likelihood (ML) *i j* , in (8) estimation of the data A following the Rayleigh distribution.

After extracting features for each image patch A, we can use the classification scheme shown in figure 2 to classify each image patch to FDS and non-FDS.

Fig. 2. Pipeline for Speckle detection.

Speckle Detection in Echocardiographic Images 75

1, 1

2 1

This solution also satisfies the constraints on fuzzy partitions. Note that equation for cluster centers gives Vj as the weighted mean of the data items that belong to a cluster, where the weights are the membership degrees. That is why the algorithm is called "c-means". The FCM algorithm computes with the standard Euclidean distance norm, which induces hyperspherical clusters. Hence it can only detect clusters with the same shape and

It is extended version of the standard fuzzy c-means algorithm by employing an adaptive distance norm, in order to detect clusters of different geometrical shapes in one data set.

 The matrices Ai are used as optimization variables in the c-means functional, thus allowing each cluster to adapt the distance norm to the local topological structure of the data. Cost function for Gustafson-Kessel clustering algorithm (GK) algorithm is defined as following:

1 1

*i k J D* 

in (17). To optimize the cluster's shape, it is proven that Ai should be as following:

1

*F*

 

*<sup>k</sup> <sup>i</sup> <sup>N</sup> <sup>m</sup>*

*C N <sup>m</sup>*

where the N x C matrix U=[µik] represents the fuzzy partitions. D means distance as defined

1/ 1 [ det( )] . *<sup>n</sup> A FF ii i i* 

where *ρi* is a fixed value for each cluster and det() means determinant of a matrix and F

1

( )( ) *<sup>N</sup> m T ik k i k i*

*X VX V*

*ik k*

,1 ,

*i N*

<sup>1</sup> ,1 ,1

(15)

(16)

(20)

*j c iN*

<sup>2</sup> ( ) ( ), 1 , 1 *<sup>T</sup> D x v Ax v i c k N ikA k i i k i* (17)

(18)

(19)

2

. *<sup>i</sup>*

*ik ikA*

1

*C ji j* 

( )

*X V X V*

*ji <sup>C</sup> i j <sup>m</sup> k i k <sup>N</sup> <sup>m</sup> ji i*

> *ji i*

*X V i N*

orientation, because the common choice of norm matrix is the identical matrix.

If ||Xi-Vj||2>0 and m>1, then (U,V) may minimize J only if

**2.3.3 Gustafson-Kessel clustering algorithm** 

Each cluster has its own norm-inducing matrix Ai:

means fuzzy covariance matrix defined by:

1

1

*<sup>i</sup> <sup>j</sup> <sup>N</sup> <sup>m</sup>*

1

#### **2.3 Unsupervised clustering for speckle classification**

Data clustering means partitioning data to fuzzy or crisp (hard) subsets. Hard clustering in a data set X means partitioning the data into a specified number of subsets of X with such a condition that an object either does or does not belong to a cluster. The number of subsets (clusters) is denoted by K. The hard partitioning is the simplest approach for data clustering, though its results are not always reliable and has numerical problems as well. However, fuzzy clustering allows objects to belong to multiple clusters in the same time, with different degrees of membership. In many real applications fuzzy clustering is more realistic than hard clustering, as objects on the boundaries between several classes are not forced to fully belong to one of the classes. In this study we used both hard (K-means and K-mediod) and fuzzy partitioning (Fuzzy C-Means, Gustafson-Kessel and Gath-Geva techniques) for speckle detection in a competitive manner.

#### **2.3.1 K-means and K-medoid**

The hard partitioning methods are simple and popular, though its results are not always reliable and these algorithms have numerical problems as well. From an Nxn dimensional data set K-means and K-medoid allocates algorithms allocates each data point to one of K clusters (in our case 2 for FDS and non-FDS) to minimize the within-cluster sum of squares (distance norm):

$$\sum\_{i=1}^{K} \sum\_{i \in A\_{\downarrow}} \left\| \mathbf{X}\_{i} - \mathbf{V}\_{j} \right\| \tag{13}$$

where Aj is a set of objects (data points) in the j-th cluster. In k-means, Vj is the mean for that points over cluster j and is called the cluster prototypes or centers. While in K-medoid clustering the cluster centers are the nearest objects to the mean of data in one cluster. Kmediod is useful when there is no continuity in the data space.

#### **2.3.2 Fuzzy C-means clustering**

The Fuzzy C-means clustering algorithm is based on the minimization of an objective function called C-means functional. It is defined by Dunn as:

$$J = \sum\_{j=1}^{C} \sum\_{i=1}^{N} \left(\mu\_{ji}\right)^{m} \left\|{X\_i} - V\_j\right\|^2 \tag{14}$$

Where µ is membership degree, m is order, and Vj is center of the cluster j, which all of them have to be determined. The minimization of the c-means functional J (cost function) represents a nonlinear optimization problem that can be solved by using a variety of available methods, ranging from grouped coordinate minimization, over simulated annealing to genetic algorithms. The most popular method, however, is a simple Picard iteration through the first-order conditions for stationary points of J, known as the fuzzy cmeans (FCM) algorithm.

The stationary points of the objective function J can be found by adjoining the following constraint for the fuzzy partitions U=[µji]:

Data clustering means partitioning data to fuzzy or crisp (hard) subsets. Hard clustering in a data set X means partitioning the data into a specified number of subsets of X with such a condition that an object either does or does not belong to a cluster. The number of subsets (clusters) is denoted by K. The hard partitioning is the simplest approach for data clustering, though its results are not always reliable and has numerical problems as well. However, fuzzy clustering allows objects to belong to multiple clusters in the same time, with different degrees of membership. In many real applications fuzzy clustering is more realistic than hard clustering, as objects on the boundaries between several classes are not forced to fully belong to one of the classes. In this study we used both hard (K-means and K-mediod) and fuzzy partitioning (Fuzzy C-Means, Gustafson-Kessel and Gath-Geva techniques) for

The hard partitioning methods are simple and popular, though its results are not always reliable and these algorithms have numerical problems as well. From an Nxn dimensional data set K-means and K-medoid allocates algorithms allocates each data point to one of K clusters (in our case 2 for FDS and non-FDS) to minimize the within-cluster sum of squares

*i j*

(13)

*X V*

where Aj is a set of objects (data points) in the j-th cluster. In k-means, Vj is the mean for that points over cluster j and is called the cluster prototypes or centers. While in K-medoid clustering the cluster centers are the nearest objects to the mean of data in one cluster. K-

The Fuzzy C-means clustering algorithm is based on the minimization of an objective

Where µ is membership degree, m is order, and Vj is center of the cluster j, which all of them have to be determined. The minimization of the c-means functional J (cost function) represents a nonlinear optimization problem that can be solved by using a variety of available methods, ranging from grouped coordinate minimization, over simulated annealing to genetic algorithms. The most popular method, however, is a simple Picard iteration through the first-order conditions for stationary points of J, known as the fuzzy c-

The stationary points of the objective function J can be found by adjoining the following

*ji i j*

*V*

2

(14)

1 *<sup>j</sup>*

*i iA*

1 1 ( ) *C N m*

*j i J X* 

*K*

mediod is useful when there is no continuity in the data space.

function called C-means functional. It is defined by Dunn as:

**2.3 Unsupervised clustering for speckle classification** 

speckle detection in a competitive manner.

**2.3.1 K-means and K-medoid** 

**2.3.2 Fuzzy C-means clustering** 

means (FCM) algorithm.

constraint for the fuzzy partitions U=[µji]:

(distance norm):

$$\sum\_{j=1}^{\mathbb{C}} \mu\_{ji} = 1, \ 1 \le i \le N \tag{15}$$

If ||Xi-Vj||2>0 and m>1, then (U,V) may minimize J only if

$$\begin{aligned} \mu\_{ji} &= \frac{1}{\sum\_{k=1}^{C} \left\| \frac{\mathbf{X}\_i - V\_j}{\|\mathbf{X}\_i - V\_k\|} \right\|\_{m=1}^2}, \quad 1 \le j \le c, \quad 1 \le i \le N\\ \sum\_{k=1}^{N} \mu\_{ji}^m \mathbf{X}\_i\\ V\_j &= \frac{i-1}{\sum\_{i=1}^{N} \mu\_{ji}^m}, \quad 1 \le i \le N \end{aligned} \tag{16}$$

This solution also satisfies the constraints on fuzzy partitions. Note that equation for cluster centers gives Vj as the weighted mean of the data items that belong to a cluster, where the weights are the membership degrees. That is why the algorithm is called "c-means". The FCM algorithm computes with the standard Euclidean distance norm, which induces hyperspherical clusters. Hence it can only detect clusters with the same shape and orientation, because the common choice of norm matrix is the identical matrix.

#### **2.3.3 Gustafson-Kessel clustering algorithm**

It is extended version of the standard fuzzy c-means algorithm by employing an adaptive distance norm, in order to detect clusters of different geometrical shapes in one data set. Each cluster has its own norm-inducing matrix Ai:

$$\left|D\_{\mathbb{kA}}\right|^2 = \left(\mathbb{x}\_{\mathbb{k}} - \mathbb{z}\_{\mathbb{i}}\right)^T A \left(\mathbb{x}\_{\mathbb{k}} - \mathbb{z}\_{\mathbb{i}}\right), \quad 1 \le i \le c, \quad 1 \le k \le N \tag{17}$$

 The matrices Ai are used as optimization variables in the c-means functional, thus allowing each cluster to adapt the distance norm to the local topological structure of the data. Cost function for Gustafson-Kessel clustering algorithm (GK) algorithm is defined as following:

$$J = \sum\_{i=1}^{C} \sum\_{k=1}^{N} \mu\_{ik}^{\;\;\;m} \, D\_{ik}^{\;\;2} \, \_{ik} \tag{18}$$

where the N x C matrix U=[µik] represents the fuzzy partitions. D means distance as defined in (17). To optimize the cluster's shape, it is proven that Ai should be as following:

$$A\_i = \left[ \rho\_i \det(F\_i) \right]^{1/n}.F\_i^{-1} \tag{19}$$

where *ρi* is a fixed value for each cluster and det() means determinant of a matrix and F means fuzzy covariance matrix defined by:

$$F\_i = \frac{\sum\_{k=1}^{N} \mu\_{ik}^{\text{ue}} (X\_k - V\_i)(X\_k - V\_i)^T}{\sum\_{k=1}^{N} \mu\_{ik}^{\text{ue}}} \tag{20}$$

Speckle Detection in Echocardiographic Images 77

(LV) cardiac images the convolutional approach that ignores the geometry of the transducer were used. As opposed to the Field II simulator ignoring geometry allows to speed up simulations but at the cost of realism of the simulated images. For the simulation it is assumed that the imaging system has a linear, space-invariant point spread function (PSF) and the transducer is linear. We also applied our speckle classification scheme on other

To evaluate our proposed speckle classification scheme and calculate performance of each unsupervised classification technique, as the ground truth, we used B-mode images simulated by ultrasound simulation programs with 100,000 scatterers and 128 RF lines (see Figure 3). After resizing of reconstructed images to the size of 1200x800 pixels, they were segmented to 12x8 image patches, where each image patch had size of 100x100 pixels. We then calculated following statistical features for image patch: R, Skewness and Kurtosis features for the k-distribution and Maximum Likelihood (ML) for the Rayleigh distribution. After calculating statistical features for each image patch, we will have 4 dimensional features for each image patch and we can classify them to FDS and non-FDS using unsupervised clustering techniques. For this purpose in this study we applied five pattern classification techniques: K-means, K-medoid, Fuzzy C-means, Gustafson-Kessel fuzzy classifier and Gath-Geva fuzzy classifier. Figures 4,5,6 and 7 show performance of the classification methods for speckle detection, respectively for cyst, fetus, LV and short axes heart phantoms. Figure 8 show the specke detection performance for a real ultrasound

Fig. 3. Simulation Examples: A) short-axis, B) left ventricle, C) cyst and D) fetus phantoms and patches. Images were segmented to 12x8 image patches, where each image patch had

phantoms: fetus and cyst phantoms generated Field II simulator.

image of a beating heart using 3D Ultrasound.

size of 100x100 pixels.

**3. Results** 

Indeed Ai is a generalized squared Mahalanobis distance norm between point Xk and the cluster mean Vi, where the covariance is weighted by the membership degrees in U=[µik]. Mahalanobis distance is a distance based on correlations between variables by which different patterns can be identified and analyzed. It is a useful way of determining similarity of an unknown sample set to a known one. It differs from Euclidean distance in that it takes into account the correlations of the data set and is scale-invariant, i.e. not dependent on the scale of measurements.

#### **2.3.4 Gath-Geva fuzzy classifier**

The Gath-Geva fuzzy clustering algorithm employs a distance norm based on the fuzzy maximum likelihood estimates (FMLE), proposed by Bezdek and Dunn :

$$D\_{lk}(X\_k, V\_l) = \frac{\sqrt{\det(F\_{wi})}}{a\_i} \exp\{0.5(X\_k - V\_l)^T F\_{wi}^{-1} (X\_k - V\_l)\} \tag{21}$$

Note that, contrary to the GK algorithm, this distance norm involves an exponential term and thus decreases faster. Fwi denotes the fuzzy covariance matrix of the i-the cluster, given by:

$$F\_{wl} = \frac{\sum\_{k=1}^{N} \mu\_{ik}^{\text{ue}} (X\_k - V\_i)(X\_k - V\_i)^T}{\sum\_{k=1}^{N} \mu\_{ik}^{\text{ue}}}, \quad 1 \le i \le c \tag{22}$$

If w = 1, Fwi is identical to the original FMLE algorithm, but if w = 2, so that the partition becomes more fuzzy to compensate the exponential term of the distance norm. The difference between the matrix Fi in GK algoritm and the Fwi in (22) is that the latter does not involve the weighting exponent m, instead of this it consists of w = 1. This is because the two weighted covariance matrices arise as generalizations of the classical covariance from two different concepts. The αi is the prior probability of selecting cluster is given by:

$$
\mu\_i = \frac{1}{N} \sum\_{k=1}^{N} \mu\_{ik} \tag{23}
$$

The membership degrees µik are interpreted as the posterior probabilities of selecting the i-th cluster given the data point xk. Gath and Geva reported that the fuzzy maximum likelihood estimates clustering algorithm is able to detect clusters of varying shapes, sizes and densities. The cluster covariance matrix is used in conjunction with an "exponential" distance, and the clusters are not constrained in volume. However, this algorithm is less robust in the sense that it needs a good initialization, since due to the exponential distance norm, it converges to a near local optimum.

#### **2.4 Echocardiographic images simulation**

In order to quantitatively compare statistical features and classification methods, as ground truth, we used two different types of ultrasound simulation programs. To simulate shortaxis echocardiographic images Field II simulator were used and to simulate left-ventricle (LV) cardiac images the convolutional approach that ignores the geometry of the transducer were used. As opposed to the Field II simulator ignoring geometry allows to speed up simulations but at the cost of realism of the simulated images. For the simulation it is assumed that the imaging system has a linear, space-invariant point spread function (PSF) and the transducer is linear. We also applied our speckle classification scheme on other phantoms: fetus and cyst phantoms generated Field II simulator.

#### **3. Results**

76 Echocardiography – New Techniques

Indeed Ai is a generalized squared Mahalanobis distance norm between point Xk and the cluster mean Vi, where the covariance is weighted by the membership degrees in U=[µik]. Mahalanobis distance is a distance based on correlations between variables by which different patterns can be identified and analyzed. It is a useful way of determining similarity of an unknown sample set to a known one. It differs from Euclidean distance in that it takes into account the correlations of the data set and is scale-invariant, i.e. not dependent on the

The Gath-Geva fuzzy clustering algorithm employs a distance norm based on the fuzzy

det( ) <sup>1</sup> ( ,) exp(0.5( ) ( )) *wi <sup>T</sup> ik k i k i wi k i*

Note that, contrary to the GK algorithm, this distance norm involves an exponential term and thus decreases faster. Fwi denotes the fuzzy covariance matrix of the i-the cluster, given

( )( )

 

*X VX V F i c*

> *ik k*

If w = 1, Fwi is identical to the original FMLE algorithm, but if w = 2, so that the partition becomes more fuzzy to compensate the exponential term of the distance norm. The difference between the matrix Fi in GK algoritm and the Fwi in (22) is that the latter does not involve the weighting exponent m, instead of this it consists of w = 1. This is because the two weighted covariance matrices arise as generalizations of the classical covariance from

> 1 1 *<sup>N</sup> <sup>i</sup> ik <sup>N</sup> <sup>k</sup>*

The membership degrees µik are interpreted as the posterior probabilities of selecting the i-th cluster given the data point xk. Gath and Geva reported that the fuzzy maximum likelihood estimates clustering algorithm is able to detect clusters of varying shapes, sizes and densities. The cluster covariance matrix is used in conjunction with an "exponential" distance, and the clusters are not constrained in volume. However, this algorithm is less robust in the sense that it needs a good initialization, since due to the exponential distance

In order to quantitatively compare statistical features and classification methods, as ground truth, we used two different types of ultrasound simulation programs. To simulate shortaxis echocardiographic images Field II simulator were used and to simulate left-ventricle

*<sup>N</sup> w T ik k i k i*

(21)

(23)

(22)

, 1

*D XV X VFX V*

1

two different concepts. The αi is the prior probability of selecting cluster is given by:

maximum likelihood estimates (FMLE), proposed by Bezdek and Dunn :

*i F*

1

norm, it converges to a near local optimum.

**2.4 Echocardiographic images simulation** 

*<sup>k</sup> wi <sup>N</sup> <sup>w</sup>*

scale of measurements.

by:

**2.3.4 Gath-Geva fuzzy classifier** 

To evaluate our proposed speckle classification scheme and calculate performance of each unsupervised classification technique, as the ground truth, we used B-mode images simulated by ultrasound simulation programs with 100,000 scatterers and 128 RF lines (see Figure 3). After resizing of reconstructed images to the size of 1200x800 pixels, they were segmented to 12x8 image patches, where each image patch had size of 100x100 pixels. We then calculated following statistical features for image patch: R, Skewness and Kurtosis features for the k-distribution and Maximum Likelihood (ML) for the Rayleigh distribution. After calculating statistical features for each image patch, we will have 4 dimensional features for each image patch and we can classify them to FDS and non-FDS using unsupervised clustering techniques. For this purpose in this study we applied five pattern classification techniques: K-means, K-medoid, Fuzzy C-means, Gustafson-Kessel fuzzy classifier and Gath-Geva fuzzy classifier. Figures 4,5,6 and 7 show performance of the classification methods for speckle detection, respectively for cyst, fetus, LV and short axes heart phantoms. Figure 8 show the specke detection performance for a real ultrasound image of a beating heart using 3D Ultrasound.

Fig. 3. Simulation Examples: A) short-axis, B) left ventricle, C) cyst and D) fetus phantoms and patches. Images were segmented to 12x8 image patches, where each image patch had size of 100x100 pixels.

Speckle Detection in Echocardiographic Images 79

Fig. 6. Simulated ultrasound image of left ventricle (A) and speckle detection results for five different unsupervised classifiers. Patches classified as fully developed speckles (FDS) are shown as black. Total number of patches (100x100 pixels) for the phantom image was 96.

Fig. 7. Simulated ultrasound image of the heart in short axis (end diastolic) (A) and speckle detection results for five different unsupervised classifiers. Patches classified as fully developed speckles (FDS) are shown as black. Total number of patches (100x100 pixels) for the phantom image was 96. Orders for statistical features respectively was 1,1 and 1.

Orders for statistical features respectively was 1,1 and 1.

Fig. 5. Simulated ultrasound image of a fetus in 12th week (A) and speckle detection results for five different unsupervised classifiers. Patches classified as fully developed speckles (FDS) are shown as black. As can be seen, FCM,K-means and K-mediod performed the same. GK and GG fuzzy classifiers were able to decrease false positives and improve accuracy of the speckle detection. Total number of patches (100x100 pixels) for the phantom image was 96. Orders for statistical features respectively was 1,1 and 0.5.

Fig. 4. Simulated ultrasound image of a cyst phantom (A) and speckle detection results for five different unsupervised classifiers. Total number of 100x100 patches for the phantom image was 24. Patches classified as fully developed speckles (FDS) are shown as black. All methods except GK-fuzzy classifier performed the same. Orders for statistical features

Fig. 5. Simulated ultrasound image of a fetus in 12th week (A) and speckle detection results for five different unsupervised classifiers. Patches classified as fully developed speckles (FDS) are shown as black. As can be seen, FCM,K-means and K-mediod performed the same. GK and GG fuzzy classifiers were able to decrease false positives and improve accuracy of the speckle detection. Total number of patches (100x100 pixels) for the phantom

image was 96. Orders for statistical features respectively was 1,1 and 0.5.

respectively were 1,1 and 0.5.

Fig. 6. Simulated ultrasound image of left ventricle (A) and speckle detection results for five different unsupervised classifiers. Patches classified as fully developed speckles (FDS) are shown as black. Total number of patches (100x100 pixels) for the phantom image was 96. Orders for statistical features respectively was 1,1 and 1.

Fig. 7. Simulated ultrasound image of the heart in short axis (end diastolic) (A) and speckle detection results for five different unsupervised classifiers. Patches classified as fully developed speckles (FDS) are shown as black. Total number of patches (100x100 pixels) for the phantom image was 96. Orders for statistical features respectively was 1,1 and 1.

Speckle Detection in Echocardiographic Images 81

Bamber, J.C., Dickinson, R.J., 1980. Ultrasonic B-scanning: a computer simulation. Phys Med

Bamber, J.C., Phelps, J.V., 1991. Real-time implementation of coherent speckle suppression

Bankman, I.N., Nizialek, T., Simon, I., Gatewood, O.B., Weinberg, I.N., Brody, W.R., 1997.

Bernard, O., D'Hooge, J., Friboulet, D., 2006. Statistics of the radio-frequency signal based on

Bernard, O., Touil, B., Gelas, A., Prost, R., Friboulet, D., 2007. Segmentation of Myocardial

Bezdek, J.C., 1975. Optimal Fuzzy Partitions: A Heuristic for Estimating the Parameters in a Mixture of Normal Distributions. IEEE Transactions on Computers. 24, 835-838. Bezdek, J.C., Hall, L.O., Clarke, L.P., 1993. Review of MR image segmentation techniques

Dutt, V., Greenleaf, J.F., 1994. Ultrasound echo envelope analysis using a homodyned K distribution signal model, Vol. 16, Dynamedia, Silver Spring, MD, ETATS-UNIS. Gath, I., 1989. Unsupervised Optimal Fuzzy Clustering. IEEE Transactions on Pattern

Gullo, F., Ponti, G., Tagarelli, A., 2008. Clustering Uncertain Data Via K-Medoids. In:

Kobayashi, T., Kida, Y., Tanaka, T., Kageyama, N., Kobayashi, H., Amemiya, Y., 1986.

with low Curie temperature. Journal of Neuro-Oncology. 4, 175-181. Matsumoto, M., Yoshimura, N., Honda, Y., Hiraoka, M., Ohura, K., 1994. Ferromagnetic

Archive for Clinical and Experimental Ophthalmology. 232, 176-181. Molthen, R.C., Shankar, P.M., Reid, J.M., 1995. Characterization of ultrasonic B-scans using non-rayleigh statistics. Ultrasound in medicine & biology. 21, 161-170. Molthen, R.C., Shankar, P.M., Reid, J.M., Forsberg, F., Halpern, E.J., Piccoli, C.W., Goldberg,

breast and liver tissue. Ultrasound in medicine & biology. 24, 93-100. Ordsmith, R.J., 1967. Abramowitz,M - Handbook Mathematical Functions with Formulas

Ossant, F., Patat, F., Lebertre, M., Teriierooiterai, M.L., Pourcelot, L., 1998. Effective density

Graphs and Mathematical Tables. Computer Journal. 10, 45-&.

moments. Ultrason Imaging. 20, 243-59.

Management. Vol., ed.^eds. Springer-Verlag, Naples, Italy, pp. 229-242. Jensen, J.A., 2004. Simulation of advanced ultrasound systems using Field II. In: Biomedical

Proceedings of the 2nd international conference on Scalable Uncertainty

Imaging: Nano to Macro, 2004. IEEE International Symposium on. Vol., ed.^eds.,

Magnetic induction hyperthermia for brain tumor using ferromagnetic implant

hyperthermia in rabbit eyes using a new glass-ceramic thermoseed. Graefe's

B.B., 1998. Comparisons of the Rayleigh and K-distribution models using in vivo

estimators based on the K distribution: interest of low and fractional order

Frequency Control, IEEE Transactions on. 53, 1689-1694.

using pattern recognition. Med Phys. 20, 1033-48.

Analysis and Machine Intelligence. 11, 773-780.

Segmentation algorithms for detecting microcalcifications in mammograms. IEEE

K distribution with application to echocardiography. Ultrasonics, Ferroelectrics and

Regions in Echocardiography Using the Statistics of the Radio-Frequency Signal. In: Functional Imaging and Modeling of the Heart. Lecture Notes in Computer Science, Vol. 4466, F. Sachse, G. Seemann, ed.^eds. Springer Berlin / Heidelberg,

Biol. 25, 463-79.

pp. 433-442.

pp. 636-639 Vol. 1.

in B-scan images. Ultrasonics. 29, 218-24.

Trans Inf Technol Biomed. 1, 141-9.

Fig. 8. Ultrasound image of the right ventricle (A) and speckle detection results for five different unsupervised classifiers. Patches classified as fully developed speckles (FDS) are shown as black. Total number of patches (30x30 pixels) for the ultrasound image was 96. Orders for statistical features respectively was 1,1 and 1.

#### **4. Summary**

In this study we reviewed the present statistical models to predict behaviour of the RF signal for different tissue types. In some ultrasonic imaging fields such as echocardiography Rayleigh distribution fits to reflect properties of reflections from blood but fails with complex structures such as myocardial tissue. However, another statistical model named K distribution have been proposed to model different kinds of tissue in ultrasound envelope imaging. To consider both statistical properties of the echocardiographic images we applied both Rayleigh and Kdistributions and for each image patch A, we computed statistical features for K-distribution and Maximum Likelihood (ML) estimation of the image patch following the Rayleigh distribution. After extracting features for each image patch, we applied the unsupervised clustering techniques to classify each image patch to FDS and non-FDS. Based on our observation, we found when we use all statistical features (R-S-K- ML) together classifiers was able to separate classes in data space better than with other combinations out of R-S-K- ML features. Based on Table 1 ranking of the fuzzy classification methods performed better. To improve specificity and sensitivity of the proposed machine-learning speckle detection scheme, feature whitening and mapping techniques such as Sammon mapping can applied to increase distance between features before applying unsupervised classification techniques.

#### **5. References**

Abonyi, J., Babuska, R., Szeifert, F., 2002. Modified Gath-Geva fuzzy clustering for identification of Takagi-Sugeno fuzzy models. Systems, Man, and Cybernetics, Part B: Cybernetics, IEEE Transactions on. 32, 612-621.

Fig. 8. Ultrasound image of the right ventricle (A) and speckle detection results for five different unsupervised classifiers. Patches classified as fully developed speckles (FDS) are shown as black. Total number of patches (30x30 pixels) for the ultrasound image was 96.

In this study we reviewed the present statistical models to predict behaviour of the RF signal for different tissue types. In some ultrasonic imaging fields such as echocardiography Rayleigh distribution fits to reflect properties of reflections from blood but fails with complex structures such as myocardial tissue. However, another statistical model named K distribution have been proposed to model different kinds of tissue in ultrasound envelope imaging. To consider both statistical properties of the echocardiographic images we applied both Rayleigh and Kdistributions and for each image patch A, we computed statistical features for K-distribution and Maximum Likelihood (ML) estimation of the image patch following the Rayleigh distribution. After extracting features for each image patch, we applied the unsupervised clustering techniques to classify each image patch to FDS and non-FDS. Based on our observation, we found when we use all statistical features (R-S-K- ML) together classifiers was able to separate classes in data space better than with other combinations out of R-S-K- ML features. Based on Table 1 ranking of the fuzzy classification methods performed better. To improve specificity and sensitivity of the proposed machine-learning speckle detection scheme, feature whitening and mapping techniques such as Sammon mapping can applied to increase distance between features before applying unsupervised classification techniques.

Abonyi, J., Babuska, R., Szeifert, F., 2002. Modified Gath-Geva fuzzy clustering for

B: Cybernetics, IEEE Transactions on. 32, 612-621.

identification of Takagi-Sugeno fuzzy models. Systems, Man, and Cybernetics, Part

Orders for statistical features respectively was 1,1 and 1.

**4. Summary** 

**5. References** 


**Part 2** 

**New Techniques in** 

**Trans-Thoracic Echocardiography** 

