An Approach for Pulmonary Vascular Extraction from Chest CT Images

Pulmonary vascular extraction from chest CT images plays an important role in the diagnosis of lung disease. To improve the accuracy rate of pulmonary vascular segmentation, a new pulmonary vascular extraction approach is proposed in this study. First, the lung tissue is extracted from chest CT images by region-growing and maximum between-class variance methods. Then the holes of the extracted region are filled by morphological operations to obtain complete lung region. Second, the points of the pulmonary vascular of the middle slice of the chest CT images are extracted as the original seed points. Finally, the seed points are spread throughout the lung region based on the fast marching method to extract the pulmonary vascular in the gradient image. Results of pulmonary vascular extraction from chest CT image datasets provided by the introduced approach are presented and discussed. Based on the ground truth pixels and the resulting quality measures, it can be concluded that the average accuracy of this approach is about 90%. Extensive experiments demonstrate that the proposed method has achieved the best performance in pulmonary vascular extraction compared with other two widely used methods.


Introduction
At present, computed tomography (CT) has become the most common imaging modality for the diagnosis of lung disease. In analysis of chest CTscans, pulmonary vascular extraction is often required before proceeding to diagnose disease and also is an important step in the treatment planning, and follow-up of lung diseases [1]. In computer-aided diagnosis of lung disease, pulmonary vascular extraction can reduce ambiguities and improve nodule detection performance [2] and is used to aid detection of other pulmonary anatomical structures such as airway tree, pulmonary tissue, lung vein, and artery [3][4][5]. Due to the complexity of the anatomical structure of the pulmonary vascular and the influence of other anatomical structures having similar intensities (e.g., tumor nodules and dense lesion) [6], there is a still more complicated task of extract pulmonary vascular from large scale chest CT images, especially for the small vascular.
For the large amount of chest CT images, manual extraction is an extremely time-consuming and tedious task for doctors. e computerized semiautomatic and automatic methods may be helpful to reduce the doctor's effort. To address these issues, various methods have been extensively studied in recent years. Xu et al. [7] propose novel stagewise convolutional networks, followed by an orientation-based region-growing method, which learn discriminative features of pulmonary vessels automatically in a stage-by-stage manner by stagewise convolutional networks. Orkisz et al. [8] reported a vascular trees segmentation method by variational region growing. is process is performed within a lung mask, where the airways and bronchial walls were previously eliminated by adaptive multiscale morphological operations. Anna [9] introduced a 3D approach for segmentation of pulmonary vascular tree from CT thoracic scans based on author's experience in airway tree segmentation. Lai et al. [10] proposes an automatic integration segmentation approach for the vascular trees in low-dose CT scans. Zhu et al. [11] proposed a vascularity-oriented level set algorithm for pulmonary vessel segmentation in imageguided intervention therapy. Zhao et al. [12] proposed a vessel segmentation method is proposed for lung images based on a random forest classifier and sparse autoencoder features. Chen et al. [13] proposed a method to simultaneously and separately segment the pulmonary nodules and blood vessels. First, a line structure enhancement (LSE) filter and a blob-like structure enhancement (BSE) filter were used to augment the initial selection of vessel regions and nodule candidates, respectively. en a front surface propagation (FSP) procedure was employed for precise segmentation of blood vessels and nodules. Buelow et al. [14] presented an automated method for the extraction of the pulmonary vessel tree from multislice CT data by a seed point-based front-propagation algorithm. Rudyanto et al. [15] presented an annotated reference dataset containing 20 CT scans and proposed nine categories to perform a comprehensive evaluation of vessel segmentation algorithms from both academia and industry. is dataset is used in the VESSEL12 challenge held at the International Symposium on Biomedical Imaging (ISBI) 2012. e currently most general and extensive vascular segmentation method can be found in [16], with a highly detailed categorization of the existing work. However, because the pulmonary vessel trees are very complex and have a huge number of branches, the overall number of works on pulmonary vascular extraction is very limited for the computer-assisted diagnosis of pulmonary disease and the major methods for pulmonary extraction are still to be developed [7,17,18].
Coming forward to meet these needs, a new approach for the extraction of pulmonary vascular from chest CT images is proposed in this work. e proposed approach is based on author's previous work in airway tree extraction [19]. e main idea of the proposed method is to extract the whole lung region, enhance the extracted lung region, and then extract the pulmonary vascular.
is paper is organized as follows. First, in Section 1 a brief review of existing approaches to segmentation of the pulmonary vascular from chest CT images is introduced. Section 2 provides a detailed description of the proposed extraction approach. e image datasets, experiment software, and results of vascular extraction from chest CT image using the introduced method are presented and discussed in Section 3. Finally, Section 4 concludes the study.

Methods
e method of this study mainly includes the following steps: extract lung region, enhance the extracted lung region, and extract the pulmonary vascular. e flowchart of the method in this study is shown in Figure 1. e specific algorithms for each step are described in the following sections.

Lung Region Extraction Method.
e anatomical structure of the lung in chest CT image is complex, as shown in Figure 2, if the pulmonary vascular is directly extracted in original CT images, it is difficult to remove the interference of other anatomical structures, such as bone, heart and muscle. erefore, to segment accurately the pulmonary vascular, it is necessary to extract the lung region to remove other anatomical structures.

Extract Lung Tissue with Region-Growing and Maximum Between-Class Variance Methods.
e chest CT image is divided into two categories of lung tissue and other anatomical structure, in which grayscale values are significantly different in the CT image. Region-growing method can segment effectively and correctly a specific set of grayscale values in an image based on the grayscale difference between pixel of that image. e region-growing method is an iterative image segmentation method with three elementary parts: seed points selection, definition of similarity, and criteria for convergence to terminate the iterative process [20]. e regions are grown from these seeds to adjacent points depending on a region membership criterion (e.g., grayscale intensity). Keep examining whether the adjacent points of seeds should be classified into the seed points until the criterion is not met any more [21].
In the previous work, we have used the region-growing method to extract airway from chest CT images [19]. e extracted airway result showed that this method is simple and effective. erefore, we also use the region-growing method to extract the lung tissue in CT images in this study. is method requires seed points of lung tissue from a particular slice image, which is the middle slice of the patient's chest CT image in this work.
First, the seed points are extracted by the method of maximum between-class variance, which can divide the original image into two parts by using threshold: foreground and background. In the chest CT image, the background is lung tissue and the foreground is the other anatomical structure. T is set to be the segmentation threshold of the foreground and background, the number of foreground points accounts for ω 0 of the image, the average grayscale value is u 0 , the number of background points accounts for ω 1 of the image, the average grayscale value is u 1 , and the total grayscale value of the image is as follows: (1) e variance of foreground and background images: Take u into formula (2) for calculating the variance of the two classes as follows: e maximum threshold T max is obtained by using formula (3). e pixels of the image less than threshold T max are extracted as the seed points of lung tissue.
Second, the similarity definition is used to determine whether the unmarked pixels of image are added to the detection region. is definition refers to the difference of image intensity between the adjacent pixels. e similarity condition is formulated [22] as follows: e unmarked voxel (x k , y k ) in the 26-adjacent N 26 can be added to the lung tissue region, if the difference between the grayscale value of I(x k , y k ) and seed voxel I(x, y) is less than the given threshold θ. And this voxel (x k , y k ) is added to the seed queue voxel for next iteration. is process is shown in Figure 3.
Since this method relies on the difference between the grayscale value of the voxels in the 26-adjacent pixels, the given threshold θ plays an important role in this process. e extraction lung tissue regions are different with different given threshold. e threshold θ is calculated according to the following formula: where max I lung tissue (x, y) and min I lung tissue (x, y) are the maximum and minimum gray value and of the pixels of the extracted seed points by maximum between-class variance method, respectively. e extracted result of the lung tissue is shown in Figure 4, which shows that there are many holes in the extracted lung region. e reasons are that the grayscale value of the vascular is different from the lung tissue.

Fill Holes of Extracted Lung Region with Morphological
Operation. In this work, morphological operation method is used to fill the holes of lung tissue on the initial extraction region.
First, the contours of the lung are extracted in the all 2-D transverse images by 8-adjacent contour tracing method.
ere may be several contours in an image. e contours and the background are set to 1 and 0, respectively. And the set A is defined to represent the contours and the background pixels of the image. We formulate the filling "holes" function as where X 0 � p, p is any point in the region, B is the 4neighbourhood structural elements, ⊕ is the dilation operator, and A c is the complement set of A. At the end of the iteration, X k � X k−1 . e union set of A and X k is the region filling the holes. e extracted region is expanded iteratively by the structural element. e number of iterations is manually based on experience. en, the extracted region is corroded by the same structural element and the number of iterations of the expanding operation. e final result of the extracted lung region after filling holes is shown in Figure 5.

Vascular Enhancement Based on Fractional Differential.
e actual chest CT image contains noise, which cause the edge of the pulmonary vascular of the image to be unclear. If the pulmonary vascular are directly extracted from the original CT image, it is easy to cause inaccurate extraction of small pulmonary vascular. erefore, the fractional order differential operator is proposed for enhancing vascular region in the extracted lung tissue(as shown in Figure 4) region before extracting pulmonary vascular.
e Grümwald-Letnikov (G-L) is used to define the numerical implementation of fractional differential of image [23][24][25]. e period of the unary signal w(t) is t ∈ [a, t], and the signal period [a, t] is equally divided by a certain time interval (such as h � 1, 2, ...), then n � [ and thus, the differential expressions with fractional differentials can be obtained: For the CT image, the differential expressions with fractional differentials can be obtained as follows: e above differential expression constructs 3 * 3 or 5 * 5 differential operators to process the image to enhance the edge of the pulmonary vascular. Since the differential operator is isotropic, it does not enhance the small pulmonary vascular information with the characteristic direction. us, a fractional differential operator template can be constructed according to formulae (8) and (9). e operator has the following effect: it has little effect on the smooth region, where the surrounding gray value changes little; and the operator plays an enhancement role for the pixels whose surrounding gray value changes greatly.

Extract lung region
Enhance the extracted lung region Extract the pulmonary vascular Pulmonary vascular Chest CT images By analysing the characteristics and distribution of pulmonary vascular in the chest CT images, it is known that the pulmonary vascular are cylindrical and their directions are mostly four diagonal directions. us, a 3 × 3-differential operator template is constructed to enhance the region of vascular, as shown in Figure 6.   Since the operator is not symmetric, the following method should be used: the coordinates of the "1" of the operator template should coincide with the coordinates of the pixel to be subjected to the fractional differential operation. e templates in these four directions are convolved with the lung tissue region, and then the convolution results in each direction are weighted, and the results in the four directions are summed. e image after pulmonary vascular enhancement can be expressed by formula (10): e lung region is processed using the above 3 * 3 fractional differential enhancement template. e different order can obtain different enhancement results. erefore, we first study the experimental results of different fractional orders for chest CT image, which is shown in Figure 7. e results show that the pulmonary vascular in the CT image are enhanced as the differential order increases. At the same time, the result of small pulmonary vascular enhancement between 0.1 and 0.2 is the best. Although the pulmonary vascular details are enhanced when the order is greater than 0.2, large pulmonary vascular are also suppressed, and the contrast between the large pulmonary vascular and other lung tissues is lowered. us, 0.2 fractional order is used to enhance the lung region of CT image in this work.

Extract Pulmonary Vascular by FMM.
e complete lung region is extracted by the method of Section 2.1 and enhanced by fractional differential method. In this section, the pulmonary vascular is extracted with the fast marching method in the lung tissue region. e fast marching method originates from the solution of the equation of the distance function. Suppose there is a curve φ moving in one direction. e FMM can be described as a family of schemes for computing the evolution of fronts.
ings become interesting when the front evolves over time. In the context of the fast marching method, the speed does not have to be the same everywhere but the speed must always be non-negative. As a given point, the motion of the front is described by the equation known as the Eikonal equation, which can be expressed as [13]: is is the general form of the equation of distance function. Godunov [26] gives a solution: In this solution, D − ij and D + ij are backward difference and forward difference operators, i and j represent the two adjacent points of the curve φ, respectively.
If ∇T is approximated as a first-order finite difference operator, then formula (12) can be written as In that case, A distance value for image must be found to estimate the length of the gradient ‖∇T‖ is equal to 1/F. e following formula (13) is proposed to solve ‖∇T‖: e D i,j , D i,k , D j,k represent the gradient value in the three directions of Cartesian coordinate system, respectively. i, j, and k represent the three adjacent points.
Considering the case of adjacent points of threedimensional CT images (as shown in Figure 8), formula (14) can be written as follows: where V A is the unknown distance value and V B , V C , V D , V E , V F , V G are the distance values at the neighbouring voxels. e distance value is the gray value of the pixels of CT images. For speeding up in homogeneous regions of the image and slowly in regions with a high value of image gradient. e F of the chest CT volume data, be considered as gradient figure, is computed by the exponentiated gradient algorithm as follows: e solution process of FMM is described as follows: Step 1. All the points are categorized into three categories: the processing points are the seed points (in the initial step, the number of the seed points is one; during the processing, the seed points may be many points); the boundary points are the points near the processing points; the pending points are all the remaining points.
Step 2. e arrival time T(x) is calculated from the initial seed point to the boundary points and sort them according to the time from small to large. Step 3.
e boundary points, less than a certain arrival time threshold, are marked as the next processing points of as the next processing points. Other boundary points are marked as the background points. e nearby points of the next processing points are marked as the new boundary points.
e arrival time T(x) of new boundary points are calculated and sorted. If the next processing points are empty, the processing is finished. Otherwise, it returns to (2). e above processing can be considered that the curve spread along an equal timeline until it arrives at the front of a wall.
In this work, the full lung region of the chest CT images is extracted by the method of Section 2.1 and enhanced by fractional differential method. en the points of pulmonary vascular are extracted with the maximum betweenclass variance method in the full lung region of the middle of slice of the chest CT images. In this case, the foreground is the pulmonary vascular. ese extracted points are selected as initial seed points. And the minimum and maximum grayscale values of these pulmonary vascular points are calculated as the limited grayscale value threshold. According to the limited grayscale value threshold and the gradient field figure F, the seed points begin to spread through the lung region of images using the fast march method as the above solution process. e points of pulmonary vascular are extracted completely when the spread processing stops.

Results and Discussion
e 3 chest CT images dataset were used in this experiment with resolution of 512 × 512, slice thickness less than 1.5 mm and slice number more than 350, which are from a hospital. In order to protect patient privacy information, the image dataset hides the hospital's and patient's name. e ground truth is manually drawn by the doctor for each database.
is method was implemented in Matlab 2015 on a PC with 4 Intel© Core is-i7-6700U CPUs 2.60 GHz, 8 GB DDR4 RAM and NVIDIA GeForce 940 Mx GPU with 2 GB video memory.
We will quantitatively analyse the results of the method for extraction of pulmonary vascular. By comparing the extraction results of the proposed method with the ground truth image, the accurate rate, the leakage rate and over rate of the proposed method are calculated to quantitatively evaluate the. e ground truth image, be manually drawn by the professional doctors on the original CT image, is provided in VESSEL12 challenge. e evaluation criteria for pulmonary vascular segmentation are as follows: where T and S is the number of pulmonary vascular pixels in the ground truth image and the extraction result by this method, respectively. In order to observe the experimental results, the comparison results of the 3 datasets between the method in this paper and the ground truth is shown in Figure 9. e results show that the method can completely extracted the pulmonary vascular, even the extraction results of small pulmonary vascular are also very good, and cause a little over segmentation in the edge of lung region. e method is compared with the level set [11] and regional growth algorithm [8], which are widely used to extract pulmonary vascular from medical image. e accurate rate, the leakage rate, and the over rate of each segmentation method are calculated by comparing with the ground truth image.
e experimental results of an image dataset of different segmentation methods are shown in Figure 10. e different places of these results are shown by yellow and blue circles in the result images. e results show that the segmentation results of these four methods are equally good, but this method is a better accurate segmentation of the small pulmonary vascular and less over segmentation than the other methods. e blue circle of Figure 10 shows that this method can extract more the small pulmonary vascular. Meanwhile, the yellow circle of Figure 10 shows that this method is less over segmentation results than the other methods. Table 1 and Figure 11 are the statistical analysis of the segmentation results of the three methods. e average of accurate rate of this method is 91.62%, which is better than the other two methods. And the average of leakage rate and over rate of this method is 8.38% and 3.43%, respectively. It is less than the other two methods.
In general, the running time of segmentation method is related to the slice thickness and the slice number. e average running time of this method is about 2 minutes,

Conclusions
is study proposes a method for segmenting pulmonary vascular from chest CT images. e lung tissue region is extracted from chest CT images by region-growing method and maximum between-class variance.
en, the holes of the extracted region are filled by morphological operation. e points of the pulmonary vascular of the middle slice of the CT images are extracted as the original seed points. Finally, the seed points are spread throughout the lung tissue based on the fast-marching method to extracted the pulmonary vascular with the gradient image. e experimental results show that this method can completely extracted the big and small pulmonary vascular and the accuracy is better than the level set method and region-growing method, and the leakage rate and over rate of this method is less than the two other methods. Furthermore, the running time of this method is about 2 minutes to meet the requirements of clinical realtime operation. is method may be used to promote the  Data Availability e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.