# Trabecular Bone Image Segmentation Using Wavelet and Marker-Controlled Watershed Transformation

Authors: Wafa Abid Fourati; Mohamed Salim Bouhlel
Journal: Chinese Journal of Engineering
 (2014)

---

## Abstract

This paper presents a new strategy for the segmentation of trabecular bone image. This kind of image is acquired with microcomputed tomography (micro-CT) to assess bone microarchitecture based  chiefly on bone mineral density (BMD) measurements to improve fracture risk prediction. Disease osteoporosis can be predicted from features of CT image where a bone region may consist of several disjoint pieces. It relies on a multiresolution representation of the image by the wavelet transform to compute the multiscale morphological gradient. The coefficients of detail found at the different scales are used to determine the markers and homogeneous regions that are extracted with the watershed algorithm. The method reduces the tendency of the watershed algorithm to oversegment and results in closed homogeneous regions. The performance of the proposed segmentation scheme is presented via experimental results obtained with a broad series of images.


---

## Body

## 1. Introduction

Osteoporosis is a metabolic bone disease [1]. The disease is defined by low bone mass and microarchitectural deterioration of bone tissues leading to enhanced bone fragility. Clinically, osteoporosis is associated with an increased risk of fractures of the vertebral bone. The exact clinical estimation of bone strength and fracture risk is important for the treatment of bone diseases such as osteoporosis. It designates a major health problem in industrialized countries. This is a silent disease without any symptoms. After the age of 50, the number of osteoporotic people clearly raises [2]. A hip fracture is costly for a system of health care. Vertebral fractures are associated with a higher risk of mortality and also having another fracture in the first months following a fracture. Segmentation of bone structures from computed tomography (CT) images (shown in Figure 1) is to extract osseous tissue on the surface of bones (cortical bone, the trabecular network) from the interior region of bones (cancellous bone, space medullar).Trabecular bone image.
(a)
(b)The trabecular network detection is a very important task to help diagnose diseases such as trabecular bone loss in osteoporotic patients or in animal models of osteoporosis [1–3].The low contrast images of the trabecular bone acquired with the sky-scan X-ray microcomputed tomography make their analysis a challenging task. Additionally, throughout acquisition, it is required to establish a high number of parameters that result in the presence of noise, nonhomogeneous illumination, and fuzzy contours. Such segmentations can be very challenging to obtain since osseous tissue does not always produce readily discernible features from soft tissue regions in CT images. Also, the trabeculae of the cancellous bone induce a certain texture pattern in the interior of the bone regions leading to large intensity variations in regions corresponding to cancellous bone.Manual segmentation of bone regions in CT images is tedious, time consuming, and subject to observer variability; as a result, there is a strong demand for automating the segmentation procedure.Watershed Transform [4] is a powerful morphological tool to segment images into regions of interest. This transform is able to characterize extremely composite objects. However, oversegmentation and sensitivity to noise are major problems of most watershed algorithms. The accomplishment of the Watershed Transform depends essentially on a good estimation of image gradients and the existence of obvious markers [5] for each of the objects of interest.The automatic determination of markers is still a difficult goal to achieve. Multiscale image analysis provides a suggestion of where visually sensible objects in an image are located and also information about their relative size or importance. Wavelet transform is a good multiresolution method for edge enhancement, noise removal, and compression [6]. Morphology is also an effective method for image processing.In this context, some methods have been proposed in the literature to analyze and simplify in different ways the image features. Although considerable research has been conducted on the segmentation problem, existing solutions to bone segmentation in CT images still remain insufficient for systematic deployment and process in a clinical setting due to drawbacks such as difficulties in algorithm initialization and parameter settings, as well as a lack of robustness to image inhomogeneity or noise. Gauch [7] has suggested an image segmentation and analysis technique employing multiscale gradient watershed hierarchies. The methodology is based on imposing a scale-based hierarchy on the watersheds associated with minima under Gaussian blurring. With the same objective, Mukhopadhyay and Chanda [8] used a hierarchy multiscale morphological operator by reconstruction for image segmentation. The result is satisfying to some extent. However, those samples of images contain a circular feature and depend on the size of the structuring element.In another study, Kim et al. [9] proposed a multiresolution based watershed image segmentation scheme, which consists of segmenting the lowest resolution image by watershed, then using a wavelet coefficient-based region merged to reduce oversegmentation, and finally projecting the segmented low resolution image onto a full resolution image. Though, the region merging process is complicated and the calculation cost is very high.Also, Jung [10] used a combined image denoising/enhancement technique to obtain an enhanced version of image gradient from the approximation image of the wavelet transform at the coarsest resolution based on an adaptive threshold. Then, the segmented image is projected up to higher resolutions. The method needs postprocessing to remove noise and some prior knowledge such as the rule of the postprocessing. O’Callaghan and Bull [11] proposed a combined morphological-spectral unsupervised image segmentation method. This used the subbands of the dual-tree complex wavelet transform in order to process both textured and nontextured objects in a meaningful fashion. Then, the watershed transform provides the segmented image.Meyer [12] introduced the leveling approach, which consists of using morphological filters to remove small details of the image. The result is satisfying to some extent which has a clear edge, a small amount of noise, and a regular shape. Yet, the micro-CT trabecular bone images present a low contrast edge and intense amount of noise.In this paper, we propose to segment the trabecular bone images using a watershed method with a multiresolution analysis to compute a gradient recognition and controlled markers. In the first step we suggest creating a reconstructed morphological gradient image by using the multiresolution representation of the image which reduces the oversegmentation of the watershed algorithm. Then, we proposed using the wavelets transform to locate the regions of interest and the morphological operations to obtain and preprocess the skeleton in order to obtain markers that allow robust segmentations in images.The paper is organized as follows: in Section2, the watershed and the wavelet transform are described; Section 3 focuses on our approach for trabecular bone segmentation and its implementation. Before concluding, the results are discussed in Section 4.


## 2. The Watershed Transform and the Wavelet Transform

### 2.1. The Watershed Transform

The watershed transformation is a morphological tool dedicated to image segmentation [4, 5, 13]. According to a flooding paradigm, this transformation considers each gray level (or color or texture) as topographic data as shown in Figure 2.Figure 2
3D profile of a gradient image.The principle is to first search out all the local minima where to perforate holes. Then, immerse the model into a lake at a uniform speed as the time. During the flooding stage, the water level slowly rises, covering an increasing amount of the image. To prevent the fusion of two catchment basins corresponding to different regional minima, a dam is erected between them at the points where they congregate each other. The watershed line corresponds to a series of this constructed dam and local extremes. For digital images, watershed constitutes closed contour of one-pix width [6]. Different approaches are employed to use the watershed principle for image segmentation. One of the most common watershed algorithms was introduced by Meyer [6]. This flooding process is performed on the gradient image. The implementation of this algorithm uses the hierarchical queue [5].This transformation is characterized by its fast computing and high accuracy in locating the weak edges of adjacent regions. But normally this will lead to an oversegmentation of the image. Therefore, many researchers [14–18] have proposed various methods for improvement, such as a preprocessing step before watershed transformation.


### 2.2. The Wavelet Transform

The wavelet transform [19] is a mathematical tool that can represent signal (image) in both spatial and frequency domains. Therefore, it is used to represent images in multiple resolutions. The wavelet transform decomposes signals into basic functions that are dilations and translations of a single prototype wavelet function:
(1)f(x)=∑a=0∞∑b=0∞CabΨa,b(x),
where Ψa,b(x)=(1/a)Ψ((x-b)/a) are obtained by translating and dilating the wavelet function Ψ(x). The discrete wavelet transform coefficient Cab is the estimation of signal components in the time-frequency plane. This coefficient is obtained by this product:
(2)Cab=〈Ψa,b,f〉.According to Mallat’s pyramid algorithm [20], the signal is decomposed using a low-pass filter “L” and a high-pass filter “H” and each filter bank is then sampled at a half rate (1/2 downsampling) of the previous frequency. By repeating this procedure, it is possible to obtain wavelet transform of any order. In the case of an image, the filtering is implemented in a separable way by filtering the lines and columns. According to this algorithm, the original image can be transformed into four subimages. A comprehensive survey of the wavelet transform is provided in [21]. This pyramidal representation is not desirable for estimation or detection problems. In order to preserve the invariability by translation, the downsampling operation must be suppressed and the decomposition obtained is called a stationary wavelet transform [21]. In practice, the structure in cascade of the filter bank does not change; the only operation of downsampling is suppressed. In this work, we choose the Haar wavelet [19, 20], because it is orthogonal, has small support (for this reason, it provides good localization in space), and requires small computational complexity (linear with respect to the size of the input image).


## 2.1. The Watershed Transform

The watershed transformation is a morphological tool dedicated to image segmentation [4, 5, 13]. According to a flooding paradigm, this transformation considers each gray level (or color or texture) as topographic data as shown in Figure 2.Figure 2
3D profile of a gradient image.The principle is to first search out all the local minima where to perforate holes. Then, immerse the model into a lake at a uniform speed as the time. During the flooding stage, the water level slowly rises, covering an increasing amount of the image. To prevent the fusion of two catchment basins corresponding to different regional minima, a dam is erected between them at the points where they congregate each other. The watershed line corresponds to a series of this constructed dam and local extremes. For digital images, watershed constitutes closed contour of one-pix width [6]. Different approaches are employed to use the watershed principle for image segmentation. One of the most common watershed algorithms was introduced by Meyer [6]. This flooding process is performed on the gradient image. The implementation of this algorithm uses the hierarchical queue [5].This transformation is characterized by its fast computing and high accuracy in locating the weak edges of adjacent regions. But normally this will lead to an oversegmentation of the image. Therefore, many researchers [14–18] have proposed various methods for improvement, such as a preprocessing step before watershed transformation.


## 2.2. The Wavelet Transform

The wavelet transform [19] is a mathematical tool that can represent signal (image) in both spatial and frequency domains. Therefore, it is used to represent images in multiple resolutions. The wavelet transform decomposes signals into basic functions that are dilations and translations of a single prototype wavelet function:
(1)f(x)=∑a=0∞∑b=0∞CabΨa,b(x),
where Ψa,b(x)=(1/a)Ψ((x-b)/a) are obtained by translating and dilating the wavelet function Ψ(x). The discrete wavelet transform coefficient Cab is the estimation of signal components in the time-frequency plane. This coefficient is obtained by this product:
(2)Cab=〈Ψa,b,f〉.According to Mallat’s pyramid algorithm [20], the signal is decomposed using a low-pass filter “L” and a high-pass filter “H” and each filter bank is then sampled at a half rate (1/2 downsampling) of the previous frequency. By repeating this procedure, it is possible to obtain wavelet transform of any order. In the case of an image, the filtering is implemented in a separable way by filtering the lines and columns. According to this algorithm, the original image can be transformed into four subimages. A comprehensive survey of the wavelet transform is provided in [21]. This pyramidal representation is not desirable for estimation or detection problems. In order to preserve the invariability by translation, the downsampling operation must be suppressed and the decomposition obtained is called a stationary wavelet transform [21]. In practice, the structure in cascade of the filter bank does not change; the only operation of downsampling is suppressed. In this work, we choose the Haar wavelet [19, 20], because it is orthogonal, has small support (for this reason, it provides good localization in space), and requires small computational complexity (linear with respect to the size of the input image).


## 3. The Proposed Method

The proposed algorithm can be summarized into the following steps.(i)
Select a desired scale2∧3, and compute the stationary wavelet transform of the image up to this scale.(ii)
Create a reconstructed gradient image by combining the morphological gradient magnitude found at each scale, as shown in Figure3.(iii)
Localize the foreground objects using stationary wavelet reconstruction of the only high frequency part and threshold the obtained image.(iv)
Clean the edges of this binary image using morphological reconstruction opening and closing operators.(v)
Calculate the skeleton of this preprocessed image to obtain the perfect suitable markers.(vi)
Compute the background markers.(vii)
Superimpose these markers in the reconstructed gradient image.(viii)
The latest obtained image is given as input to the watershed transform.Figure 3
Overall scheme of the reconstructed gradient image.These steps are detailed next.


### 3.1. Multiresolution Morphological Gradient

The gradient magnitude [22] is used often to preprocess a gray-scale image prior to using watershed transform. The combination of the morphological gradients in different resolutions [23] is insensitive not only to noise but also to extract various fineness of the edges, because from high-resolution images, small details can be detected, while from low-resolution images, larger structures are detected. Thus, to obtain a reconstructed gradient image, the original image is transformed into a subimage using stationary wavelet transform. We retained only the low frequency part of this subimage. The morphological gradient magnitude is calculated for each resolution as
(3)G(I)=(I⊕B)-(IΘB),
where ⊕ and Θ stand, respectively, for dilation and erosion. “B” is the structuring element (3*4 (width*height) rectangular SE). The reconstructed gradient is obtained by combining the gradient magnitude images found at each scale such as
(4)RG(I)=Sup⁡(Gi(Ii))i=1…n.n is the coarsest resolution.


### 3.2. Creating Markers

The oversegmentation of the watershed algorithm can be reduced by finding markers for each object. However, automatic marker selection is a difficult process and in most applications human interaction is required. In addition, if inappropriate markers are selected, the final result changes dramatically. Based on the analysis of the characteristics of the trabecular bone image, the traditional function that detects the regional minimum as markers is not efficient in the trabecular bone. This method identifies many minimums in the same region, so we obtain an over-segmentation in this region after the segmentation. Therefore, a novel combination method is proposed to detect markers. Firstly, wavelet transform and morphology are used to locate the elements and get the subimages including the elements. The major processes of this step are described as follows.(i)
Get the input image.(ii)
Adopt “Haar” wavelet for the wavelet transform of the input image.(iii)
Set the low frequency component after wavelet decomposition as 0; the high frequency part does not change.(iv)
Wavelet reconstruction: get the reconstructed image.(v)
Threshold the reconstructed image using Otsu’s method [14] and binary image is obtained.Figure4 shows the localization of the object of interest.(a) Result of the reconstructed image after decomposition; (b) binary image after wavelet transform.
(a)
(b)Based on the investigation of the characteristics of the micro-CT trabecular bone image, we suggest computing the skeleton of the objects of interest as markers, provided that each object only contains one marker. In order to create foreground and background markers, a binary image enhancement is needed. So, firstly, we use the morphological opening by reconstruction followed by the morphological closing reconstruction to remove small particles due to noise. These operators [6] are listed as follows: Let I,J be the two images defined on the same domain and J≤I. B (circular SE) is the flat structuring element.(i)
Opening by reconstruction: the reconstruction ofI from J is obtained by iterating elementary geodesic erosions of J above I until stability is reached and denoted as
(5)ΦI(rec)(J)=∩n≥1E(n)(J),
where E(I)(J) can be obtained by iterating “n” elementary geodesic erosion and the geodesic erosion is defined as
(6)E(I)(J)=(JΘB)∪I∪standsforpointwisemaximum.(ii)
Closing by reconstruction: the reconstruction of “I” from “J” is obtained by iterating elementary geodesic dilations of “J” under “I” until stability is reached and denoted as
(7)γI(rec)(J)=∪n≥1D(n)(J),
where D(n)(J) can be obtained by iterating “n” elementary geodesic dilation and the geodesic dilation is defined as
(8)D(n)(J)=(J⊕B)∩I∩:standsforpointwiseminimum.Then, to get smooth edge foreground objects of the binary image, we use some morphology operations such as closing and dilating. These morphology operations fill the gaps between pixels and smooth their outer edges. The size of the gaps between pixels must be the maximum of structuring element size.Figure5 shows the result of binary image enhancement with morphology filters.(a) Result of the morphological opening closing by reconstruction image; (b) preprocessed binary image with morphological operations (circular SE).
(a)
(b)The foreground and background markers are created by calculating the morphology skeleton of the preprocessed binary image. In order to find the skeleton of trabecular network or medullar space in trabecular bone image, we use the algorithm of Gonzales and Woods [24] for thinning binary regions. Region points are assumed to have value 1, and background points value 0. Thinning is a morphological operation, somewhat like erosion or opening. It has-been designed as an iterative edge deleting items on the area. By considering the constraints that deletion of these points [24, 25] does not remove end points nor break connectivity nor cause excessive erosion of the region.The thinning [24] of an image I by a structuring element B is related to the hit-and-miss transform and it is defined as
(9)Thin(I,B)=I-(hit-and-miss(I,B)),
where the subtraction is a logical subtraction defined by
(10)X-Y=X∩notY.The hit-and-miss transform of an imageI by structuring elements B1 (circular SE, size 2) and B2 (circular SE, size 5) (they should not be having any overlapping element) is defined as
(11)I⊗(B1,B2)=(IΘB1)-(IcΘB2),
where Θ and Ic stand, respectively, for erosion and the complement of “I.”Figure6 shows the obtained foreground and background markers.(a) Foreground markers object; (b) background objects.
(a)
(b)


### 3.3. Watershed Segmentation

Before we apply the watershed transformation to the reconstructed morphological gradient, we impose foreground and background markers as minimums. Then, the flood will come off only from these marker pixels. The modified gradient with imposed minima is given by(12)Mod-grad=IMMin(grad(I),IMbIMf).The result image is given as an input to the watershed algorithm. Let us represent the watershed segmentation operator byWt(·); the output image is defined as
(13)J=Wt(Mod-grad).


## 3.1. Multiresolution Morphological Gradient

The gradient magnitude [22] is used often to preprocess a gray-scale image prior to using watershed transform. The combination of the morphological gradients in different resolutions [23] is insensitive not only to noise but also to extract various fineness of the edges, because from high-resolution images, small details can be detected, while from low-resolution images, larger structures are detected. Thus, to obtain a reconstructed gradient image, the original image is transformed into a subimage using stationary wavelet transform. We retained only the low frequency part of this subimage. The morphological gradient magnitude is calculated for each resolution as
(3)G(I)=(I⊕B)-(IΘB),
where ⊕ and Θ stand, respectively, for dilation and erosion. “B” is the structuring element (3*4 (width*height) rectangular SE). The reconstructed gradient is obtained by combining the gradient magnitude images found at each scale such as
(4)RG(I)=Sup⁡(Gi(Ii))i=1…n.n is the coarsest resolution.


## 3.2. Creating Markers

The oversegmentation of the watershed algorithm can be reduced by finding markers for each object. However, automatic marker selection is a difficult process and in most applications human interaction is required. In addition, if inappropriate markers are selected, the final result changes dramatically. Based on the analysis of the characteristics of the trabecular bone image, the traditional function that detects the regional minimum as markers is not efficient in the trabecular bone. This method identifies many minimums in the same region, so we obtain an over-segmentation in this region after the segmentation. Therefore, a novel combination method is proposed to detect markers. Firstly, wavelet transform and morphology are used to locate the elements and get the subimages including the elements. The major processes of this step are described as follows.(i)
Get the input image.(ii)
Adopt “Haar” wavelet for the wavelet transform of the input image.(iii)
Set the low frequency component after wavelet decomposition as 0; the high frequency part does not change.(iv)
Wavelet reconstruction: get the reconstructed image.(v)
Threshold the reconstructed image using Otsu’s method [14] and binary image is obtained.Figure4 shows the localization of the object of interest.(a) Result of the reconstructed image after decomposition; (b) binary image after wavelet transform.
(a)
(b)Based on the investigation of the characteristics of the micro-CT trabecular bone image, we suggest computing the skeleton of the objects of interest as markers, provided that each object only contains one marker. In order to create foreground and background markers, a binary image enhancement is needed. So, firstly, we use the morphological opening by reconstruction followed by the morphological closing reconstruction to remove small particles due to noise. These operators [6] are listed as follows: Let I,J be the two images defined on the same domain and J≤I. B (circular SE) is the flat structuring element.(i)
Opening by reconstruction: the reconstruction ofI from J is obtained by iterating elementary geodesic erosions of J above I until stability is reached and denoted as
(5)ΦI(rec)(J)=∩n≥1E(n)(J),
where E(I)(J) can be obtained by iterating “n” elementary geodesic erosion and the geodesic erosion is defined as
(6)E(I)(J)=(JΘB)∪I∪standsforpointwisemaximum.(ii)
Closing by reconstruction: the reconstruction of “I” from “J” is obtained by iterating elementary geodesic dilations of “J” under “I” until stability is reached and denoted as
(7)γI(rec)(J)=∪n≥1D(n)(J),
where D(n)(J) can be obtained by iterating “n” elementary geodesic dilation and the geodesic dilation is defined as
(8)D(n)(J)=(J⊕B)∩I∩:standsforpointwiseminimum.Then, to get smooth edge foreground objects of the binary image, we use some morphology operations such as closing and dilating. These morphology operations fill the gaps between pixels and smooth their outer edges. The size of the gaps between pixels must be the maximum of structuring element size.Figure5 shows the result of binary image enhancement with morphology filters.(a) Result of the morphological opening closing by reconstruction image; (b) preprocessed binary image with morphological operations (circular SE).
(a)
(b)The foreground and background markers are created by calculating the morphology skeleton of the preprocessed binary image. In order to find the skeleton of trabecular network or medullar space in trabecular bone image, we use the algorithm of Gonzales and Woods [24] for thinning binary regions. Region points are assumed to have value 1, and background points value 0. Thinning is a morphological operation, somewhat like erosion or opening. It has-been designed as an iterative edge deleting items on the area. By considering the constraints that deletion of these points [24, 25] does not remove end points nor break connectivity nor cause excessive erosion of the region.The thinning [24] of an image I by a structuring element B is related to the hit-and-miss transform and it is defined as
(9)Thin(I,B)=I-(hit-and-miss(I,B)),
where the subtraction is a logical subtraction defined by
(10)X-Y=X∩notY.The hit-and-miss transform of an imageI by structuring elements B1 (circular SE, size 2) and B2 (circular SE, size 5) (they should not be having any overlapping element) is defined as
(11)I⊗(B1,B2)=(IΘB1)-(IcΘB2),
where Θ and Ic stand, respectively, for erosion and the complement of “I.”Figure6 shows the obtained foreground and background markers.(a) Foreground markers object; (b) background objects.
(a)
(b)


## 3.3. Watershed Segmentation

Before we apply the watershed transformation to the reconstructed morphological gradient, we impose foreground and background markers as minimums. Then, the flood will come off only from these marker pixels. The modified gradient with imposed minima is given by(12)Mod-grad=IMMin(grad(I),IMbIMf).The result image is given as an input to the watershed algorithm. Let us represent the watershed segmentation operator byWt(·); the output image is defined as
(13)J=Wt(Mod-grad).


## 4. Experimental Results

The proposed segmentation algorithm is applied to a wide number of micro-CT trabecular bone images of human bone biopsy, taken from the femoral neck [2]. The results are subjectively analyzed. The trabecular network is effectively segmented from the medullar space. The markers are selected adaptively, depending on the characteristics of the trabecular bone image. They are employed to join adjacent regions because they should be corresponding to the same object and should eliminate the split of the wanted objects. Although the proposed algorithm prevents oversegmentation, it preserves the subfeatures of the objects required for getting morphological calculation results and reconstructed 2D slices. These segmentation results can be used in the feature extraction and classification stages to distinguish a normal bone from an osteoporotic bone. The segmentation results as produced by the proposed algorithm are shown in Figure 7.(a) Image markers; (b) watershed segmentation results; (c) overlay results on the original image.
(a)
(b)
(c)


### 4.1. Comparison with Other Techniques

To evaluate the performance of the presented method, we compare the proposed technique with other preprocessing techniques known to improve the watershed segmentation, namely, Marker-Controlled Watershed Segmentation [24], combined morphological-spectral unsupervised image segmentation method (VPSN) [11], and watersheds segmentation method using wavelets to enhance the image [16–22]. Results were analyzed in both quantitative metric and visual inspection.The segmentation results using the classical Marker-Controlled Watershed, the method combining the watershed with the dual-tree complex wavelet transform (VPSN), and the method using the wavelet to enhance the image after applying the watershed are shown in Figures8(c), 8(d) and 8(e). Comparing with Figure 8(b), it is obvious that our proposed strategy produces a much more accurate segmented image. Also, we applied our strategy to other low contrast images such as Figures 9 and 10.(a) Trabecular bone image. (b) Segmentation using our proposed technique. (c) Segmentation using classical Marker-Controlled Watershed. (d) Segmentation using combined morphological-spectral method (VPSN). (e) Segmentation using watershed and wavelet.
(a)
(b)
(c)
(d)
(e)(a) Inherently noisy spine image. (b) Segmentation using our proposed technique. (c) Segmentation using classical Marker-Controlled Watershed. (d) Segmentation using combined morphological-spectral method (VPSN). (e) Segmentation using watershed and wavelet.
(a)
(b)
(c)
(d)
(e)(a) Industrial circuit image. (b) Segmentation using our proposed technique. (c) Segmentation using classical Marker-Controlled Watershed. (d) Segmentation using combined morphological-spectral method (VPSN). (e) Segmentation using watershed and wavelet.
(a)
(b)
(c)
(d)
(e)Figures9(b) to 9(e) and Figures 10(b) to 10(e) show the segmentation results obtained by the proposed technique and the other proposed techniques known to improve the watershed segmentation cited above. In fact, visual inspection shows that the proposed technique produced the best result.Figure11 compares computational time (in seconds) of the four segmentation methods of the industrial circuit image.Figure 11
Computational time (in seconds) of the enhanced image by wavelet, VPSN, Marker-Controlled Watershed, and proposed methods of the circuit image.The proposed technique using the multiresolution to estimate the reconstructed morphological gradient image and the high frequency subbands to locate and compute the markers is significantly faster than the other methods.Furthermore, results were evaluated in quantitative points of view. The chosen quantitative metric was a percentage error obtained by computing the difference between the segmented image and the reference image according to the equation below:(14)E=Np((Bgt∪Bseg)-(Bgt∩Bseg))Np(Bgt),
where “Bgt” represents a binary image of the reference image, “Bseg” represents a binary image of the segmented image, and the “Np” represents the number of the pixels of the corresponding binary image.For example, a quantitative analysis (according to (14)) indicates that segmentation errors (in %) corresponding to Figures 8(b) to 8(e) are, respectively, 12.73, 17.56, 65.74, and 22.56% (reference image was obtained by manual delineation). The proposed technique produced the smallest segmentation error. Also the segmentation errors (in %) corresponding to Figures 9(b) to 9(e) are, respectively, 20.47, 27.50, 25.22, and 84.27%. Once more, the proposed technique produced the smallest error.


## 4.1. Comparison with Other Techniques

To evaluate the performance of the presented method, we compare the proposed technique with other preprocessing techniques known to improve the watershed segmentation, namely, Marker-Controlled Watershed Segmentation [24], combined morphological-spectral unsupervised image segmentation method (VPSN) [11], and watersheds segmentation method using wavelets to enhance the image [16–22]. Results were analyzed in both quantitative metric and visual inspection.The segmentation results using the classical Marker-Controlled Watershed, the method combining the watershed with the dual-tree complex wavelet transform (VPSN), and the method using the wavelet to enhance the image after applying the watershed are shown in Figures8(c), 8(d) and 8(e). Comparing with Figure 8(b), it is obvious that our proposed strategy produces a much more accurate segmented image. Also, we applied our strategy to other low contrast images such as Figures 9 and 10.(a) Trabecular bone image. (b) Segmentation using our proposed technique. (c) Segmentation using classical Marker-Controlled Watershed. (d) Segmentation using combined morphological-spectral method (VPSN). (e) Segmentation using watershed and wavelet.
(a)
(b)
(c)
(d)
(e)(a) Inherently noisy spine image. (b) Segmentation using our proposed technique. (c) Segmentation using classical Marker-Controlled Watershed. (d) Segmentation using combined morphological-spectral method (VPSN). (e) Segmentation using watershed and wavelet.
(a)
(b)
(c)
(d)
(e)(a) Industrial circuit image. (b) Segmentation using our proposed technique. (c) Segmentation using classical Marker-Controlled Watershed. (d) Segmentation using combined morphological-spectral method (VPSN). (e) Segmentation using watershed and wavelet.
(a)
(b)
(c)
(d)
(e)Figures9(b) to 9(e) and Figures 10(b) to 10(e) show the segmentation results obtained by the proposed technique and the other proposed techniques known to improve the watershed segmentation cited above. In fact, visual inspection shows that the proposed technique produced the best result.Figure11 compares computational time (in seconds) of the four segmentation methods of the industrial circuit image.Figure 11
Computational time (in seconds) of the enhanced image by wavelet, VPSN, Marker-Controlled Watershed, and proposed methods of the circuit image.The proposed technique using the multiresolution to estimate the reconstructed morphological gradient image and the high frequency subbands to locate and compute the markers is significantly faster than the other methods.Furthermore, results were evaluated in quantitative points of view. The chosen quantitative metric was a percentage error obtained by computing the difference between the segmented image and the reference image according to the equation below:(14)E=Np((Bgt∪Bseg)-(Bgt∩Bseg))Np(Bgt),
where “Bgt” represents a binary image of the reference image, “Bseg” represents a binary image of the segmented image, and the “Np” represents the number of the pixels of the corresponding binary image.For example, a quantitative analysis (according to (14)) indicates that segmentation errors (in %) corresponding to Figures 8(b) to 8(e) are, respectively, 12.73, 17.56, 65.74, and 22.56% (reference image was obtained by manual delineation). The proposed technique produced the smallest segmentation error. Also the segmentation errors (in %) corresponding to Figures 9(b) to 9(e) are, respectively, 20.47, 27.50, 25.22, and 84.27%. Once more, the proposed technique produced the smallest error.


## 5.  Conclusions

A novel method using morphological operations, wavelet, and watershed segmentation is proposed to automatically detect the trabecular network from micro-CT trabecular bone images. The proposed method does not require human intervention. The experiment shows that the proposed algorithm produces accurate boundaries of the region of interest. In addition to that, the algorithm is implemented automatically leading to efficient computation work. We should also notice that the segmentation result depends on the quality of the gradient image and the suitable markers.



---

*Source: 891950-2014-07-13.xml*