J. Biomedical Science and Engineering, 2010, 3, 1078-1084 JBiSE
doi:10.4236/jbise.2010.311140 Published Online November 2010 (http://www.SciRP.org/journal/jbise/).
Published Online November 2010 in SciRes. http://www.scirp.org/journal/jbise
Medical ultrasound image segmentation by modified local
histogram range image method
Ali Kermani, Ahmad Ayatollahi, Ahmad Mirzaei, Mohammad Barekatain
School of Electrical Engineering, Iran University of Science and Technology, Tehran, Iran.
Received 8 September 2010; revised 8 October 2010; accepted 13 October 2010.
ABSTRACT
Fast and satisfied medical ultrasound segmentation is
known to be difficult due to speckle noises and other
artificial effects. Since speckle noise is formed from
random signals which are emitted by an ultrasound
system, we can’t encounter the same way as other
image noises. Lack of information in ultrasound im-
ages is another problem. Thus, segmentation results
may not be accurate enough by means of customary
image segmentation methods. Those methods that
can specify undesirable effects and segment them by
eliminating artificial effects, should be chosen. It
seems to be a complicated work with high computa-
tional load. The current study presents a different
approach to ultrasound image segmentation that re-
lies mainly on local evaluation, named as local histo-
gram range image method which is modified by
means of discrete wavelet transform. Thus, a signifi-
cant decrease in computational load is then achieved.
The results show that it is possible for tissues to be
segmented correctly.
Keywords: Segmentation; Local Histogram; Ultrasound
Image; Morphological Image Processing; Discrete Wave-
let Transform
1. INTRODUCTION
Fast and reliable ultrasound image segmentation is a
complicated process with particular difficulties. Because
of the presence of speckle in these images and their poor
contrast, using common methods for segmentation is not
possible. General approaches to segmentation can be
grouped into three categories: 1) pixel-based, 2) continu-
ity, and 3) edge based methods. Many segmentation
methods that are presented, include clustering, watershed,
active contours, and those methods that are based on
statistical models [1-4]. Pixel-based methods are the
easiest to implement, since they apply to one element
while, they are particularly susceptible to noise. Conti-
nuity-based and edge-based methods encounter the seg-
mentation problem from opposing sides: continuity-
based methods search for similarities while edge-based
methods search for differences [5].
On the other hand, almost these segmentation methods
consider the entire image as a unit, which contains a
large area of gray value, and uses spatial or frequency
for segmenting; however, it can’t recognize speckle
noises accurately. Local evaluation methods can be a
solution for ultrasound images, since speckle pixels will
be recognized precisely. In [6], local estimation of
Rayleigh parameter was proposed to identify different
tissues. However, its complexity and computational load
are high. In [7], researchers defined “local histogram
range image” (LHRI) based on histogram distribution
which reduced statistic complexity. They applied a clas-
sification method to recognize edge or border of organ
while, speckle noises remained unchanged. LHRI was
proposed for segmentation in [8]; although, its computa-
tional load was too high.
In this study, the combined edge-based and continuity-
based method are presented by modified LHRI. Discrete
wavelet transform (DWT) is applied to reduce dimen-
sion of input images of LHRI algorithm; whereas, the
image energy compresses in few number of wavelet's
coefficients. Hence, we expect to have proper segmented
images and to decrease computational load by means of
reducing input image size.
We will introduce the LHRI method for an ultrasound
image in Part 2. Section 3 describes DWT. Combination
of LHRI method with DWT is introduced in Section 4.
In Section 5 the simulation results are presented.
2. LOCAL HISTOGRAM RANGE IMAGE
FOR SEGMENTATION
2.1. Definition
The recognition of ultrasound signal characteristics is
the primary requirement. Speckle has a random nature as

A. Kermani et al. / J. Biomedical Science and Engineering 3 (2010) 1078-1084
1079
it is formed from signals which are produced by scatter-
ers in the medium under observation. This analysis of
speckles is one of the most important studies. In addition,
today's advanced scanners apply a log-compression to
the envelope data to enhance the lower signal values.
Considering that log-compression is a nonlinear opera-
tion, statistics of envelope signal are different before and
after this compression.
Because of different scatterer’s distribution and den-
sity in different tissues, researchers have studied several
probability density models of ultrasound envelope sig-
nals. We divide image areas into two general groups to
evaluate an ultrasound image: 1-fully developed speckle
(a large number of randomly located scatterers with
small scatterer spacing were compared to the wavelength
of ultrasound), 2-areas that have low scatterers. Several
distributions are presented for each mentioned area that
we only introduce two usual distributions. Statistics
demonstrate Rayleigh model in fully developed speckle
areas [9]. In the second area that the number of scatterers
is low and its spatial locations are dependent variables,
Nakagami distribution is the most flexible model among
the others [10]. Of course, these distributions are due to
previous logarithmic compression step. However, refer-
ence [7] has shown that each distribution can be
achieved after applying logarithmic compression. For
the first area (fully developed speckle), Raleigh distribu-
tion convert into double exponential distribution that its
mean and variance are obtained as follows:
(
)
( )
11
12
ln2
ln,
22
LogRayleigh
nn
Meannn
γ
σ
−=+−+ (1)
( )
22
2
1
24
LogRayleigh
n
VAREyy π
−
=−=
 (2)
where
γ
=0.5772 is Euler constant,
1
n
is image dy-
namic range, and
2
n
is gain setting (ratio of minimum
to maximum output signal). It can be concluded that the
mean depends on three factors: initial Raleigh variance,
1
n
and
2
n
, while the variance only depends on dy-
namic range
1
n
.
Also, log-Nakagami [9] was derived for latter area
that its variance is as follows:
()( )
2
21
,
4
LogNakagami n
VAREyym
ξ
−
=−=
 (3)
where
ξ
is defined as follow:
( )( )
22
11
...,
1
mmm ξ=++
+ (4)
Log-Nakagami distribution tends to Rician in the case
of m = 1.5, which leads to the variance below:
( )
22
2
1
42.23
LogRice
n
VAREyy π
−
=−=
 . (5)
Since both distributions depend on one parameter
1
()
n
, their comparison shows that the variance of log-
Rayleigh exceeds the variance of log-Nakagami. Hence,
variance of log-Rayleigh is more than log Nakagami’s
variance.
After introducing the distribution, now it's the time to
find out the relationship between borders and mentioned
distributions. The borders of ultrasound images are to-
tally divided into two groups: 1) the edge between two
fully developed speckle areas, and 2) the Edge between
speckle area and specular scattered areas. The experi-
ments on ultrasound images show that although histo-
grams of fully developed speckle areas have similar
structure; however, their mean values could be different
for these two different organs. Also, equal variances and
different means are concluded from double exponential
equations. In addition, local histogram range in borders
between two fully developed speckle areas is more than
the range of local histogram in fully developed areas. We
can see from (5) that the variance of Rician distribution
after log compression will be a constant, and its value is
smaller than the one in Rayleigh distribution. Therefore,
valid range of the histogram will not be very large, while
its mean is larger than double exponential mean. Local
histogram between two specular scatterer and speckle
has the largest range among other areas.
We can recognize edges, and segment image by means
of proper function via prior analysis and relation be-
tween distributions and valid local histogram range,
without high computation load of statistic calculations.
Therefore, we use local histogram range image defini-
tion [7]: at first we should choose proper size for moving
window that suits our approach. The experiment on sev-
eral ultrasound images shows that the best window size
for segmentation is
44
×
to
1515
×
. This window
moves around original image and an
NM
×
matrix is
obtained as below:
,
, if u different gray values
exist in the moving window
0, if no different gray values
exist in the moving window,
ij
u
k




=




(6)
where
,
ij
k
is the pixel value of LHRI at the position a(i,
j).
2.2. Classifier Function
Here, a function is described that separates different im-
age areas after computation of the LHRI matrix. In [7],

A. Kermani et al. / J. Biomedical Science and Engineering 3 (2010) 1078-1084
1080
authors have used a linear adaptive function to enhance
edges; we modified it for segmentation. Our purpose is
clustering pixels in two levels, foreground and back-
ground. The mean of whole image is a general margin
for separating two classes, but it is a static parameter;
therefore, it doesn ’t work properly for all ultrasound
images. We need dynamic parameter to classify pixels
truly. For this reason, we define classifier function below
with dynamic parameter
,
ij
k
obtained from LHRI ma-
trix:
(
)
,,,,ijijijij
fgkgg
α=+− (7)
where
,
ij
f
is the enhanced value,
,
ij
g
is the old one,
g
is the mean of whole image, and
α
is the coeffi-
cient value. Distance of each pixel’s gray value from
mean is computed at first in this equation. This distance
can be positive or negative. Output pixel would be
brighter if this value is positive, and darker if it is nega-
tive. The proportion of this change depends on
,
ij
k
. It
can be explained as local adaptive classifying parameter,
whereas more
,
ij
k
causes more margin parameter near
the whole image mean.
Figure 1 shows this function for 8 bit gray level in-
puts and several k (0 to 9) and
127
g
=
. It can be seen
that the defined function classifies and ranks gray level
values of pixels in terms of LHRI pixel value and whole
image mean.
Another problem with this approach is the
α
value.
Implementation of this algorithm indicates that proper
results would be attained by coefficient values between
0.1 to 1 which depend on gray value range of original
image.
2.3. Morphological Image Processing
Since LHRI algorithm leads to the existence of small
holes in the obtained matrix, it can’t classify tissues
properly. For this reason, we apply morphological image
processing to achieve perfect and smooth tissues. Con-
sequently, dilation and erosion are selected. At first, we
apply dilation to eliminate holes and fill the contour gaps
[12]. The dilation is defined as:
(
)
{
}
z
ABzBA
ϕ
⊕=∩≠
)
(8)
where... is a median filtered image, B is a circle or a
square structuring element and denotes reflection of B.
Erosion is applied to eliminate the effects caused by di-
lation in the tissue size. It is defined as:
(
)
{
}
z
ABzBA
Θ=⊆
)
(9)
The key problem with this explanation is the size of
structuring element, because small structuring doesn ’t
fill holes correctly. Furthermore, large size of them
makes artificial edges. Experiments show that
33
×
element could fill holes properly.
3. DISCRETE WAVELET TRANSFORM
The concept of “wavelet” was first introduced in 1984
Figure 1. Classifier function for k = {0, 2, 3, ..., 9}.

A. Kermani et al. / J. Biomedical Science and Engineering 3 (2010) 1078-1084
1081
by A. Grossmann and J. Morlet [13]. By wavelet trans-
form (WT), each signal included in 2
()
LR
is presented
as weighted summation of wavelet basis functions. In
WT, the basis functions are obtained from a single func-
tion named “wavelet basis functions” by time translation
(b) and dilation (a). Eq.10 presents wavelet function:
()
,1 ba
xb
x
a
a
ψψ
−

=×


(10)
If a, and b, in continuous wavelet transform, are used
in binary form, the DWT is yielded [14]:
2,2
mm
abn
==×
(11)
Thus, DWT is defined as:
( )()
()
2
2,222,
m
mmm
f
Wnxnfxdx
ψ
+∞
−−
−∞
=−
∫ (12)
Inverse discrete wavelet transform (IDWT) is ob-
tained by following equation:
()
( )
()
2,2
mm
fmn
mn
fxWnx
ψ
+∞+∞
=−∞=−∞
=∑∑ (13)
Each wavelet basis function has nonzero value in spe-
cial frequency interval (band). Therefore, it concludes
signal information in the special frequency interval. Also,
basis functions are selected in a way that they are or-
thogonal. This means that there are no overlaps between
the special frequency intervals of various basis functions
[14]. This special frequency interval of wavelet basis is
introduced as follows:
(
)
{
}
2.
m
m
WspanxnmZ
ψ−
≡−∈
(14)
Another set of functions named scaling function, is
introduced in a way that its special frequency band
comes from union the special frequency interval of
wavelet basis. If its special frequency band of scaling
function is defined as
m
V
, we have [14]:
( )
12
-
, V, V
m
mi
i
VWLR
ϕ
−
∞+∞
=−∞
===
U (15)
Scaling functions that are orthogonal basis for
m
V
,
are obtained thorough time translation and dilation of
special function as follows:
()
( )
2
,
22 nZ
mm
nm xxnϕϕ
−
−
=−∈
(16)
Therefore, for survey of the signal in all frequency
bands, DWT should be calculated for
[,)
ma
∈+∞
along with signal decomposition with
(1)
an
ϕ− [15]. In
such a way signal decomposition is as follows:
()
( )( )
0
00
2
222
j
j j
jkjk
kjjk
fxdxkcxk
ϕψ
+∞
=
=−+−
∑∑∑
(17)
where
jk
c
and
0
jk
d
are approximations and details co-
efficients of WT respectively, and they are computed by:
() ()
00
jkjk
dfxxdx
ϕ
+∞
−∞
=∫ (18)
() ()
jkjk
cfxxdx
ψ
+∞
−∞
=∫ (19)
These coefficients for discrete signal are calculated by
[16]:
(
)
1,
2
jnjk
k
chknc+ =−
∑
(20)
(
)
1,
2
jnjk
k
dgknc
+
=−
∑
(21)
where h(k) low pass filter coefficient and g(k) is high
pass filter coefficient.
We need to generalize mentioned subjects into two
dimensions under to apply DWT on images. Four fun-
damental functions must be used for two dimensional
DWT, which elicit low frequency information and high
frequency information in three directions ‘x’, ’y’ and
diagonal. These functions are assumed to be separable
for simplicity. Therefore, we have from one dimensional
wavelet and scaling functions [17]:
)()(),(
)()(),(
)()(),(
)()(),(
yxyx
yxyx
yxyx
yxyx
d
h
V
ψψψ
ϕψψ
ψϕψ
ϕϕϕ
=
=
=
=
(22)
where these functions compute four types of coefficients:
approximation, vertical, horizontal and diagonal coeffi-
cients. Two dimensional DWT computations is con-
verted to double calculation of one dimensional DWT by
using separable mentioned functions.
Input image is divided into four images which consists
of approximation image, vertical, horizontal and diago-
nal edges for each decomposition level in DWT. Ap-
proximation part will be a shrank image that contains a
large amount of primary image energy. So we can apply
some image processing algorithms (such as segmenta-
tion) on approximation coefficients to decrease compu-
tation load.
4. COMBINATION OF LHRI METHOD
WITH DWT
As mentioned in previous section, for decreasing execu-
tion time of algorithm, we can calculate approximation
coefficients from DWT, which have a large amount of
image energy, and then use them in LHRI algorithm.
It should be noted that using DWT and its inverse in
proper portions of algorithm is extremely important. In

A. Kermani et al. / J. Biomedical Science and Engineering 3 (2010) 1078-1084
1082
general we can do it in two forms. We apply them and
choose the best one by comparing them.
In the first method, we determine DWT of original
image, and then approximated wavelet coefficients are
used as LHRI input. We calculate IDWT before applying
local histogram range image (LHRI) on primary image.
Subsequently, the obtained matrix is defined as LHRI
matrix in next algorithm processes.
The difference between the second method and former
is in using LHRI matrix coefficients on approximation
coefficients of prior image, which is obtained with DWT.
We finally use inverse transform for conclusion, to make
similar output size with input image. We applied two
introduced methods on almost fifty ultrasound images,
and came to the conclusion that the former method seg-
ments images with more details. Figure 2 shows block
diagram of proposed combined method. Step by Step
implementation of this algorithm on an example image is
demonstrated in Figure 3.
More specifically, if we apply DWT for more than one,
the execution time of algorithm decreases considerably,
and algorithm accuracy is reduced. Therefore, there is a
trade-off between them. The experimental results show
that two levels of decomposition have the best efficiency.
Input image size become quarter for adding each de-
composition level in DWT, whereas wavelet computa-
tion time is low. Therefore, total computation time is
reduced in similar ratio approximately. As a result, the
execute time decreases in
116
ratio for applying two
decomposition levels.
5. SIMULATION RESULTS
We applied our algorithm on two images of carotid and
breast lesion to evaluate the results. At first, we get
DWT from input images. We should choose optimum
filter which used in this transform. As [18] proved, we
choose biorthogonal filter. Figure 4(a) shows original
ultrasound breast lesion. Figure 4(b) illustrates seg-
mented image by LHRI alone. As shown in Figure 4(c),
Figure 2. Block diagram of combined method.
Figure 3. Step by step implementation of main algorithm.

A. Kermani et al. / J. Biomedical Science and Engineering 3 (2010) 1078-1084
1083
(a) (b) (c)
Figure 4. Segmentation of the ultrasound breast lesion image. (a) original image; (b) LHRI
method result; (c) the method used in this study.
(a) (b) (c)
Figure 5. Segmentation of the ultrasound carotid image. (a) original image; (b)LHRI method
result; (c) the method used in this study.
(a) (b) (c)
(d) (e) (f)
Figure 6. Breast lesion ultrasound image segmentation, (a, c, and e) original images, (b, d, and f)
segmented results respectively.
segmentation result by LHRI & wavelet method is simi-
lar to the former, whereas execution time is reduced con-
siderably. We can prove these results for Figure 5, too.
Table 1 shows the amount of simulation time for two
algorithm implementations with 1.86 GHz CPU by Mat-
lab software for Figure 4. These times are illustrated for
Figure 5 in Table 2. Figure 6 illustrates three breast
lesions that are segmented by proposed method.

A. Kermani et al. / J. Biomedical Science and Engineering 3 (2010) 1078-1084
1084
Table 1. Execution time of simulation for Breast lesion.
LHRI DWT IWT Total
Algorithm
LHRI
Algorithm 1.054 sec - - 1.054 sec
LHRI +
Wavelet 0.127 sec 0.156 sec
0.016 sec 0.297 sec
Decrement
percent 87.95% - - 71.82%
Table 2. Execution time of simulation for carotid.
LHRI DWT IWT Total
Algorithm
LHRI
Algorithm 6.227 sec
- - 6.227 sec
LHRI +
Wavelet 0.549 sec
0.235 sec 0.078 sec 0.862 sec
Decrement
percent 91.18% - - 86.15%
We can see that the modified method by Wavelet seg-
mented images with details and reduced implementation
time about 90% of primary algorithm.
6. CONCLUSIONS
LHRI method segments envelope ultrasound images
effectively. However, also its computation load is prac-
tically high. We show that DWT can decrease execution
time without any change in segmentation result. Conse-
quently, the required memory for algorithm implementa-
tion, is reduced. To obtain more accurate results while
decreasing the computation load, we will use the modi-
fied wavelet and more advanced morphological proc
REFERENCES
[1] Chan, T. and Vese, L. (2001) Active contours without
edges. IEEE Transaction on Image Processing, 10, 266-
272.
[2] Grau, V., Mewes, A.U.J., Alcaniz, M., Kikinis, R. and
Warfield, S.K. (2004) Improved watershed transform for
medical image segmentation using prior information
IEEE Transactions on Medical Imaging, 23(4), 447-458.
[3] Wu, J. and Chung, A.C.S. (2007) “A segmentation model
using compound Markov random fields based on a boun-
dary model,” IEEE Transaction on Image Processing, 16,
241-252.
[4] Makrogiannis, S., Economou, G., Fotopoulos, S. and
Bourbakis, N.G. (2005) Segmentation of colour images
using multiscale clustering and graph theoretic region
synthesis. IEEE Transaction on System Man and Cyber-
netics, 35(2), 224-238.
[5] Semmlow, J.L. (2004) Biosignal and biomedical image
processing. CRC Press, Boca Raton.
[6] Aja-Fern´andez, S., Martin-Fern´andez, M. and Alber-
ola-L´opez, C. (2007) Tissue identification in ultrasound
images using Rayleigh local parameter estimation. Pro-
ceedings of Bioinformatics and Bioengineering, Boston
[7] Wang, B. and Liu, D.C. (2008) A novel edge enhance-
ment method for ultrasound imaging. The 2nd Interna-
tional Conference on Bioinformatics and Biomedical En-
gineering, 3, 645-649.
[8] kermani, A., Ayatollahi, A. and Talebi, M. (2010) Seg-
mentation medical ultrasound image based on Local his-
togram range image. The 3rd International Conference
on BioMedical Engineering and Informatics, in Press.
[9] Shankar, P.M. (2000) A general statistical model for ul-
trasonic backscattering from tissues. Ultrasonics, Ferro-
electrics and Frequency Control, 47(3), 727-736.
[10] Wagner, R.F., Smith, S.W., Sandrik, J.M. and Lopez, H.
(1983) Statistics of speckle in ultrasound B-scans. IEEE
Transactions on Sonics and Ultrasonics, 30(3), 156-163.
[11] Ghofrani, S., Jahed-Motlagh, M.R. and Ayatollahi, A.
(2001) An adaptive speckle suppression filter based on
Nakagami distribution. Proceeding of IEEE International
Conference on Trends in Communications, Bratislava, 1,
84-87.
[12] Gonzalez, R.C. and Woods, R.E. (2002) Digital image
processing. 2nd Edition, Publishing House of Electronics
Industry, Beijing, 523-527.
[13] Topiwala, P.N., Ed., (2007) Wavelet image and video
compression. The Springer International Series in Engi-
neering and Computer Science, Cambridge.
[14] Vetterli, M. H. (1992) Wavelets and filter banks: Theory
and design. IEEE Transactions on Signal Processing, 40
(9), 2207-2232.
[15] Zhenzhu, Y., Yong, Y., Zhenxi, C. and Kongyang, Z.
(2009) Study of the de-noising method based on wavelet
and fractal. Second International Workshop on Computer
Science and Engineering, WCSE’09, 1, 15-19.
[16] Guo, F., Derong, T.L. and Can, Z. (2005) A new com-
pression algorithm for medical images using wavelet
transform. IEEE Conferences on Networking, Sensing
and Control, Tucson.
[17] Gupta, P.K. and Kanhirodan, R. (2006) Wavelet trans-
form based error concealment approach for image de-
noising. 1st IEEE Conference on Industrial Electronics
and Applications, Singapore.
[18] Xiaojuan, L., Guangshu, H. and Shangkai, G. (1999) De-
sign and implementation of a novel compression method
in a tele-ultrasound system. IEEE Transactions on In-
formation Technology in Biomedicine, 3(3), pp. 205-213.