Extraction of pulmonary Trachea by dynamic tubular edge contour algorithm

Introduction

Lung cancer is one of the most common malignant tumors, and the number of patients has been increasing year by year (1). At present, the method used to accurately judge malignant lung tumors is pathological examination, which requires tissue samples of tumor lesions obtained through a puncture operation. The puncture needle needs to penetrate muscle and lung tissue to reach the tumor, and avoid the pulmonary trachea and blood vessels as much as possible, otherwise serious complications such as pneumothorax, pulmonary hemorrhage and hemothorax can occur, resulting in the death of patients (2,3). At present, lung puncture surgery is mainly performed manually by doctors referring to CT images, resulting in high levels of complications and death. Data surveys show that pneumothorax and pulmonary hemorrhage are the main complications, with the morbidity of pneumothorax, pulmonary hemorrhage and puncture death reaching 10–40%, 26–33%, and 1% respectively. With the development of computer technology and image technology, lung CT images can be used to extract pulmonary mask using random forest algorithm, pulmonary trachea and vessels to create a three-dimensional image as a visual reference for lung puncture procedures, guiding the puncture needle to actively avoid the pulmonary trachea and blood vessels and reducing the morbidity and mortality rates associated with lung puncture (4).

The commonly used tracheal extraction methods can be divided into five categories: (I) morphological method; (II) regular tracheal detection method; (III) template matching method; (IV) machine-learning method; and (V) shape analysis method (5). Due to volume effects and image noise (artifact) and other reasons, all of the commonly used extraction methods generally have the phenomenon of bronchial loss. Of them, the three-dimensional region growth algorithm based on morphology (I) is one of the more mature methods, which can extract the grade 7 lung bronchus, but there is still the phenomenon of trachea loss, and the selection of seed points needs to be high, otherwise leakage will occur. Therefore, there is great practical and research significance to being able to extract the complete lung and trachea structure, including the high-level lung and trachea, based on CT images, making it a hot topic in the field of computer vision and image processing.

Based on the characteristic “tubular structure” of the pulmonary trachea in CT images, we present an algorithm based on dynamic edge contour for use in axial, coronal and sagittal CT images. The point clouds of lung and trachea with tubular features (ring or semi-closed ring) are extracted respectively, and the point clouds in three directions are fused to obtain a complete three-dimensional model of lung and trachea, which effectively solves the problems of tracheal loss, leakage and the like, and can extract higher numbers of lung bronchi.

We present the following article in accordance with the MDAR reporting checklist (available at http://dx.doi.org/10.21037/atm-20-7300).

Methods

Expression of pulmonary trachea in CT images

Figure 1A is an axial CT image with 512×512 pixels. CT images contain information on the structure and distribution of pulmonary muscle tissue, bone, trachea, blood vessels, and air noise. In order to reduce the workload of data processing, only the CT image of the lung parenchyma needs to be preserved and the air noise in the pulmonary parenchyma needs to be removed, as shown in Figure 1B. Although the CT pixel values of the blood and the wall of the pulmonary vessels are close to each other, the CT value of the air in the lung and trachea is quite different from that of the trachea wall. Therefore, those features with an annular structure (ring, elliptical ring and long ring) can be determined as lung and trachea in CT images,, as shown in Figure 1C, while those with a solid structure (round, elliptical and long strip) can be determined as pulmonary blood vessels, as shown in Figure 1D, and the probability of a non-closed annular structure, especially a C-shaped structure, is relatively high. To extract the pulmonary trachea, it is necessary to distinguish it from other tissues and organs in the pulmonary parenchyma according to its different scale characteristics (“tubular features”) in CT images, to effectively extract the point cloud data of the pulmonary trachea.

Figure 1 Lung CT image, lung tissue, and lung vessel and airway characteristics. (A) Axial CT image; (B) lung tissue image; (C) lung vessel characteristics and (D) lung airways characteristics.

It should be specially noted that in the CT image in one orientation, the trachea presents as a ring and an elliptical ring, while in the CT image in the other two orientations, the trachea presents as a long ring or an elliptical ring. In the CT image in one orientation, it presents as a circular or oval solid structure, and in the CT image in the other two orientations, it presents as a solid strip or an oval solid structure, as shown in Figure 2.

Figure 2 Shapes of a lung vessel and airway in axial, sagittal, and coronal CT images. Lung airway in (A) axial CT image; (B) sagittal CT image and (C) coronal CT image. Lung vessel in (D) axial CT image, (E) sagittal CT image and (F) coronal CT image.

Extraction of pulmonary trachea and pulmonary vessels by connected domain

In order to avoid the loss of suspected pulmonary trachea and further reduce the workload of image processing, it is necessary to first extract all of the trachea and vessels in the pulmonary parenchyma. The connected region algorithm is used to extract the connected region of lung airways and blood vessels.

The connected region refers to adjacent image regions having the same pixel values (6). As shown in Figure 3, in the 8-neighborhood (A, B, C, D, F, G, H, I) of pixel point E, if the pixel values of points B, C, D and E are the same and adjacent, then the pixel points B, C, D, E are considered to be the same connected domain, and likewise the pixel point A and the pixel points F, G, H, I are considered to be different connected domains. The connected region algorithm calculates the connected region of each pixel point, and merges the connected regions of the pixel points to find and mark each connected region in the image. The pulmonary trachea and blood vessels extracted by the connected region algorithm are shown in Figure 4, where (A) is the blood vessel connected region, (B) is the suspected trachea connected region, and (C) is the trachea connected region.

Figure 3 Connected domain judgment and connected domain algorithm. (A) Connected domain graph and (B) connected domain judgment principle.

Figure 4 Pulmonary trachea and vessels extracted by the connected method. (A) Blood vessel connected region; (B) suspected trachea connected region; (C) trachea connected region.

Extraction of annular pulmonary trachea by dynamic tubular edge contour algorithm

In CT images, the annular structure is the pulmonary trachea, and the annular trachea is extracted by the edge contour algorithm. Edge contour extraction must first obtain contour points, assuming that white areas are the target pixel values and black is the background, as shown in Figure 5. Firstly, it is judged whether the pixel values of the target pixel 8-neighborhood are the same. If they are, the target pixel is an internal point; otherwise it is a contour point. As shown in Figure 5A, dots 1, 2, 3, 4, 5, 6, 7, 8, 9 are all internal dots, while corresponding dots 10, 11, 12, 13, 14, 20, 21, 22, 24, 25 are all contour dots (7). Secondly, the contour line is determined according to the contour line algorithm, as shown in Figure 5C, and three contour boundary lines are searched.

Figure 5 Edge contour extraction algorithm. (A) Original image, (B) outline point, and (C) contour boundary line.

Because an annular structure has two contour lines compared with other solid structures, it can be judged as an annular or non-annular structure according to the number of contour lines. If it is an annular structure, it can be regarded as trachea as shown in Figure 6A. Figure 6B and C are solid structures that need further determination.

Figure 6 Contours of pulmonary vessels and airways: (A) closed airway contour; (B) closed lung vessel contour; and (C) unclosed airway contour.

Determination of unclosed annular pulmonary trachea

As shown in Figure 1, an annular or elongated annular structure is determined as trachea, while a solid structure may be a blood vessel or unclosed trachea due to CT imaging quality problems. The geometric features of the unclosed trachea are obviously different from the structure of solid blood vessels. The curvature of bone line of the unclosed trachea is obviously smaller than that of the bone line of other solid structures. The unclosed trachea is extracted by using this feature.

Extraction of skeleton lines in connected domains of solid structures by template elimination method is as follows. Set the target pixel value as 1, the background point as 0, and mark the 8-neighborhood of point P as shown in Figure 7 (8). Eqs. [1]–[6] are the judgment condition. If conditions [1], [2], [3], and [4] are met, pixel points are marked as points deleted. Next, cycle the pixel points in each solid structure. If conditions [1], [2], [5] and [6] are met simultaneously, the pixel points are also marked as points deleted. At the end of a loop, all removable points are deleted. By repeatedly running these steps, the curve formed by the remaining pixel points is the bone line of the connected domain of the solid structure. Figure 8A shows the connected region of unclosed annular lung and trachea, Figure 8B shows the skeleton line of the connected region, and the end points of the skeletonization curve are connected, as shown in Figure 8C. Figure 8D is the connected domain of a pulmonary vascular solid structure, Figure 8E is the skeleton line of the connected domain, Figure 8F is the closed skeletonization curve.

Figure 7 Skeleton template.

Figure 8 Extract and complete unclosed pulmonary airway and vessel bone lines. (A) Unclosed pulmonary airway; (B) airway bone lines; (C) complete unclosed pulmonary airway; (D) pulmonary vessel; (E) pulmonary vessel bone lines; and (F) complete pulmonary vessel bone lines.

A(P) represents the number of “01” in the P1, P2, …P8 sequence, and B(P) is the number of “1” in the 8-domain window of P.

Ellipse fitting is used to ellipse the closed skeleton curve shown in Eq. [7] (9):

$$f\left(x,y\right)=x^2+Axy+B{y}^2+Cx+Dy+E=0$$[7]

$q\left(x_i,y_i\right)\left(i=1,2,\mathrm{..},N\right)$ Set to any point on the boundary of the closed curve, which is the sum of squares of errors of the fitting points, as shown in Eq. [8]:

$$F\left(A,B,C,D,E\right)={\underset{i=1}{\overset{N}{{\displaystyle \sum \left(f\left(x_i,y_i\right)\right)}}}}^2$$[8]

According to the extreme value theorem, when Eq. [9] is satisfied, the sum of error squares F of fitting points is the smallest, and coefficients A, B, C, D and E can be calculated by the Gaussian elimination method, thus fitting the elliptic equation of the bone curve.

$$\frac{\partial F}{\partial A}=\frac{\partial F}{\partial B}=\frac{\partial F}{\partial C}=\frac{\partial F}{\partial D}=\frac{\partial F}{\partial E}=0$$[9]

The elliptic long axis a and the elliptic short axis b fitting the ellipse are shown in Eqs. [10] and [11], and the ratio of the long axis to the short axis is t.

$$a={\left(\frac{2\left(ACD-B{C}^2-D^2+4BE-A^{2}E\right)}{\left(A^2-4B\right)\left(B-\sqrt{A^2+\left(1-B^2\right)+1}\right)}\right)}^{\frac{1}{2}}$$[10]

$$b={\left(\frac{2\left(ACD-B{C}^2-D^2+4BE-A^{2}E\right)}{\left(A^2-4B\right)\left(B+\sqrt{A^2+\left(1-B^2\right)+1}\right)}\right)}^{\frac{1}{2}}$$[11]

Finally, the pulmonary trachea is determined according to the proportional relationship between the long and short axes of the ratio. Because image quality results in an unclosed annular trachea, the ratio of the long and short axes of the fitting ellipse is small, as shown in Figure 9. For solid blood vessels, whether elliptical, circular or elongated, the ratio of the long and short axes fitting the ellipse is generally large. According to experimental statistical analysis, the relationship between the ratio value range and the extraction effect (the number of tracheal pixel points extracted) is shown in Figure 10, and t=20 (pixel) is more appropriate.

Figure 9 Elliptical fitting curve for (A) unclosed pulmonary airway bone line and (B) pulmonary vessel line.

Figure 10 Relationship between the ratio of long and short axes and number of extracted pixels.

Of course, in the very few cases of the difference between the long and short axes being greater than 20, it is still possible to be trachea. In this case, if ratio is greater than 20 in the axial CT images and less than 20 in the other coronal or sagittal CT images, this feature can be used to ensure that all of the lungs and trachea are not missed.

Statistical analysis

EXACT09 index of PICASSO was used to evaluate the advantages and disadvantages of the algorithm. Firstly, the ability of dynamic tubular edge contour algorithm to reduce leakage and maintain small branches was demonstrated by processing the small trachea part of pulmonary trachea. Meanwhile, the standard data set in EXACT09 was used to test its indicators. The selected data sets were Case 21, Case 22, Case 23, Case 24, Case 25, Case 26, Case 27, Case 28 and Case 29.

Results

The experimental materials were 680 CT images obtained by 360-degree scanning of the chest of two patients with a SOMATOM Definition AS+16 spiral CT scanner in the Department of Respiratory Medicine of West China Hospital. The pixel size of each image was 512×512, and the slice distance between images was 0.6 mm. This study passed the ethical review of West China Hospital and encrypted the patients’ identities. The flow chart of extracting the lung and tracheal tree is shown in Figure 11.

Figure 11 Dynamic edge extraction algorithm flow chart.

Figure 12 is a pulmonary tracheal tree extracted by the proposed method. The pulmonary tracheal level reaches 10–11, which is two more levels than the existing methods, thus better solving the problem of lung and trachea loss and leakage.

Figure 12 Pulmonary airway extraction image.

In order to further compare the effectiveness of our novel “dynamic tubular edge contour algorithm” in extracting pulmonary trachea, extraction experiments were carried out on the first nine CT data sets (Cases 21–29) in the “EXACT09 Airway Segmentation Challenge” (http://image.diku.dk/extract/) using this method (10,11). The extraction effect is shown in Figure 13. In the airway segmentation algorithm of the EXACT09 Airway Segmentation Challenge, 20 types of CT scan are evaluated, and the participating methods are evaluated through seven performance indicators. In order to compare fairly, we evaluated our method by the same protocol, and the performance parameters are shown in Table 1. Among them, Case 21 in the standard data set is compared. Figure 14 shows the effect of extracting pulmonary trachea using the dynamic tubular edge contour algorithm, it can be seen that the dynamic tubular edge contour algorithm extracts more bifurcation points and branches of the bronchus at the marked position, extract higher-level bronchial trees, and has richer extraction details. In particular, the blue circled area in the figure shows the identification and reconstruction of the small tracheacan.

Figure 13 Extraction of three-dimensional model of pulmonary airways from standard data sets.

Table 1 Quantitative evaluation of segmentation method
Full table

Figure 14 Comparison of pulmonary airway extraction effects in a standard CT library.

Table 1 shows that our method has good reproducibility. It achieved high sensitivity (detected tree length =81.4%) with a low false-positive rate (4.42%). Figure 14 shows the results for all eight cases. For a single chest CT scan, on a standard PC, our airway segmentation takes 1–2 minutes. Compared with the other two algorithms, the proposed algorithm has advantages in extracting high-level bronchus and can effectively solve the problems of tracheal loss and omission.

Discussion

The main algorithm of extracting high-level pulmonary airways using CT images of patients’ lungs is generally missing and leaking, and it is not possible to extract more advanced numbers of airways. Based on the characteristics of “tube structure” in CT image, an algorithm based on dynamic tube wall edge contour is proposed, a complete three-dimensional model of pulmonary trachea can be obtained, which effectively solves the problems of trachea loss and leakage, and can extract 10–11 pulmonary trachea in standard CT data set. In the dynamic tube wall edge contour extraction lung tracheal algorithm, the length axis ratio of the ellipse fitted to the pulmonary tracheal skeletal line will affect the recognition rate and misreading rate of the extraction lung tracheal tube. In the Figure 10, We find that when we set t as 20, we get the best result, but it is an empirical value, which changes in the face of different data sets, and its value has certain limitations. At the same time, we find that different CT resolutions have an important influence on the accuracy of extraction, and the resolution of CT images has a linear positive relationship with the accuracy of tracheal extraction. In addition, our results show that in the process of medical reconstruction, the organ tissue in CT image is described and identified according to the pre-shaped shape of organs and tissues, which is characterized by high accuracy compared to other reconstruction algorithms. These studies in this paper provide a new reconstruction strategy for organ and tissue reconstruction based on medical images.

Conclusions

Our novel dynamic tubular edge contour algorithm uses the patient’s chest CT image to extract the lung parenchyma according to threshold segmentation and morphological principles, converts the axial scanning CT image into coronal and sagittal “CT images” of the human body, and extracts the tracheal and blood vessel information in three directions (axial, sagittal and coronal). The annular pulmonary trachea is extracted by the edge contour algorithm, the skeleton line in the connected domain of the solid structure is extracted by the template elimination method, the closed skeleton line is elliptically fitted by ellipse fitting, and the unclosed annular pulmonary trachea is extracted according to the ratio of the long and short axes of the ellipse. The point cloud data of the pulmonary trachea in three directions are fused to obtain a complete three-dimensional model. The method can effectively extract pulmonary bronchioles to level 10–11, with a low false-positive rate (4.42%) and high sensitivity (detected tree length =81.4%). It can effectively solve the problems of tracheal loss and omission, and the extraction effect has certain advantages in bronchial level and integrity. During puncture operation, the puncture path can be selected according to the reconstructed pulmonary trachea to ensure that the large trachea is avoided during puncture, which has important theoretical and practical significance for reducing complications and the mortality rate.

Acknowledgments

Funding: Lung Puncture Path Planning and Puncture Tracking and Positioning Device Research (19ZDYF1402); Research on Path Planning and Puncture Tracking Location Methods of Pulmonary Puncture (2018SCUH0041).

Footnote

Reporting Checklist: The authors have completed the MDAR reporting checklist. Available at http://dx.doi.org/10.21037/atm-20-7300

Data Sharing Statement: Available at http://dx.doi.org/10.21037/atm-20-7300

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at http://dx.doi.org/10.21037/atm-20-7300). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved.

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Bray F, Ferlay J, Soerjomataram I, et al. Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin 2018;68:394-424. [Crossref] [PubMed]
2. Azfar Ali H, Lippmann M, Mundathaje U, et al. Spontaneous hemothorax: a comprehensive review. Chest 2008;134:1056-65. [Crossref] [PubMed]
3. Shalhub S, Neptune E, Sanchez DE, et al. Spontaneous pneumothorax and hemothorax frequently precede the arterial and intestinal complications of vascular Ehlers-Danlos syndrome. Am J Med Genet A 2019;179:797-802. [Crossref] [PubMed]
4. Liu C, Zhao R, Pang M. A fully automatic segmentation algorithm for CT lung images based on random forest. Med Phys 2020;47:518-29. [Crossref] [PubMed]
5. Graham MW, Gibbs JD, Cornish DC, et al. Robust 3-D airway tree segmentation for image-guided peripheral bronchoscopy. IEEE Trans Med Imaging 2010;29:982-97. [Crossref] [PubMed]
6. Gharasuie MM, Gaffari A. editors. An efficient run-based method for connected component labeling. 2015 9th Iranian Conference on Machine Vision and Image Processing (MVIP); 2015.
7. Singer HM, Bilgram JH. Extracting contours of crystals. Journal of Crystal Growth 2004;261:122-34. [Crossref]
8. Jinhe W, Zongying O, Xiaodong X. Discuss of Skeletonization and Recognition for Engineering Drawings. Computer Engineering 1999.
9. Teng Z, Kim JH, Kang DJ. Ellipse detection: a simple and precise method based on randomized Hough transform. Optical Engineering 2012;51:7203. [Crossref]
10. Ukil S, Sonka M, Reinhardt JM. Automatic segmentation of pulmonary fissures in X-ray CT images using anatomic guidance. Proceedings of SPIE - The International Society for Optical Engineering 2006;6144:213-23.
11. Fabijańska A. Segmentation of pulmonary vascular tree from 3D CT thorax scans. Biocybernetics and Biomedical Engineering 2015. [Crossref]

(English Language Editor: K. Brown)