Abstract
Purpose. Most applications in the field of medical image processing require precise estimation. To improve the accuracy of segmentation, this study aimed to propose a novel segmentation method for coronary arteries to allow for the automatic and accurate detection of coronary pathologies. Methods. The proposed segmentation method included 2 parts. First, 3D region growing was applied to give the initial segmentation of coronary arteries. Next, the location of vessel information, HHH subband coefficients of the 3D DWT, was detected by the proposed vessel-texture discrimination algorithm. Based on the initial segmentation, 3D DWT integrated with the 3D neutrosophic transformation could accurately detect the coronary arteries. Results. Each subbranch of the segmented coronary arteries was segmented correctly by the proposed method. The obtained results are compared with those ground truth values obtained from the commercial software from GE Healthcare and the level-set method proposed by Yang et al., 2007. Results indicate that the proposed method is better in terms of efficiency analyzed. Conclusion. Based on the initial segmentation of coronary arteries obtained from 3D region growing, one-level 3D DWT and 3D neutrosophic transformation can be applied to detect coronary pathologies accurately.
1. Introduction
Efficient and automatic image segmentation methods are useful for the isolation and visualization of vessels in computed tomographic angiography (CTA). There are many proposed methods for the segmentation of vessels [1–14]. A vessel filter [1] can be used to enhance tubular structure; however, it cannot address the problem of the image force and veins, which can lead to a narrowed or broken segmentation of vessels. Parametric shape models [2–5] do not directly allow for the detection of topological changes, and they usually obtain a seriously narrowed segmentation in the neighborhood of a branch point in the vessel. Level-set approaches [6–13] are computationally expensive. They also suffer from leakage at places where the intensity gradients of the edges are relatively weak and are very sensitive to the placement of the initial contour of the propagating front. Metz et al. [14] used the minimum cost path of the specified start and end points in vessel to detect the coronary arteries centerline. This is not an automatic method; detecting branches is difficult. Friman [15] proposed multiple hypothesis template tracking, which follows the direction of centerline obtained in advance. However, it is difficult to detect small branches and vessels by using this method.
In this study, we propose a new method for automatically and correctly segmenting coronary arteries from CTA data sets. In image preprocessing, we detected the aorta automatically by using methods proposed in the literature [5, 16]. The proposed coronary arteries segmentation method is summarized as follows. First of all, we automatically obtained the seed point of a 3D region growing by the difference between the two adjacent slices due to the small changes of the aorta between two adjacent slices. Next, 3D region growing was applied to initially search for the probable location of coronary arteries, which was then dilated by 3 voxels. Based on the dilation of the probable location, we detected the coronary arteries accurately by applying the 3D discrete wavelet transformation (DWT) and 3D neutrosophic transformation to the CTA volume. The location of vessel information, in HHH subband coefficients, was detected by the proposed vessel-texture discrimination algorithm. Accordingly, HHH subband coefficients were used, which were characterized and classified by - operation and -means clustering. Finally, the proposed method was tested on several CTA data sets, and the experimental results indicated that the proposed method had a good performance.
The rest of this study is organized as follows. Section 2 reviews some preliminaries, and Section 3 uses 3D region growing and 3D DWT to propose a new method for segmenting coronary arteries. Section 4 contains the experiments and discussion, and the conclusions are drawn in Section 5.
2. Preliminaries
In this section, we will briefly introduce the concepts of DWT and provide an overview of some fundamental mathematical concepts that are used in this study.
2.1. Region Growing
Region growing is a simple, well-developed, region-based image segmentation technique [17]. It postulates that neighboring voxels within the same region have similar intensity values and is also classified as a voxel-based image segmentation method since it involves the selection of initial seed points. In other words, this method of segmentation examines neighboring voxels of initial seed points and determines whether neighboring voxels should be added to the region. Consequentially, the general concept of region growing is to group voxels with the same or similar intensities to one region according to the given seed points and a homogeneity criterion.
2.2. Discrete Wavelet Transform
Wavelet transform is obtained by a single prototype function which is regulated with scaling and shift parameters. To construct , a scaling function is determined. The discrete normalized scaling and wavelet basis functions are defined as where and are the dilation and translation parameters. Orthogonal wavelet basis functions not only provide a simple method to calculate coefficient expansion but also span in signal processing. As a result, signal can be expressed as a series expansion of orthogonal scaling functions and wavelets. More specifically, where and are the low-pass and high-pass coefficients, respectively; is an integer to define an interval on which is a piecewise constant. The two-scale equations for scaling and wavelet basis function are given as follows: where . The coefficient in (3) has to meet several conditions for the set of the wavelet basis function to be unique and orthonormal and have a certain degree of regularity.
The coefficients and play a very crucial role in a given DWT. Performing the wavelet transformation does not require the explicit forms of and but only depends on and . The final output of the wavelet decomposition includes a set of - wavelet coefficients. One method to implement DWT is to use a filter bank that provides perfect reconstruction. DWT involves local analysis of frequency in space and time domains, and it provides multiscale image details step by step. If the scale becomes smaller, every part becomes more accurate, and ultimately all imaging details can be focalized accurately. If DWT is applied to a volume, it will produce the highest-frequency, middle-frequency, and lowest-frequency parts. Figure 1 shows the results of applying 3D DWT to a volume, which includes eight parts: LLL, LLH, LHL, LHH, HLL, HLH, HHL, and HHH. The lowest-frequency and highest-frequency parts are LLL and HHH, respectively [16, 17].
2.3. -Means Clustering
-means clustering is a method of cluster analysis which aims to partition observations into clusters in which each observation belongs to the cluster with the nearest mean. Given a set of observations , where each observation is a -dimensional real vector, -means clustering aims to partition the observations into sets so as to minimize the within-cluster sum of squares: where is the mean of points in .
3. The Proposed Segmentation Method
In order to segment coronary arteries accurately from CTA data sets, 3D region growing was initially applied to search for the probable location of the coronary arteries. Next, we used the 3D DWT and 3D neutrosophic transformation to accurately detect the coronary arteries.
3.1. Initial Segmentation of Coronary Arteries
This section discusses the initial segmentation of coronary arteries using 3D region growing. In image preprocessing, we found the aorta automatically by using methods proposed in the literature [5, 16]. The selection of the seed point of 3D region growing was initially made to check which slice began the information of the coronary arteries. Due to the small changes in aorta area between the two adjacent slices, we automatically obtained the seed points by the difference between the two adjacent slices.
As shown in Figure 2, we use the difference between the two adjacent slices (a) and (b) to automatically obtain the seed points bounded by the blue line in (c) which indicates the boundary of coronary arteries. The boundary of coronary arteries denotes the high-frequency subband in (2) when comparing vessel lumen and background.
(a)
(b)
(c)
Since coronary arteries do not exhibit abrupt intensity changes along their centerline [4], a rough tubular mask of coronary arteries can be easily constructed by 3D region growing. We chose a 26-connected neighborhood for our adjacent pixel relationship, and then the 3D region growing method was applied with a set of prespecified seed voxel(s) and grown from these seeds by merging neighboring voxels whose properties were most similar to the premerged region. The homogeneity criterion was defined as the difference between the intensity of the candidate voxel and the average intensity of the premerged region. The selection of the seed point was initially intended to check which slice began the information of the coronary arteries. Next, the homogeneity criterion was applied to group voxels with the same or similar intensities into one region. If the homogeneity criterion was satisfied, the candidate voxel was merged with the premerged region. The process was repeated until no more voxels were assigned to the region, and then the number of all merged voxels was calculated. In order to avoid leakage, the total number of merged voxels was limited to 12000; otherwise, 3D region growing was restarted by automatically using an improved homogeneity criterion. Finally, the initial segmentation of the coronary arteries in a volume was completed.
3.2. Vessel-Texture Discrimination
According to Parseval’s theorem, the energy in a signal is given as follows [16, 17]: This equation implies that the energy of a signal is the summation of low-frequency and high-frequency coefficients. DWT is a good analytic tool for image texture analysis or line-based patterns [18–21]. Since a vessel is a type of 3D line-based pattern in CT volume, we used - of these wavelet coefficients to find the energy of line-based patterns. - was defined as where the vector was the wavelet coefficients of a frequency channel. The searching algorithm is summarized as follows.
Algorithm 1. (1) Transform a given vessel volume into frequency channels by a specified number of decomposition levels. We usually set the number to one in the first search.
(2) Use (7) to calculate the average - of each channel and maximum of for the vessel volume.
(3) If the maximum of was significantly greater than another channel’s , the search was stopped. Otherwise, the number of decomposition levels was increased followed by a repeat of step 1.
By using the above algorithm, we observed that the most significant information of the vessel texture often appeared in the high frequency channels. Thus, we used the subband HHH to detect vessels in this study.
3.3. Accurate Detection of Coronary Arteries
The initial segmentation of coronary arteries in a volume was completed by 3D region growing as described in Section 3.1. Since region growing is a simple region-based image segmentation method, it was only used to search for the initial location of the coronary arteries. We then accurately detected the coronary arteries by applying DWT to each slice in the volume, as described in this subsection.
First, the initial location of the coronary arteries obtained from 3D region growing was dilated by 3 voxels. Next, we used the Haar wavelet bases in (1) to transform the host images into the orthogonal DWT domain by one-level decomposition. Only HHH subbands were employed for further processes, because most of the information on the coronary arteries and boundaries were in the HHH subbands. We calculated the mean energy using coefficients of HHH subbands in a local window as follows: where Next, the subbands HHH were characterized by 3 membership sets , , and . Consider where That is, a pixel could be represented as a neutrosophic domain which means the pixel is true, false, and uncertain, where varies in , varies in , and varies in . In order to reduce the uncertainty , - operation was employed as follows: where the parameter is a positive number and The inverse DWT was then applied to obtain a new volume which possessed the true subset. Finally, we applied -means clustering () in (2) to differentiate vessel lumen, vessel boundary (true subset), and background. The true subset was retained, respectively.
4. Experiments and Discussion
To test the proposed method, CTA volumes obtained from a CT system were segmented for coronary arteries. The slice thickness was 0.625 mm and the volume was in different data sets. The window size was set to 3 which was enough to capture the local texture characteristics. The parameter was set to 0.2. We tested 20 data sets, most of which were segmented successfully except for a few small branches that were lost in 2 of the data sets due to the local failure in region growing. To evaluate the performance of our segmented coronary arteries, we compared our results with that obtained from the ground truth values obtained from the commercial software from GE Healthcare and the level-set method. We used the overlapping metric (OM) and Hausdorff distance () to analyze the efficiency of each method.
4.1. The Segmenting Efficiency of 2D Imaging
Figure 3(b) shows the results of 6 slices in 1 CTA volume obtained using the proposed method. The areas bounded by the red line are coronary arteries. The segmenting efficiency was compared with the manually delineated ground truth data in Figure 3(a) by using OM which was defined as where indicates the pixels/voxels of the segmented coronary arteries. The OM was close to 1 when the segmentation was well matched to the reference ground truth and approached zero when the results had no similarity to the reference.
(a) Ground truth data
(b) The segmented results
In the 6 slices in Figure 3(b), the segmentation results showed that the proposed method detected coronary arteries accurately. As shown in Table 1, the average OM of the proposed method was 0.96.
4.2. The Segmenting Efficiency on a 3D Volume
The first focus of the comparison was the correctness of the 4 main branches: the right coronary artery (RCA), the left anterior descending artery (LAD), the circumflex (CRX), and the first diagonal artery (DA). Figures 4 and 5 show the detected coronary arteries in the CTA volume obtained from the GE Healthcare and the proposed method, and the proposed method segmented these 4 main branches correctly compared to the 4 main branches in Figure 5.
Another test of performance is the correctness of the remaining branches. Due to the multiresolution of the DWT, each subbranch in Figures 3 and 5 was correctly segmented. Furthermore, the coronary arteries obtained from the proposed method were much better than those obtained from the level-set method, which had several leakages as seen in Figure 6.
In addition to the OM, the difference between the segmented vessel surface and the manually delineated ground truth data was measured by the Hausdorff distance [22] which was defined as follows: where and are the vertices of the mesh surfaces of the arteries corresponding to the segmentation results and the ground truth and measure the Euclidean distance between points and belonging to vertices and . Table 2 lists the mean OM and mean Hausdorff distance for the proposed method and the method of Yang et al. [11]. The results show that the proposed method was much better than that of Yang et al. [11] in terms of both OM and Hausdorff distance.
4.3. Diameter Measurement
In this subsection, we computed the diameters of segmented coronary arteries using the proposed method and the method of Yang et al. Many efficient algorithms have been proposed to extract the tube centerline. We applied the algorithm proposed by Lee et al. [23] to extract the centerline of the segmented coronary arteries. Using these extracted centerlines, we obtained the cross sections of the segmented coronary arteries, as shown in Figure 7. By computing the area of each cross section, the diameter was estimated as follows:
(a)
(b)
Table 3 shows the estimated diameter of the same cross section for the proposed method and the method of Yang et al. [11]. The diameter of the proposed method was closer to the diameter in ground truth data at the same cross section than to that obtained from the method of Yang et al.
4.4. Experimental Environment and Execution Time
The proposed method was implemented in MATLAB (R2011a) on a standard specification PC with a 3.2 GHz CPU and 12 GB RAM. The average execution time was 58 seconds to extract the entire coronary tree, compared to approximately 47 seconds for the method by Yang et al. for the same process.
5. Conclusions
Accurate extraction of coronary arteries is important to assess artery lesions in clinical practice. In this study, we propose a novel method to segment coronary arteries automatically. Based on the initial segmentation obtained from 3D region growing, one-level 3D DWT and 3D neutrosophic transformation were applied to detect coronary arteries accurately. The location of vessel information, in HHH subband coefficients of DWT, was successfully detected by the proposed vessel-texture discrimination algorithm. Accordingly, the HHH subband coefficients were used and characterized and classified by 3D neutrosophic transformation and -means clustering. The experimental results verify the efficiency of the proposed method.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This study was supported by the National Taiwan University Hospital and the National Science Council, Taiwan, under the NSC Grant: NSC 98-2221-E-002-098-MY3.