CT Scanners Based on Airway Modeling and Seed Prediction

Proceedings of the 7th IFAC Symposium on Modelling and Control in Biomedical Systems Aalborg, Denmark, August 12-14, 2009

Airway Segmentation for Low-Contrast CT Images from Combined PET/CT Scanners Based on Airway Modeling and Seed Prediction Chaoqun Fang*, Xiu Ying Wang*,**, David Dagan Feng*,**,*** *

Biomedical and Multimedia Information Technology (BMIT) Research Group School of Information Technologies, J12 The University of Sydney, NSW 2006, Australia ** Center for Multimedia Signal Processing Department of Electronic and Information Engineering Hong Kong Polytechnic University, Hong Kong

***

Med-X Research Institute, Shanghai Jiao Tong University, Shanghai, China

Abstract: The combination of positron emission tomography (PET) and computed tomography (CT) scanning provides a superior access to both functional information and anatomical structures of the airway tree. However, airway tree segmentation from such low-dose and low-contrast CT images is a challenging task due to the limitation of the image resolutions. Complex anatomical structure of airway tree and partial volume effect pose other difficulties in airway segmentation. Conventional airway segmentation algorithms often produce less than satisfying results. In this paper, we propose a novel method for fully automatic airway tree segmentation for CT images from combined PET/CT scanners. In our method, airway modeling is used in seed extraction and prediction, and a new strategy is devised for identifying potential airway branches that are not detectable by conventional 3D region growing. In comparison with traditional 3D region growing segmentation algorithm, our method outperforms with not only retrieving considerably larger number of branches, but also providing more accurate geometric information. Keywords: biomedical system, image segmentation, combined PET/CT scanner, airway segmentation, region growing

1. INTRODUCTION Airways are important parts of the respiratory system, which are structured as a tree with branches that become narrower, shorter, and more numerous as they penetrate deeper into the lungs (Aykac et al., 2003). There are plenty of diseases related to airway and its surrounding areas such as bronchitis, lung cancer, emphysema and tuberculosis. With remarkable progress in the use of CT scanners in the medical field, an essential component of these new devices is the design of sensitive and reliable methods for assessing alterations in regional airway structure and function. Segmentation of pulmonary CT images occupies a significant step for most of these airway image analysis applications. Traditional ways to estimate regional volumes in the thorax include manually tracing boundaries, which is both laborious and time-consuming. Airway segmentation also plays an important role in airway construction. Once the airway tree is extracted, quantitative analysis can be performed to evaluate the tree structure and function. Based on the constructed image, many physiological and pathological conditions, such as stenosis and tumors, can be detected and analyzed (Aykac et al., 2003). Airway tree segmentation in CT image is a challenging problem because of its complex anatomic structure. Intensity

978-3-902661-49-4/09/$20.00 © 2009 IFAC

homogeneity between airway bronchus and lung tissues makes the airway segmentation even more difficult. Partial volume effect (PVE) poses a significant challenge in identifying small bronchus at the lower layer of the airway tree (Tschirren et al., 2005). If the small bronchus axis is parallel to the slice, the density values in the airway lumen will change and make the branch difficult to detect based on the variance of the intensity values (Mori et al., 1996). Moreover, while large airways can be well distinguished from the airway walls and adjacent vessel, small airways often have gray levels similar to that of the surrounding lung tissues. This often leads to either under- or over-segmentation errors (or namely “leaks”). A variety of airway segmentation techniques have been developed for segmenting the 3D airway tree. These methods can be categorized into four classes (Sluimer et al., 2006): (1) knowledge-based techniques (Sonka et al., 1996, Park et al., 1998); (2) region growing based methods (Mori et al., 1996, Mayer et al., 2004, Singh et al., 2004); (3) central-axis analysis based methods (Tschirren, 2002, Swift et al., 2002); (4) mathematical morphology based methods (Kiraly et al., 2002, Aykac et al., 2003, Francoise et al., 1999). In contrast, the theory of fuzzy connectivity was used by J.Tschirren (Tschirren et al., 2005). In his method, based on the comparison between intensity values of the input image and the values of two seeds (trachea, airway wall), fuzzy

198

10.3182/20090812-3-DK-2006.0051

7th IFAC MCBMS (MCBMS'09) Aalborg, Denmark, August 12-14, 2009

membership values (affinity), were assigned to the voxels. It allowed two regions to compete against each other and decide the classes of voxels. Also recently, Aykac (Aykac et al., 2003) used 2D gray-scale reconstruction of the airway (grayscale closing with increasing kernel size), and then applied iterative closing and dilation to get all the potential airways. This process was followed by three dimensional reconstructions that contained forward pass and backward pass to extract the final airway tree.

(EC) is used to determine the threshold. An initial airway volume is extracted at this stage. 2.3 Iterative 2D Operation

However, the above algorithms were not specially designed for low-dose and low-contrast CT images and hence few of them have conclusively proved adequate for low-contrast 3D CT images due to limitation of the image resolutions. In the case of low-dose scans, existing segmentation algorithms may fail by stopping early or leaking into the lung regions. Furthermore, reported algorithms usually require manual adjustment of parameters to achieve an optimal result. Developing an automatic and robust airway segmentation method to assist physicians to access regional abnormalities is thus in high demand. 2. METHOD 2.1 Overview Three dimensional region growing, taking advantage of speed and simplicity, is widely used as a method for airway segmentation from High Resolution CT images (HRCT) (Rizi et al., 2008). However, as mentioned above, this method frequently misses airway regions in low-contrast CT images. To improve segmentation accuracy, in this paper, a new prediction strategy is proposed to identify airway regions that are not detectable by conventional 3D region growing. Our method consists of two main stages (Fig. 1): automatic 3D region growing and iterative 2D operation. In the first stage, automatic 3D region growing is performed to acquire an initial airway segmentation result. This method has the drawback of under-segmentation, due to PVE and the homogeneity of the intensity value between airways and its surrounding pulmonary lung. In the second stage, iterative 2D airway segmentation is initiated slice-by-slice to identify potential airway locations in the 2D cross-sectional images. This stage is primarily divided into three steps. Firstly, according to the topological modeling of bronchus region, seeds are extracted from the previously segmented regions; secondly, based on the seeds selected, potential seeds on the following slice are predicted. Thirdly, 2D region growing process is applied on each slice to extract the detailed airway branches. This 2D operation will be repeated until there are no new regions being identified. 2.2 Automatic 3D Region Growing In this stage, 3D painting region growing proposed by Mori (Mori et al., 1996) is applied to extract the bronchus area by gradually increasing a region growing threshold until the segmented region exploded into the lung. Explosion Control

Fig. 1. Flow chart of the automatic airway tree segmentation. To obtain a more complete airway tree, new seeds are predicted on the basis of the extracted regions from 3D region growing. This 2D operation step is to identify potential airway locations in the 2D cross section images. Because we need to extract and predict airway tree branches in both upper and lower parts of the lung tissue, bi-directional (forward and backward) detections are processed sequentially. (1) Seed localization from the extracted regions Each extracted region may contain one or multiple seeds (as shown in Fig. 2(a) and Fig. 2(b), each seed represents one bronchus region). Locating all the seeds from these regions is important not only for completing airway branches but also for predicting new seeds at the following stage. Seeds from identified regions are extracted based on geographic models of bronchi regions (Fig. 2(c)). We start by searching the local minimum points, which are called valleys, in each separated region. Among these valleys, the one with the lowest intensity value is selected as the first seed. To find out the other seeds, the intensity values of all pixels along the route between this seed and the testing valley are assessed. If these two points do not fall into one bronchus region, pixels along the route will have much higher gray level values than the existing seed. Utilizing this feature, the other seeds are selected based on the scheme as below. We

define (1)

(2)

ℜ as a set of elements along route

ψ = c ,c ,..., c ( n )

between the testing valley and the

existing seed. As Fig. 2(c) presents, point c(i) with the highest

199

7th IFAC MCBMS (MCBMS'09) Aalborg, Denmark, August 12-14, 2009

intensity value along route ψ is identified as boundary point cB (1). The valley with boundary point cB which satisfies the criterion described as (2) is selected as another seed, where the Threshold value is predefined (i.e. 150). If two or more seeds have been located, new valley points need to be

In this paper, a new airway model is designed for the procedure of seed prediction. As shown in Fig. 3, while shapes of cylinders and circles represent bronchi and seeds of the airway tree respectively, branch direction is formed by the last two seeds in the same branch. Based on this airway

valley valley

Seed

Seed

New seed

Seed

New seed

(a)

(b)

valley

(a)

valley

cB

Fig. 3. The procedure of candidate seed prediction. Cylinders represent the imaginary airway branches. Gray ellipses represent extracted regions. White circles represent localized seeds and black circles represent newly predicted seeds. (a) Newly found seed is discarded as it belongs to one bronchus region with existing seed. (b) Newly found seed is kept.

cB c(1) c(2) ….

(n)

c(n-1)c

ψ valley

seed

ψ seed seed

(c)

Get all the valley points from regions in previous slice For each branch in the old tree If the branch hasn’t been terminated From the valleys, find the one that is closest to the predicted point If the valley point doesn’t belong to any existing region Add the valley to new branch Add the branch to the tree Else if the branch has no child Terminate the branch Add the branch to new tree End if End if End for

seed (d)

Fig. 2. Seed extraction from valleys. (a) Airway region with one seed.(b) Airway region with two seeds. (c) The procedure of assessing boundary point between the seed and the valley point. (d) Valleys need to be tested with each existing seed. examined to ensure they do not fall into the same bronchus region with any existing seeds (Fig. 2(d)). Otherwise, these valley points are removed. To trace the 2D operation of airway tree, we employ list structures to simulate the growing of the tree, with each list representing a simulated branch. Every simulated branch is comprised by the seeds selected from sequential slices, indicating the growing path of the tree. As there is only one seed for each of the first two slices, we add both seeds to the same branch. c B = max ( Intensity ( c ( i ) )) (i )

(1)

Intensity ( c B ) − Intensity ( ExistSeed ) > Threshold

(2)

c

(b)

∈ℜ

(2) Potential seed prediction For airway prediction, the utilization of adaptive region of interest (ROI) was proposed by Kitasaka et al. (Kitasada et al., 2003) and Tschirren et al. (Tschirren et al., 2005). A similar prediction strategy was also used by Aykac et al. (Aykac et al., 2003) when identifying the locations of the trachea and carina.

Fig. 4. Pseudo code of seed prediction model, the procedure of seed prediction is summarized as follows: Candidate region, which restricts the area of valley candidates, is the projection of the selected region at the prior slice. Predicted point is chosen as the linear extension of the last two seeds in one simulated branch. Finally, seed candidate is the valley candidate that is closest to the predicted point in candidate regions. Since some branches combine as one bifurcation point, newly predicted seeds may fall into one bronchus region with other existing seeds. In such a case, this newly predicted seed is discarded (Fig. 3(a)), otherwise, the whole branch after attaching newly seed is added to the new tree (Fig. 3(b)). The pseudo code of the seed prediction strategy is given in Fig. 4.

200

(3) Seed validation

7th IFAC MCBMS (MCBMS'09) Aalborg, Denmark, August 12-14, 2009

In order to filter falsely detected seeds, 2D region growing procedure is applied to each seed candidate while local EC is employed to control the growth area, after which a set of rules are defined to ascertain the final seeds. The airways, which are surrounded by the airway walls, exhibit high contrast borders (Sonka et al., 1996). For the points that belong to the airway, the process of region growing would end with an optimal threshold value. In contrast, because there are no airway walls encircling those points inhabiting inside lungs, the regions would explode dramatically and be deleted after validation.

3.1 Sensitivity and negative predictive value (NPV) Reliability of the segmentation algorithm was statistically measured by the sensitivity and the negative predictive value (NPV), which are defined as follows (10-11).

Table 1. Sensitivity and NPV Sensitivity comparison of five clinical cases (%) (proposed method/3D region growing)

Gen 1

To evaluate the airway areas, protocols (3-9) are defined in this paper to monitor all extracted airway pixels by verifying the size and seed number of each bronchus candidate. These protocols also serve as the control strategy to prevent the problem of large amounts of leaks.

0 1

The thresholds for large, middle and small regions are defined in three ranges: [0, n2],,[ n2, n1], [n1, N1] respectively . Since airway branches of the first two airway tree generations usually have larger regions and more accurate than smaller branches, where partial volume effect makes a significant influence on the result, we define different local EC value (35) and various policies as validation rules (6-7). SizeSon represents the size of bronchus region where the tested seed inhabits, whereas SizePA represents that of the its parent point. By examining the relationship of SizeSon and SizePA, we allow airway seeds to grow sufficiently while prohibit leaks. In addition, for those regions with only one seed inside, there is a possibility that it may bifurcate to two branches. As a result, the evaluating parameters are determined by the number of seeds in the region of its parent seed (8-9). At the end of this procedure, all the seeds failing to pass the rules are deleted from the simulated branch and that tree would be terminated. Otherwise, the size of bronchus region where testing seed stays would be recorded. SizePA > n1 ⇒ EC = N 1

n 2 < SizePA ≤ n1 ⇒ EC = α × SizePA

SizePA ≤ n2 ⇒ EC = N 2 SizePA > n1 ⇒ SizeSon < β * SizePA SizePA ≤ n1 ⇒ SizeSon < δ × SizePA NumOfSeeds= 1 ⇒ β = 2, δ = 3

NumOfSeeds≠ 1 ⇒ β = 1.4, δ = 3

2

3

4

5

All subjects

100.0/100. 100.0/100.0 0 100.0/100. 100.0/100.0 100.0/100.0 100.0/100.0 100.0/100.0 100.0/100.0 0 100.0/100.0 100.0/100.0 100.0/100.0 100.0/100.0

2

100.0/40.0 100.0/25.0 100.0/60.0

100.0/60

100.0/60.0 100.0/54.2

3

87.5/25.0

71.4/14.3

75.0/25.0

81.9/45.5

87.5/25.0

80.4/28.3

4

100.0/75.0

50.0/0.0

50.0/30.0

75.0/50.0

85.7/14.3

66.7/24.2

5

50.0/0.0

null

100.0/0.0

100.0/0.0

100.0/0.0

85.0/0.0

6

null

null

100.0/0.0

100.0/0.0

null

100.0/0.0

70.0/25.0

73.7/31.6

87.9/48.5

93.1/31.0

82.4/34.5

26.3

15.2

19.4

28.9

18.8

Total 86.4/31.8 NPV (%)

5

Sensitivity = AutoD / ManuD

(10)

NPV = FalseD /(FalseD+AutoD)

(11)

The values of AutoD, ManuD and FalseD in the above equations are defined as below: AutoD: the number of automatic correctly detected airway branches ManuD: the number of manually detected airway branches FalseD: the number of automatic falsely detected airway branches The branch segmentation starts from generation zero (main trachea) and increases by one as it is bifurcated. For the cases that one branch splits into three bronchi, the two bronchi with similar diameter are classified as sister bronchi and are assigned as the next generation of the third bronchus.

(3) (4) (5) (6) (7) (8) (9)

To evaluate the results, our proposed method is compared with the conventional 3D region growing. The comparison results of the airway branch detection procedures in all the generations are summarized in Table. 1.

where n = 100mm 2 , n = 60mm 2 , N = 600mm 2 , N = 60mm 2 , α = 3 1 2 1 2

3. EXPERIMENTAL RESULTS Our automated airway tree detection method was evaluated by using low-contrast volumetric CT images of the human lungs from a combined PET/CT scanner. Five normal subjects were scanned at full inspiration (Royal Prince Alfred Hospital, Sydney, Australia).

Our newly proposed segmentation method identified all the bronchi branches up to the second generation (main trachea represents generation zero). The sensitivity was 80.4% for the third generation and 66.7% for the fourth generation; while 3D region growing method had the sensitivity of 28.3% and 24.2% respectively. Beyond the fifth generation, the 3D region growing failed to detect any new branches, while our proposed method was able to identify the smaller bronchi. Among all five clinical cases, there were a total of 142 branches detected manually and 117 branches detected by our

201

7th IFAC MCBMS (MCBMS'09) Aalborg, Denmark, August 12-14, 2009

approach, indicating an overall detection sensitivity of approximate 82.4%, compared with 34.5% by 3D region growing. The false detecting percentage in the five data sets was an overall of 18.8%. It is noticeable that with the growth of the sensitivity of the airway detection by our method, the false detecting percentage increased at the same time. In the worst case, the false detecting number reached 28.9% of the manually detecting numbers. 3.2 Visual comparison of Reconstructed Trees (a)

(b)

(c)

(d)

As is shown in Fig. 5, besides identifying more airway tree regions, our method provides more accurate geographic information than conventional 3D region growing method at cross section CT images. Fig. 6 shows a surface-rendered display of airway trees reconstructed from the segmented image datasets. All the left images were outputs of 3D region growing method while the right images were the result of our proposed method. As shown in Fig. 6(a), the right branch of the airway tree is completely missed and the left branch is rather thin and incomplete. In contrast, it is noticeable that our method extracts more branches than traditional 3D region growing (Fig. 6(b)). For the case with better image quality, while 3D region growing detects more detailed branches (Fig.6(c)), our proposed method improves the number of correctly detected branches for the reconstructed tree (Fig. 6(d)). Compared with 3D region growing method, our airway segmentation method not only extracts more airway regions, but also reconstructs a more robust airway tree. The utilization of local threshold, which varied among different parts of the tree, is the main reason for this advantage of proposed method over 3D region growing.

incomplete

missing

(a)

Fig. 6. Comparison of the visualized airway tree structures between 3D region growing process and our method. 4. CONCLUSION We have developed a novel automatic airway segmentation method for 3D low-contrast CT thoracic images from combined PET/CT scanners. The proposed method takes advantage of priori anatomical knowledge of pulmonary airway structure to guide airway modeling and introduces a new prediction scheme to extract the airway tree. Experiments on clinic cases demonstrate that this proposed method is able to segment more complete airway trees for the low-contrast CT images, although minor leaks may still happen in the cases when small bronchi were extracted, due to the limit of image resolutions.

(b)

ACKNOWLEDGEMENT

missing

This research was supported by ARC and PolyU grants. REFERENCES incomplete

incomplete

(c)

(d)

Fig. 5. Comparisons of 3D region growing method and our proposed method. (a)(c) Fused image of chest volume and extracted airways after being applied by 3D region growing. (b)(d) Fused image of chest volume and extracted airways after being applied by proposed method.

Aykac, D., Hoffman, E. A., Mclennan, G. & Reinhardt, J. M. (2003) Segmentation and analysis of the human airway tree from three-dimensional X-ray CT images. IEEE Transactions on Medical Imaging, volume (22), 940-950. Francoise, P., Catalin, I. F., Andre, C. & Philippe, G. (1999). Modeling, segmentation, and caliber estimation of bronchi in high resolution computerized tomography. Journal of Electronic Imaging, volume (8), 36-45.

202

7th IFAC MCBMS (MCBMS'09) Aalborg, Denmark, August 12-14, 2009

Kiraly, A. P., Higgins, W. E., Mclennan, G., Hoffman, E. A. & Reinhardt, J. M. (2002). Three-dimensional human airway segmentation methods for clinical virtual bronchoscopy. Academic Radiology, volume (9), 1153-1168. Kitasaka, T., Mori, K., Hasegawa, J.-I., Suenaga, Y. & Toriwaki, J.-I. (2003). Extraction of bronchus regions from 3D chest X-ray CT images by using structural features of bronchus. In: Pros. of 17th Int. Congr. Computer Assisted Radiology and Surgery (CARS), volume (1256), 240-245 Mayer, D., Bartz, D., Fischer, J., Ley, S., Del R O, A., Thust, S., Kauczor, H.-U. & Heussel, C. P. (2004). Hybrid segmentation and virtual bronchoscopy based on CT Images. Academic Radiology, volume (11), 551-565. Mori, K., Hasegawa, J., Toriwaki, J., Anno, H. & Katada, K. (1996). Recognition of bronchus in threedimensional X-ray CT images with applications to virtualized bronchoscopy system. In: Proc. of 13th Int. Conf of Pattern Recognition, volume (3), 528– 532. Park, W., Hoffman, E. A. & Sonka, M. (1998). Segmentation of intrathoracic airway trees: a fuzzy logic approach. IEEE Transactions on Medical Imaging, volume (17), 489-497. Rizi, F. Y., Ahmadian, A., Sahba, N., Tavakoli, V., Alirezaie, J., Fatemizadeh, E. & Rezaie, N. (2008). A hybrid fuzzy based algorithm for 3D human airway segmentation. In: Proc. of 2nd Int. Conf on Bioinformatics and Biomedical Engineering, 22952298. Singh, H., Crawford, M., Curtin, J. & Zwiggelaar, R. (2004). Automated 3D segmentation of the lung airway Tree using gain-based region growing approach. In: Proc. of 7th Int. Conf on Medical Image Computing and Computer-Assisted Intervention, volume (3217), 975-982 Sluimer, I., Schilham, A., Prokop, M. & Van Ginneken, B. (2006). Computer analysis of computed tomography scans of the lung: a survey. IEEE Transactions on Medical Imaging, volume (25), 385-405. Sonka, M., Wonkyu, P. & Hoffman, E. A. (1996). Rule-based detection of intrathoracic airway trees. IEEE Transactions on Medical Imaging, volume (15), 314326. Swift, R. D., Kiraly, A. P., Sherbondy, A. J., Austin, A. L., Hoffman, E. A., Mclennan, G. & Higgins, W. E. (2002). Automatic axis generation for virtual bronchoscopic assessment of major airway obstructions. Computerized Medical Imaging and Graphics, volume (26), 103-118. Tschirren, J. H. (2002). Segmentation, Skeletonization, and Branchpoint Matching – A Fully Automated Quantitative Evaluation of Human Intrathoracic Airway Trees. In: Proc of 5th Int. Conf of Medical Image Computing and Computer Assisted Intervention, volume (2489), 12-19 Tschirren, J., Hoffman, E. A., Mclennan, G. & Sonka, M. (2005). Intrathoracic airway trees: segmentation and airway morphology analysis from low-dose CT

203

scans. IEEE Transactions on Medical Imaging, volume (24), 1529-1539.