Automatic pulmonary artery-vein separation in CT images using a twin-pipe network and topology reconstruction

Introduction

Lung cancer is a significant global health concern with high morbidity and mortality rates (Siegel et al., 2022). In 2020, an estimated 19.3 million new cancer cases and nearly 10 million cancer deaths are expected worldwide, and lung cancer remains the leading cause of cancer deaths, accounting for around 18% of all cancer deaths. Computed tomography (CT), with its high resolution and ability to show the fine structures and density differences of lung tissues, has become a major tool for thoracic surgeons to screen, diagnose and treat diseases. Currently, in the clinic, surgeons generally manually reconstruct lung anatomy by using commercial medical software such as IQQA (Xu et al., 2019) and Mimics, which can then be used to analyze lung disease and the anatomical relationships among pulmonary arteries, veins, airways (Guo et al., 2022), and nodules (Zheng et al., 2021). However, manual reconstruction of arteries and veins is time-consuming and visually demanding due to the large quantity of CT data, which puts a heavy burden on surgeons. Moreover, the accurate separation of lung veins and arteries plays a crucial role in diagnosing and evaluating lung diseases such as lung cancer. Therefore, automatic and rapid vessel segmentation and separation of arteries and veins through CT scans is crucial for clinical diagnosis and can provide better survival prediction for lung cancer patients (Haq et al., 2022).

At present, the automatic separation of pulmonary arteries and veins from CT scans has become a popular topic in research. However, this task exist several difficulties, which can be attributed to the following reasons. Firstly, arteries and veins are indistinguishable due to their similar intensity values in noncontrast CT images. Secondly, the vessel tree structure is extremely complex and dense, with arteries and veins close to each other and intertwined. And finally, artifacts, partial volume effects, and patient-specific vessel tree structural abnormalities cause difficulty in A/V separation.

Recently, deep learning-based methods have demonstrated powerful feature extraction and analysis capabilities, and have achieved remarkable results in medical image analysis (Abedalla et al., 2021; Nguyen et al., 2021). With the development of deep learning, lung vessel segmentation has been widely explored (Nam et al., 2021). Although these studies have had a certain theoretical basis and achieved satisfactory results, there are still some intractable problems in the results after the segmentation of vessels. For example, the pulmonary vessel is a typical tree-like structure, with thick and thin vessels connected and numerous bifurcations. Existing deep learning methods rarely consider the topological information of the vessel structure, which can result in the separated arteries and veins being discontinuous or with breaks in the vessel branches. Additionally, tubular vessels in the lungs are small and have thin walls, making them difficult to detect and classify accurately using deep learning methods. These vessels are often mistaken for noise or artifacts in the image, leading to misclassification or missed detection. Furthermore, due to the similarity in size, shape, and appearance of arteries and veins, they may be difficult to distinguish in medical images and can be confused with each other after the A/V separation.

To address these concerns, we propose a pulmonary A/V separation method using topology information and a twin-pipe network. Our contributions are summarized as follows:

• 1) To address the discontinuity problem in A/V trees after separation, a vessel tree topology is constructed by combining scale-space particles and multi-stencils fast marching (MSFM) to ensure the continuity and authenticity of the topology reconstruction.

• 2) To address the poor separation effect of tubule vessels, a twin-pipe network is proposed to learn the multiscale differences between arteries and veins, as well as the characteristics of small arteries that are closely accompanied by bronchi, given that arteries and bronchi tend to accompany each other.

• 3) To address the intertwining of arteries and veins, a topology optimizer is designed that considers both interbranch and intrabranch topological relationships to optimize the A/V classification results.

Our method is evaluated on our private database and on the public dataset CARVE14 (Charbonnier et al., 2015). Furthermore, we evaluated our method on another device and a CT scan of an arterial enhancement model to demonstrate its generalizability across different devices.

The rest of this article is organized as follows. In “Methods” shows the proposed approach to solve the A/V separation problem in noncontrast CT images. “Experiments” mainly introduces the experiment methods, including data sources, experimental setup details, and evaluation metrics. Then, in “Results”, the experimental results, the ablation experiment, and the generalization experiment are presented. “Discussion” elaborates on the discussion. Finally, “Conclusion” concludes this article.

Related works

Although pulmonary A/V separation is a difficult problem, many scholars have proposed methods in the past decade. Studies have shown that the vasculature within the lungs is highly variable, but some inherent anatomical properties are usually present. One such property is that arteries in the lungs typically follow the bronchial tree, while veins tend to run in the interstitial spaces between the branches. As shown in Fig. 1, the arteries with accompanying bronchi are not evident when the vessels are near the hilum of the lung. As the arteries move away from the hilum, the bronchi begin to follow the arteries closely. Some methods rely on bronchial features for A/V separation. For instance, Tozaki et al. (2001) used information about distances between vessel segments and bronchi to separate arteries and veins, and Nakamura et al. (2005) classified pulmonary arteries and veins based on the distance from the bronchi region to the vessel segment and the distance from the nearest interlobar to the vessel. Similarly, Bülow et al. (2005) based upon the fact that the pulmonary artery tree accompanies the bronchial tree, designed a method of “arterialness” to classify each vessel segment. These pulmonary A/V separation methods rely on the quality of airway segmentation, while the use of bronchial features to aid A/V separation is rarely applied in deep learning. Our method takes this into account.

Another type of A/V separation study focuses on the anatomical priors of the vessels themselves, using information on their connectivity to separate arteries and veins. For example, Saha et al. (2010) proposed a method to separate A/V of pulmonary using morphological features of vessels. The approach involves selecting seed points and utilizing fuzzy distance transformation for tracking to ensure vessel connectivity. Finally, arteries and veins points are separated through multiscale iterative growth. Wala et al. (2011) designed an automated trace-based separation scheme that tracked arterial seed points from the basal pulmonary artery region and detected bifurcations to separate the artery from other nearby isometric structures. Park et al. (2013) proposed a method to find voxel-based trees with minimal construction energy inside vessel segmentation and divided them into a group of subtrees that share the same A/V classification. Additionally, Kitamura et al. (2016) designed a method to classify pulmonary artery and vein voxels based on vascular connection information and energy minimization of high-order potentials. However, these methods are usually semi-automatic or dependent on specific CT quality. Besides, the use of the anatomical priors for the pulmonary A/V separation method is susceptible to limitations of the anatomical structures themselves, such as bronchial leakage, vessel adhesions and discontinuity.

Moreover, several new ideas for automatic A/V separation have been proposed. Charbonnier et al. (2015) proposed a method to convert vessel segmentation into a geometric graph representation. Payer et al. (2016) proposed an automated algorithm for A/V separation based on two integer programming. Through integer programming, parameter tuning was carried out to extract subtrees, and A/V separation was performed. The method used local information to detect the attached arteries and veins in the geometric graph to generate subgraphs. Then, the volume difference between arteries and veins was used for classification. Park, Bajaj & Gladish (2018) developed a method that utilizes sphere inflation tracking from the endpoint to the heart. This approach aims to differentiate between bifurcation and crossover points and automatically detect the pulmonary trunk. Jimenez-Carretero et al. (2019) proposed a new pulmonary A/V separation scheme, in which a specialized graph-cut method was designed to ensure the connectivity and consistency of the vessel subtrees. Yu et al. (2021) proposed a combined algorithm for pulmonary A/V separation based on a subtree relationship to separate the adhesion points of the vessel tree structure. Generally, these methods are effective on the image with noise-free or uniform intensities but require complex parameter adjustment and parameter sensitive.

In recent years, many researchers have attempted to solve the problem of A/V separation using deep learning techniques, which offer better robustness and stability. Nardelli et al. (2018) proposed a convolutional neural network for A/V classification of pulmonary lobar vessel particles, an initial A/V classification is performed for each particle in the extracted vessel tree using CNN. The final classification result was obtained by a graph-cut strategy that combined connectivity and pre-classification information. Zhai et al. (2019) designed a network linking the CNN and graph convolutional network (GCN) to combine local image information with graph connection information to train the constructed graph, and the A/V separation results were obtained. Although U-Net is a common method in medical image segmentation, few works have been done in recent years for A/V separation. Qin et al. (2021) proposed learning tube-sensitive U-Net for lung A/V classification and anatomy prior of lung context map and distance transform map is combined to detect more arterioles and veins. Heitz et al. (2021) relies on three-paths 2.5D U-Net networks along axial, coronal, and sagittal slices to extract the tubular structures of the lung. These networks have partly solved the problems in pulmonary A/V separation, however, none of them consider the complex structures and topological connectivity of the lung vessels. Confusion and discontinuity of the pulmonary arteries and veins remained, especially in the hilar, adhesions, and terminal vessel regions.

Methods

The overall framework of the pulmonary A/V separation method in this article is shown in Fig. 2, including vessel tree topology extraction module, twin-pipe network, and topology optimizer. Portions of this text were previously published as part of a preprint (https://arxiv.org/abs/2103.11736).

In the vessel tree topology extraction module, the vessel tree is segmented by 3D U-Net from the chest CT images, and the topology is extracted from the vessel tree. Then, the distance transform is used to guide and compensate for the missing points. In the twin-pipe network, 3D patches are taken from each particle on the topological tree as the centers. The vessel segment and the tubule vessel are then trained separately to derive the preliminary classification results of the vessel particles. In the final optimizer module, the topology of the subtrees and their branches are extracted from the topological tree, the branch confidence is calculated, and subtrees are pruned to optimize A/V classification using topological connectivity. Finally, the arteries and veins of pulmonary topology are reconstructed, and the reconstructed A/V and separately segmented arteries and veins near the hilum are fused to obtain the final result.

Vessel tree topology module

The proposed method begins with vessel tree segmentation. Vessel trees are extracted by using 3D U-Net with the patch of [256, 256, 16]. Then, a complete topological vessel tree is constructed by topological extraction of the vessel tree. Due to the applicability of the skeleton algorithm and the characteristics of the vessel tree itself, it is often not available to obtain accurate results of skeleton structure. To obtain a better topological veesel tree, we flexibly combine scale-space particles (Estépar et al., 2012) and the multi-stencils fast marching (MSFM) (Liu et al., 2018), superior to some existing special methods such as thinning, geometry, shortest path, etc., the specific steps are as follows.

Firstly, based on the characteristics of the vessel tubular structure, scale-space particles reconstruction is used to extract the initial skeleton of the vessel and obtain the local information of the vessel mask. The scale-space particle sampling method exploits the theory of linear scale-space to localize the features of the image described by the Hessian. The vessel tree is represented as a set of particles, each of which contains vessel scale, orientation, and intensity information. Therefore, particles are represented as X = { $x_i$}, where $x_i$ represents a particle point. However, after vessel tree reconstruction, the vessels are discontinuous due to the inability of this method to identify non-tubular structures, such as the bifurcation of vessels, resulting in local point loss. Additionally, there is no parent-child relationship between the particles.

Secondly, the skeleton extraction algorithm based on distance transform can maintain the connectivity of vessel tree. The MSFM is used to extract information about global connectivity from the vessel mask. The potential terminal points and the root node are obtained by calculating the distance map. Then, the vessel trajectory is obtained by iteratively tracing from the terminal points to the root node. During this process, the branch online confidence score is calculated in the time map to determine whether the trace iteration should be updated or stopped. In this article, we take the 3D vessel mask as input, aim to output the vessel skeleton tree G, and assign a 3D spatial coordinate and radius to each vessel skeleton tree G node. The degree of each vessel tree node can be between 1 and 3. However, the radius of each node of the obtained vessel tree is only an approximation and does not really reflect the vessel tubular shape. There are some differences between the reconstructed results and the real vessel tree.

Finally, in order to ensure the continuity and authenticity of the topology reconstruction, we flexibly combined the global and local information of the vessel skeleton. In this article, a particle-based 26-neighborhood search method is adopted to extract the topological structure of the vessels by guiding the time map to make up for the lost vessel particles. vessel particles can be classified into three categories: terminal points (or false-positive terminal points), branching points, and bifurcating points. We indicate the number of points in the 26 neighborhoods of vessel particles $x_i$ as ${\mathrm{\Omega }}_{26}$( $x_i$). The values of ${\mathrm{\Omega }}_{26}(x_i)$ is different based on the type of particle $x_i$. For the terminal points, ${\mathrm{\Omega }}_{26}(x_i)=1$. Branching points are assigned ${\mathrm{\Omega }}_{26}(x_i)=2$. Bifurcating points have ${\mathrm{\Omega }}_{26}(x_i)>2$. Further discrimination is needed due to the presence of branch point loss or branch fracture resulting in false positive terminal points. Travel in the direction of the terminal vessel points, and if the terminal point is still on the vessel mask after traveling, it is considered a false positive terminal point, where the travel distance is between one and two scales. Otherwise, it is considered the true terminal point. For the false positive terminal points, the trajectory of the lost vessel points is obtained by MSFM in the time map, in which the time map is calculated from the 3D distance map. Finally, the complete topology tree is obtained. Based on the extracted topology tree, the information from each particle is used to construct the graph structure we need. Through the category of particle points and the parent-child node relationship between particles, we redefined a graph T = $\{X,\epsilon \}$ composed of nodes X = { $x_i$} and edges $\epsilon$ = { ${\epsilon }_{i_j}$}. Node are defined as ${\epsilon }_{i_j}$ = 1, when $\mathrm{\Omega }$( $x_i$) = 2. The process of the vessel tree topology extraction method is shown in Fig. 3. And compared with the scale-space particles, the advantages of the vessel tree topology extraction in this article are shown in Fig. 4.

Twin-pipe network

The twin-pipe network is designed to improve the preliminary classification accuracy by learning the differences in A/V characteristics caused by different scales. One pipeline network is trained on the full vessel graph to learn image information and global connectivity, while the other pipeline network is specifically trained on the tubule vessel graph. Taking into account the physiological feature of small arteries accompanying the small trachea, the twin-pipe network incorporates the CT original images and the vessel enhancement images as input patches. This enables the network to learn additional distinguishing features of pulmonary arteries and veins. Finally, the preliminary classification results are obtained through the mutual correction module.

The network structure is shown in Fig. 5. Twin-pipe network is made up of two non-local CNN-GCN classifiers. Traditional CNN networks have limitations in capturing long-range dependencies that extract the global understanding of visual scenes. In this article, we consider both image information and connectivity information by connecting Non-local CNN and graph convolutional network (GCN). The Non-local CNN (Sun et al., 2020) considers the global context information of the image, while the GCN module (Selvan et al., 2020) can learn the local and graph connection information. Combining these two modules together is useful for analyzing the vessel tree. In order to connect the Non-local CNN with GCN, this article adapts the end-to-end method of CNN connecting GCN proposed by Zhai et al. (2019). The definition of the GCN layer is as follows:

(1) $$H^{(l+1)}=\sigma (W{H}^{(l)}{\mathrm{\Theta }}^{(l)})$$

The input of layer (L + 1) is the output of layer (L), the input size of $H^{(l+1)}$ is $N\times F^{(l+1)}$, and the size of output $H^{(l)}$ is $N\times F^{(l)}$, where W is the weight matrix of nodes X. ${\mathrm{\Theta }}^{(l)}$ is a hierarchical training parameter matrix, $\sigma (\cdot )$ is an activation function, such as ReLU function. Image feature matrix Y, the output of Non-local CNN network, is taken as the input of GCN layer.

(2) $$H^{(l+1)}=\sigma (WY{\mathrm{\Theta }}^{(l)}),whereY=\mathrm{\Phi }(P|\theta )$$

P is the local image patches centered around each particle, oriented perpendicularly to the blood vessel. Each image patch with the size of S = [32, 32, 3] is labeled as artery or vein according to the center voxel, i.e., $y_i$ $\in$ {1, 2}. $\mathrm{\Phi }()$ is a general Non-local CNN layer. $\theta$ is the Non-local CNN network training parameter.

To train the twin-pipe network, a set of n nodes and their corresponding image patches are randomly selected from graph T. These selected nodes and patches are represented as N with a size of n * S. The image patches of their neighbors are also required as we use GCN networks. Neighbor patches are denoted by MN which has the size m $\times$ n $\times$ S, where m is the number of neighbors. N and the neighborhood MN are the input of the network. We can extract each branch subgraph from graph T, and each branch subgraph is either an artery or a vein. The tubule vessel graph is the set of tubule branches in the full vessel graph, we perform vessel enhancement and normalization on the original image to enhance the differences between arterioles and small trachea. The extraction process is shown in Fig. 6. The twin-pipe network is trained and predicts whether each central voxel is an artery or a vein.

Experiments

Results

The evaluation structure of the method in the article is as follows. First, this article presents the results of our method on the Siemens-Non-16 dataset, as well as the results of the A/V separation methods in recent years. Second, the ablation experiments are performed on the Siemens-Non-16 dataset to validate each component of the proposed method. In addition, the generalization experiments are carried out to verify the proposed method performance from the GE-Non-8 dataset, the GE-A-8 dataset, and the CARVE14 dataset.

Discussion

The proposed end-to-end framework constructs topological trees by extracting vessels from CT and obtains preliminary classification results through a twin-pipe network. Subsequently, a topological optimizer is used to refine the A/V classification results and reconstructs pulmonary A/V through the scale information at last. An objective inter-work comparison could be difficult due to the difference in medical data sets and implementation, so our experiments are mainly conducted on our own CT data. The experimental results in Table 1 indicate that the average accuracy of the proposed method can reach 96.2% compared with that of manual reference annotation. The experimental results in Table 4 show that our method is also applicable to CT images from different devices and different modes, thus the proposed method has good generalization.

Pulmonary A/V separation is a challenging problem in medical image analysis. The separation ought to overcome the complexity of the pulmonary structure as well as the relatively limited resolution of CT images. Our evaluation method for A/V separation is mainly based on topological particles rather than voxel classification. The reason is that in clinical practice, surgeons focus more on the structural branching direction of vessels, consistent with the evaluation system based on topological particles. A/V separation provides effective information for surgical planning and navigation, while proximal vascular extension aids surgeons in locating vessels more quickly. Therefore, the final A/V separation results include vessels close to the hilum of the lung, as shown in Fig. 7. Our evaluation system for the A/V separation method does not include vessels near the hilum, mainly because the vessels near the hilum are abnormally large and non-tubular, and the vessel topology could not be extracted.

First, we extract vessel topology for the training of the classification network. Currently, some topology methods such as thinning (Palágyi & Németh, 2017), geodesics and shortest paths (Chen, Mirebeau & Cohen, 2016), and some other special methods (Liu et al., 2017) are used to construct vessel topology. However, validating these methods with accuracy is difficult, mainly because (1) it is difficult for these methods to obtain manual annotations based on high-resolution images, and (2) these methods pay more attention to vessel branch direction, bifurcation pixel redundancy, loopback, fracture, etc., We propose a vessel tree topology method that fully uses the advantages to solve the topological fracture problem caused by particle loss. As shown in Fig. 8, topological skeleton nodes are basically located on the central axis of vessels, and there is no obvious pixel redundancy or fracture at the bifurcation of vessels.

Second, we design a twin-pipe network for the preliminary classification of A/V and to improve the classification accuracy of tubule vessels. As shown in Table 2, the classification accuracy of tubule vessels in the baseline network is lower than that of the full vessels, while the sensitivity of the twin-pipe network in the tubule vessels is higher than that of the baseline network. The experimental results show that the error frequency is lower in the tubule vessels in the twin-pipe network. That is, the twin-pipe network can effectively improve the tubule vessel classification accuracy. The experiment proves that: (1) the Non-local CNN-GCN network can effectively combine local information, global image information, and graph connection information. Differences in A/V characteristics can be automatically learned, eliminating the need for complex parameter tuning. (2) The twin-pipe network effectively captures the differences in A/V characteristics at various scales and learns the close association between the tubule bronchus and the artery. It is superior to the baseline Non-local CNN-GCN network and does not depend on the airway segmentation results.

Finally, the topology optimizer extracts the topology of the subtrees and their branches’ refinement results using the method in “Topology Optimizer”. Then, we use topology subtrees and topology branches for postprocessing. As shown in Table 3, the proposed topology optimizer accuracy is superior to that of the subtree-based and branch-based topology optimization. Fig. 9 shows the reconstructed results of different topology strategy optimizers. The results in the first row show that when the number of points on the branch is small, the branch-based topology optimization method is prone to prediction errors. This is owing to the branch-based refinement strategy focusing on the relationships within the branches and ignoring the topological relationships between branches. The second row shows that topology optimization based on the subtree strategy is prone to prediction errors in the case of A/V intersection. Due to the complex intertwining of arteries and veins trees, the points where arteries and veins intersect are often mistaken for arterial subtree bifurcation points. As a result, venous branches are misclassified as arterial sub-trees, leading to classification errors.

Figure 10 shows the final separation result of our method at the intersections of the arteries and veins. The proposed topology optimizer considers the interbranch and intrabranch relationships. The extracted topological subtree is used to maintain the spatial connectivity of topological particles. The branch confidence calculated using topological branches is used to correct the topological subtree. To some extent, our method can solve the misclassification problem caused by interlaced arteries and veins. This method made some achievements in A/V separation, but there are still some unavoidable shortcomings. The method is based on topological particle reconstruction of A/V classification, and some errors are found between the reconstructed vessel size and the actual situation. Although it is possible to classify the A/V at the intersections, the scale of the intersections after reconstruction is abnormal and often results in fractures or abnormal expansion. This condition occurs because topological particles at non-tubular intersections of arteries and veins cannot be accurately extracted during topology construction, thereby leading to topological particle scale mutation or deletion. To apply it to clinical practice, we can further improve the classification results and performance by optimizing the topology extraction algorithm and detecting the intersection situation, as well as increase the clinical applicability through semi-automated graphical user interface interaction operations.

Conclusion

In this article, we propose a novel method for automatically separating pulmonary arteries and veins. This method is applicable for noncontrast CT scans that lack vessel edge information, and by incorporating vessel topology information, it can resolve adhesions of arteries and veins and misclassification of small vessels. We extract the vessel skeleton through the vessel tree topology method, obtain the preliminary classification results after twin-pipe network training, and use the topology structure information for postprocessing. Our approach has been tested on the public dataset CARVE14 and our private dataset. Experimental results show that the proposed method has good performance and generalizability to CT images from different devices and different modes. The method proposed in this article will better assist physicians in the planning of surgery for lung disease.

Supplemental Information