# Pulmonary Lobe Segmentation in CT Images Based on Lung Anatomy Knowledge

Authors: Yuanyuan Peng; Hualan Zhong; Zheng Xu; Hongbin Tu; Xiong Li; Lan Peng
Journal: Mathematical Problems in Engineering
 (2021)

---

## Abstract

In computed tomography (CT) images, pulmonary lobe segmentation is an arduous task due to its complex structures. To remedy the problem, we introduce a new framework based on lung anatomy knowledge for lung lobe segmentation. Firstly, the priori knowledge of lung anatomy is used to identify the fissure region of interest. Then, an oriented derivative of stick filter is applied to isolate plate-like structures from clutters for lobar fissure verification. Finally, a surface fitting model is employed to complete the incomplete fissure surface for lung lobe segmentation. Compared with manually segmented fissure references, the designed approach obtained a high medianF1-score of 0.8865 in the left lung and obtained a high median F1-score of 0.9200 in the right lung. The average percentages of the segmented lung lobes in the lung lobe ground truth are 0.960, 0.989, 0.973, 0.920, and 0.985 for the left upper, left lower, right upper, right middle, and right lower lobes, respectively. The perfect performance of the proposed scheme is tested by visual inspection and quantitative evaluation.


---

## Body

## 1. Introduction

In medical practice, recognition of lung lobes is useful for clinical diagnosis and lung disease assessment [1]. Hence, accurate segmentation of lung lobes is urgently needed. Anatomically, the left lung is composed of two lobes, whereas the right lung is composed of three lobes [2]. Clinically, lung lobe segmentation is extremely important in surgical treatment of lung diseases [3, 4]. However, manual segmentation of lung lobes is impracticable in CT images due to the large number of images. Therefore, automatic segmentation methods for lung lobes are necessary.A schematic diagram is shown in Figure1; the human lungs are composed of five lobes. The right lung is subdivided into the right upper lobe (RUL), the right middle lobe (RML), and the right lower lobe (RLL), separated by the right oblique fissure and the right horizontal fissure, whereas the left lung is subdivided into the left upper lobe (LUL) and the left lower lobe (LLL), separated by the left oblique fissure. Lung lobes are served by separate bronchial trees and their corresponding pulmonary arteries and are anatomically independent. Unfortunately, lung lobe segmentation is an arduous task due to complex lung anatomy.Figure 1
Schematic diagram of the human lungs.Pulmonary lobe segmentation algorithms can be divided into three categories. The first approaches consist of methods that search for the fissure surface based on shape features and gray-level information. The second approaches consist of methods that use this information from airways, vessels, and fissures to estimate lobar boundaries. The third approaches consist of methods that construct pulmonary atlases to guide lung lobe segmentation. Although many methods have been presented, no standard method has yet emerged [5].Considering the fact that fissures are physical boundaries between different lobes, Pu et al. [6] combined a threshold with a marching cube algorithm to identify fissures and utilized the direction difference for fissure decomposition; then, an implicit surface fitting model was applied for pulmonary lobe segmentation. Shamonin et al. [7] designed a fissure enhancement filter to distinguish between fissures and noise; then, an iterative B-spline fitting model was applied to segment pulmonary lobes. However, their algorithms [6, 7] may not be able to accurately segment pulmonary lobes in some cases due to their fissure decomposition error and confounding structures. For the same purpose, Ross et al. [8] employed a particle system and ridge surface sampling to identify lobar fissures; then, a thin plate spline was used to segment pulmonary lobes. The approach fails if the particle signal is extremely weak. Using a different strategy, Lo et al. [9] used a fissure confidence function to verify pulmonary fissures; then, a surface evolution algorithm was applied for pulmonary lobe segmentation. Lassen-Schmidt et al. [10] firstly selected a number of representative lobar boundary points; then, a plane was constructed by the given points to segment pulmonary lobes, but their algorithms [9, 10] need manual interaction. Recently, Wei and Hu [11] used a texture analysis method based on the neural network and a stack projection to verify lobar fissures; then, a fissure surface model was used for pulmonary lobe segmentation. Qi et al. [12] used a line enhancement filter followed by a uniform cost search to verify pulmonary fissures; then, an implicit surface model was used to segment pulmonary lobes, but their methods [11, 12] may not be able to accurately verify pulmonary fissures. The reason is that the orientation of fissures can vary too much in some cases [7].Motivated by the assumption that the information of airways and blood vessels may contribute to segment pulmonary lobes, Ukil and Reinhardt [13] used airways and vascular trees to construct a region of interest (ROI) encompassing fissures and applied a ridgeness measure to enhance fissures; then, a fissure decomposition technique and a minimum curvature biharmonic spline interpolation algorithm were employed to segment pulmonary lobes. In a similar way, Doel et al. [14] applied a similar approach to identify fissure ROI and developed a fissureness filter to highlight the fissure representation; then, a multilevel B-spline method was used for lung lobe segmentation. Their algorithms [13, 14] fail if pulmonary veins cross lobar fissures. Using a different strategy, Lassen et al. [15] constructed a cost image based on fissures, vessels, and bronchi for the watershed transformation to segment pulmonary lobes. Its improved version was presented by Bragman et al. [16] using a groupwise fissure prior approach for pulmonary lobe segmentation. Later, Bauer et al. [17] presented a subject-specific prior method for lung lobe separation. Recently, Giuliani et al. [18] proposed a lung lobe segmentation method using alpha-expansion. However, their methods [15–18] heavily depend on airway segmentation.In addition, a number of studies constructed pulmonary atlases to guide lung lobe segmentation. Following this strategy, Zhang et al. [19] combined the ridgeness measure and an atlas to verify pulmonary fissures; then, a fissure surface was constructed for pulmonary lobe segmentation. Later, van Rikxoort et al. [20] improved the method and presented a multiatlas approach to segment pulmonary lobes. Recently, Schmidt-Richberg et al. [21] combined the level sets with atlases to segment pulmonary lobes, but their methods [19–21] have a longer computation time.The proposed scheme aims at automatically segmenting pulmonary lobes. In order to identify fissure ROI, many studies [13, 14, 22, 23] used bronchial trees and pulmonary vessels to approximately construct lobar shapes. Their algorithms [13, 14, 22, 23] fail in some cases because of pulmonary veins may cross fissures. Different from the algorithms [13, 14, 22, 23] which have been previously proposed, we present a novel and effective framework to automatically identify pulmonary fissure ROI. For subsequent segmentation, we improved and extended our previous works [24–26] to accurately verify pulmonary fissures. Finally, a surface fitting model [6, 11] is used to complete the incomplete fissure surface for lung lobe segmentation. Experimental results show that the proposed method has a good performance in pulmonary fissure segmentation and lung lobe segmentation.


## 2. Materials and Methods

### 2.1. Data and References

Fifteen CT images were used to evaluate the validation of our scheme and acquired from the LObe and Lung Analysis 2011 (LOLA11) challenge [27], including various pathologic abnormalities. To better evaluate the segmented right fissures and left fissures, we divided CT images into left lungs and right lungs.A manual segmentation of lung lobes in 9 coronal slices was made by two human observers. Unfortunately, the ground truth is not publicly available. To make full use of the publicly available data for experimental evaluation, we manually segmented 15 lung lobes and regarded them as the lung lobe ground truth. In addition, lobar fissure segmentation is different from lung lobe segmentation [28]. Therefore, two internationally renowned experts verified the LOLA11 reference [25] and removed parts of the invisible reference and added parts of the visible reference that were lost in the LOLA11 reference. In this study, we regarded the fissure reference [25] as the fissure ground truth for experimental evaluation.


### 2.2. Overview of the Proposed Method

Our scheme (developed in MATLAB and C++) uses a new framework for lung lobe segmentation. The designed algorithm employs three stages: (1) fissure ROI identification, (2) lobar fissure verification, and (3) lung lobe segmentation. Figure2 depicts the flowchart of this computerized scheme. The following context presents these three stages.Figure 2
The flowchart of lung lobe segmentation.


### 2.3. Fissure ROI Identification

In this paper, we present a new computerized scheme for fissure ROI identification. Firstly, an effective method based on the graphic processing unit (GPU) is used to extract the centerline of airway trees [29], and a graph search approach is employed to obtain the subtrees in airway trees. Then, a shape prior and graph cut method [30] is applied to separate pulmonary arteries and veins. Furthermore, the priori knowledge of airways and pulmonary arteries is fused to segment pulmonary arteries. Finally, an alpha-shape algorithm [31, 32] is applied to identify fissure ROI.


#### 2.3.1. Airway Segmentation

To accurately segment airways, a useful and valuable method [29] based on GPU is applied for tube extraction. Firstly, Gaussian smoothing with the parameter σ = 0.5 is used to blur the dataset I. Subsequently, a gradient vector field V⟶ is designed:(1)V⟶=−∇I∗σ.Normalize(2)V⟶nv⟶=V⟶v⟶V⟶v⟶ifV⟶nv⟶≥Fmax,V⟶v⟶Fmax,else.Here,Fmax controls the normalization. To prevent different tubular structures diffusing together, a gradient vector flow (GVF) method [33] is used:(3)EV⟶=∫μ∇V⟶v⟶2+V⟶0v⟶V⟶v⟶−V⟶0v⟶2dv⟶2,where V⟶0 and μ denote the initial gradient vector field and the weighting constant. Xu and Prince [33] proposed an iterative approach to solve the Euler equation:(4)μ∇2V⟶−V⟶−V⟶0V⟶02=0.The purpose of the GVF model is to make the following tube detection filter (TDF) [34, 35] invariant to contrast of airways [29]. Subsequently, the TDF is applied to detect tubular structures. Firstly, three eigenvectors of the Hessian matrix corresponding to the tube are estimated: one of three eigenvectors e⟶1 is related to the direction of the tube, and the other two eigenvectors e⟶2 and e⟶3 are associated with the cross-sectional plane of the tube. The main idea is to create a circle with a radius r in the cross-sectional plane, and N is the number of spaced points i on the circle. Therefore, the direction from the center to each point i is(5)d⟶i=e⟶2sin2πiN+e⟶3cos2πiN.The average dot product is(6)Tv⟶,r=1N∑i=0N−1v⟶+rd⟶i⋅−d⟶i.Whenr is increased, the average dot product is calculated again. In this paper, r ranges from 0.5 to 25. Lastly, the airway centerline is extracted by traversing a ridge in the TDF results. All values in TDF results above a threshold are extracted and treated as candidate starting points, starting from the candidate point with the positive tangent direction t0 and negative −t0. For a current centerline point xi and the corresponding tangent direction ti, all local neighbors xin are considered:(7)xixin⟶⋅ti>0.The neighbor point with the highest TDF response is chosen as the next centerline pointxi+1. Finally, the airway centerline is acquired and marked with S.Different subtrees in the airway tree belong to different lobes [15]. Based on this theory, a graph search method is used to search for major bifurcations for bronchial tree separation in appropriate orientations in Figure 3. In the first step, a centerline extraction method [29] is used to transfer the airway tree segmentation to its skeleton representation. Then, a directed graph is modeled from the airway tree with the trachea as the root. To extract the bifurcations, a simple but effective criterion [37, 38] is used:(8)N26p∩S∪p≥4,where p and N26 denote the pixel in the centerline of the airway tree S and its 26-neighborhood area, respectively. If there are 4 or more points within a 3 ∗ 3 ∗ 3 cube (p and its 26 neighborhood), we consider the pixel p as a bifurcation. Finally, a graph search method is used to remove major bifurcations labeled with yellow dots in Figure 3(b) and calculate for each subtree in Figure 3(c) along with the particular separation direction n⟶ in Figure 3(b), in which EndRMB and LMB denote the ends in the right main bronchus and the left main bronchus, respectively.Figure 3
Subtrees in an airway tree. (a) An airway tree. (b) Graphical representation of a typical skeletonized airway tree (figure taken from [36]). (c) Skeleton of an airway tree; subtrees belonging to different lobes are displayed in different colors.
(a)(b)(c)


#### 2.3.2. Pulmonary Artery and Vein Separation

Airways are tubular structures and less adhered with each other. Furthermore, the GVF model has the advantage that it is feature preserving. Therefore, a method based on the GVF model is used to segment the airway tree, but the approach is not applicable for pulmonary artery and vein separation. As shown in Figure4, pulmonary arteries may adhere with pulmonary veins, and pulmonary veins may cross lobar fissures (marked with yellow arrows). To remedy the problem, a shape prior and graph cut method [30] is applied to separate and segment multiple interwoven vessel trees.Figure 4
Pulmonary arteries and veins. (a) Pulmonary arteries always adhered with pulmonary veins, and pulmonary veins may cross lobar fissures. (b) The local image corresponding to (a). As shown in the region labeled with yellow arrows, pulmonary arteries and veins are joined together. Furthermore, pulmonary veins may cross pulmonary fissures.Pulmonary artery and vein separation is a challenging task due to its complex structure. To solve the problem, the shape prior and graph cut method [30] is applied to segment the interwoven three-dimensional (3D) tubular trees. The first step is to identify all tubular objects by TDF. Three eigenvectors of the Hessian matrix corresponding to tubular objects are estimated: one of the three eigenvectors e⟶1 is related to the direction of the tube, and the other two eigenvectors e⟶2 and e⟶3 are associated with the cross-sectional plane of the tubular objects. The main idea is to create a circle with a radius r in the cross-sectional plane, and N is the number of spaced points i on the circle. Therefore, the direction from the center to each point i is(9)d⟶i=e⟶2sin2πiN+e⟶3cos2πiN.To segment the interwoven 3D tubular trees, a medialness method [30] is applied to separate pulmonary arteries and veins. The initial boundariness contribution and averaging boundariness contribution are(10)bi=σbV0v⟶+rd⟶id⟶i,b¯σb=1N∑i=0N−1σbV0v⟶+rd⟶id⟶i,where N represents the number of spaced points on the circle, r is the circle radius, di is the tube direction, and the gradient vector field V0 is with σb=rη0≤η≤1. The variance can be calculated as follows:(11)s2v⟶,r=1N∑i=1N−1bi−b¯2.Thus, the response of the final offset medialness is(12)R+v⟶,r=1N∑i=0N−1bi1−s2v⟶,rb¯2.To suppress the response away from the tube center, Bauer et al. [30] proposed a method based on an adaptive threshold:(13)R−v⟶,r=σbV0v⟶.The multiscale medialness response is(14)Rmultiv⟶=maxrmin≤r≤rmaxmaxR+v⟶,r−R−v⟶,r,0.Here,rmin=1 and rmax=5. After acquiring the multiscale medialness response, we set a threshold r=3 and take the medialness response larger than the threshold as seed points. Finally, a tree construction algorithm is used to separate pulmonary arteries and veins. The main idea is to extract the centerline of the multiscale medialness response by traversing a ridge map.Pulmonary arteries may adhere with pulmonary veins in Figure5(a). To separate the multiple interwoven vessel trees, we use the medialness method to extract the centerline of vessel trees, and the results are shown in Figure 5(b). As observed in the regions marked with yellow arrows in Figure 5, pulmonary arteries and veins have been separated.Figure 5
Pulmonary artery and vein separation. (a) Pulmonary arteries adhered with pulmonary veins, and the adhered regions are marked with yellow arrows. (b) Pulmonary arteries and veins are separated.
(a)(b)


#### 2.3.3. Fissure ROI Identification

Based on the reality that pulmonary arteries are close to bronchial trees [39–41], we combine the bronchial trees (BT) with vessel subtrees (VST) to remove mounts of pulmonary veins. We dilate BT and use a connected-component analysis method (CCAM) to select 5 biggest objects for pulmonary artery segmentation. To give a mathematical expression,(15)PA=max1≤i≤5CCAMVST+DBT∩VST,where D denotes the morphological dilation.To illustrate the validation of the proposed scheme, we integrate with bronchial trees, pulmonary arteries, and lobar fissures to describe the distribution map from different angles of view. Bronchial trees and lobar fissures are labeled with yellow color and green color in Figure6, respectively. It can be seen that pulmonary arteries are close to bronchial trees [40, 41], and pulmonary arteries cannot cross lobar fissures [15, 39].Figure 6
The distribution map among bronchial trees, pulmonary arteries, and lobar fissures is displayed in different angles of view, in which bronchial trees, pulmonary arteries, and lobar fissures are marked with different colors.
(a)(b)(c)As a result, we consider different bronchial trees and the corresponding pulmonary arteries belonging to different pulmonary lobes as a whole in Figure7(a). Then, an alpha-shape algorithm [31, 32] is used to construct lobe shapes. The alpha-shape algorithm was presented by Edelsbrunner and Mücke [31] that constructs an entire family of shapes for a given set of points. We adopt for the implementation of alpha-shape algorithm to construct the geometric shape of different bronchial trees and their corresponding pulmonary arteries. Finally, surface fitting processing and region filling processing are used to approximately estimate the lobe distribution map for fissure ROI identification in Figure 7(b).Figure 7
Lobe shapes. (a) Bronchial trees and the corresponding pulmonary arteries belonging to different lobes are displayed in different colors. (b) Lobe shapes corresponding to (a). (c) The relationships among lobe shapes, fissures, and bronchial trees. (d, e) The relationships between lobe shapes and fissures.
(a)(b)(c)(d)(e)Fissure ROI identification is useful in lung lobe segmentation. As shown in Figure8(a), the estimated right upper, middle, and lower lobe distribution maps are marked with ST1, ST2, and ST3, respectively. We just regard ST2 and ST3 as an entirety and ST1 as another part; then, the alpha-shape algorithm is applied to identify the right oblique fissure ROI, and its corresponding results are shown in Figure 8(b).Figure 8
Lobe shapes. (a) Right fissure ROI, in which the estimated right upper, middle, and lower lobe distribution maps are marked withST1, ST2, and ST3, respectively. (b) Right oblique fissure ROI.
(a)(b)


### 2.4. Lobar Fissure Verification

To effectively verify pulmonary fissures, an oriented derivative of stick (ODoS) filter [26] is applied to enhance pulmonary fissures. Then, an orientation partition scheme is used to separate pulmonary fissures and adhering clutters. Next, a connected-component analysis approach is adopted to remove noise for plate-like structure extraction. Finally, the identified fissure ROI is employed to isolate lobar fissures from adhering clutters.


#### 2.4.1. Retrospection of the DoS Filter and ODoS Filter

To efficiently detect the weak and thin pulmonary fissures, Xiao et al. [25] presented a derivative of stick (DoS) filter for pulmonary fissure enhancement. The filtering kernel consists of three parallel sticks (LS, MS, and RS), and θ is the stick orientation. Two 2D line strength measures lmaxθ and lminθ are(16)lmaxθx=max1L∑j=1LIjM−1L∑j=1LIjL,1L∑j=1LIjM−1L∑j=1LIjR−k⋅EIjM2−EIjM2,lminθx=min1L∑j=1LIjM−1L∑j=1LIjL,1L∑j=1LIjM−1L∑j=1LIjR−k⋅EIjM2−EIjM2,where IjL, IjM, and IjR denote the intensity of the jth pixel along the left, middle, and right stick, k=0.7 represents a positive coefficient, and E is an expected value operator. To suppress noise and other structures, Xiao et al. [25] integrated the line strength measures lmaxθ and lminθ in multidirections:(17)Fmaxx=maxmax0<θ≤πlmaxθ,0,Fminx=maxmax0<θ≤πlminθ,0.The DoS filtering response can be formulated as(18)Fο=FmaxοFmin.ο represents the cascading operator. To integrate the DoS filtering response in multidirections, Xiao et al. defined a shape-tuned fissureness function:(19)F3Dx=FοA+FοS+FοC⋅medianFοA,FοS,FοCmaxFοA,FοS,FοC.Here,FοA, FοS, and FοC denote the DoS filtering response Fο across the axial, sagittal, and coronal sections. However, the original DoS filter cannot efficiently suppress adhering clutters. To solve the problem, Peng and Xiao [26] improved and extended it by constructing its orientation response. The orientation response is(20)θmax=argmax0<θ≤πlmaxθiflmaxθ>0,NaN,else.A vector representation is(21)V⟶max=cosθmax,sinθmaxiflmaxθ>0,0,else.Therefore, the response of the ODoS filter can be described as follows:(22)V⟶F=F3D⋅V⟶max.


#### 2.4.2. Plate-Like Structure Extraction

To suppress noise and other structures, Peng and Xiao [26] firstly used a threshold T to modulate the ODoS filtering response:(23)V⟶T=V⟶FifV⟶F>T,0,else,where the threshold T we set is equal to 1. Based on the fact that pulmonary fissures appear as plate-like structures and have similar directions, Peng and Xiao [26] presented an orientation partition scheme to separate lobar fissures and noise in Figure 8. They divided the enhanced 3D CT image into n overlapped subregions, and the subregions are marked with R1,R2,…Rn−1 and Rn. To give a mathematical expression,(24)FRi=1ifθmax∈RiandV⟶T≠0,0,else.In 3D partitioned imagesFRi, clutters such as airways and vessels usually appear as small objects. As a result, they used a connected-component analysis method to select the candidate fissures, and the results are displayed in Figure 9.Figure 9
Plate-like structure extraction in the right lung. After using an orientation partition scheme, the enhanced images are divided into a number of planar structures.


#### 2.4.3. Oblique Fissure and Horizontal Fissure Verification

From the above processing, some plate-like structures are easily acquired, including lobar fissures and accessory fissures. As is shown in Figure10, accessory fissures have bad influence on fissure segmentation, but few methods have been proposed by using the priori knowledge of pulmonary structures to remove them. In this study, the identified fissure ROI is used to limit the search area and remove those unrelated accessory fissures for right fissure verification. Similarly, we can use the previously identified fissure ROI to verify right oblique fissures and left oblique fissures.Figure 10
Lobar fissure verification in the right lung, in which fissure ROI is used to limit the search area for fissure verification.After acquiring the right fissure (RF) and the right oblique fissure (ROF), we use a simple but effective criterion to verify the right horizontal fissure (RHF):(25)RHF=RF−ROF.As shown in Figure11, it can be seen that the proposed method has a good performance in lobar fissure verification.Figure 11
Right oblique fissure and horizontal fissure separation. (a) Right fissure (RF). (b) Right oblique fissure (ROF). (c) Right horizontal fissure (RHF).
(a)(b)(c)


### 2.5. Lung Lobe Segmentation

Lung lobe segmentation is an arduous task due to incomplete, deformed, and disrupted fissures. To remedy the problem, a surface fitting model [6, 11, 12] is applied to construct a complete surface. Subsequently, the complete surface is employed to segment lung lobes.


#### 2.5.1. Surface Fitting Model

Once oblique fissures and horizontal fissures are verified, we modeled a fissure surface in the form of(26)z=fx,y,where x, y, and z are the coordinate points in the 3D space. Imaging a thin and flexible plate for a given set of points by elastic bands is a nontrivial task. Each of these elastic bands draws the plate towards the given set of points, stretching only in the Z-direction. The given set of points was to fall on a simple planar surface [6, 11]:(27)z=w0+w1x1+w2y2+e.A matrix expression is(28)Z=AW+e.At every point of the surface, the surface fitting model attempts to force the partial derivatives of the surface in neighboring cells to be equal. A mathematical expression in the second set of linear equations is of the form(29)BW=0.To construct a smooth fissure surface, a Laplacian function is treated as the regularization term. It can be described as(30)B=α∂2fx,y∂2x∂2fx,y∂2y,where α is a scale factor that controls the effect of regularization. Then, we solve for the vector W by performing the following minimization:(31)Wmin=argminAW−Z2+BW2.Therefore, the surface fitting model is employed to construct a complete fissure surface, as shown in Figure12.Figure 12
Surface fitting model. (a) Right oblique fissure. (b) The surface fitting results corresponding to (a). (c) Right horizontal fissure. (d) The surface fitting results corresponding to (c).
(a)(b)(c)(d)


#### 2.5.2. Cutting Model

To segment pulmonary lobes, we present a cutting model to subdivide human lungs into five lobes. As shown in Figure13, the right complete oblique fissure surface is used to divide the right lung into the special part and the lower lobe, and then the right complete horizontal fissure surface is employed to separate the special part into the upper and middle lobes. In the same way, we use the left complete fissure surface to divide the left lung into left upper and lower lobes. As a result, the human lung is subdivided into five lobes.Figure 13
Pulmonary lobe separation. (a) Right complete oblique fissure surface divides the right lung into two parts. (b) Right complete horizontal fissure surface divides the special part into two parts. (c) The flowchart of the right lung lobe segmentation.
(a)(b)(c)


## 2.1. Data and References

Fifteen CT images were used to evaluate the validation of our scheme and acquired from the LObe and Lung Analysis 2011 (LOLA11) challenge [27], including various pathologic abnormalities. To better evaluate the segmented right fissures and left fissures, we divided CT images into left lungs and right lungs.A manual segmentation of lung lobes in 9 coronal slices was made by two human observers. Unfortunately, the ground truth is not publicly available. To make full use of the publicly available data for experimental evaluation, we manually segmented 15 lung lobes and regarded them as the lung lobe ground truth. In addition, lobar fissure segmentation is different from lung lobe segmentation [28]. Therefore, two internationally renowned experts verified the LOLA11 reference [25] and removed parts of the invisible reference and added parts of the visible reference that were lost in the LOLA11 reference. In this study, we regarded the fissure reference [25] as the fissure ground truth for experimental evaluation.


## 2.2. Overview of the Proposed Method

Our scheme (developed in MATLAB and C++) uses a new framework for lung lobe segmentation. The designed algorithm employs three stages: (1) fissure ROI identification, (2) lobar fissure verification, and (3) lung lobe segmentation. Figure2 depicts the flowchart of this computerized scheme. The following context presents these three stages.Figure 2
The flowchart of lung lobe segmentation.


## 2.3. Fissure ROI Identification

In this paper, we present a new computerized scheme for fissure ROI identification. Firstly, an effective method based on the graphic processing unit (GPU) is used to extract the centerline of airway trees [29], and a graph search approach is employed to obtain the subtrees in airway trees. Then, a shape prior and graph cut method [30] is applied to separate pulmonary arteries and veins. Furthermore, the priori knowledge of airways and pulmonary arteries is fused to segment pulmonary arteries. Finally, an alpha-shape algorithm [31, 32] is applied to identify fissure ROI.


### 2.3.1. Airway Segmentation

To accurately segment airways, a useful and valuable method [29] based on GPU is applied for tube extraction. Firstly, Gaussian smoothing with the parameter σ = 0.5 is used to blur the dataset I. Subsequently, a gradient vector field V⟶ is designed:(1)V⟶=−∇I∗σ.Normalize(2)V⟶nv⟶=V⟶v⟶V⟶v⟶ifV⟶nv⟶≥Fmax,V⟶v⟶Fmax,else.Here,Fmax controls the normalization. To prevent different tubular structures diffusing together, a gradient vector flow (GVF) method [33] is used:(3)EV⟶=∫μ∇V⟶v⟶2+V⟶0v⟶V⟶v⟶−V⟶0v⟶2dv⟶2,where V⟶0 and μ denote the initial gradient vector field and the weighting constant. Xu and Prince [33] proposed an iterative approach to solve the Euler equation:(4)μ∇2V⟶−V⟶−V⟶0V⟶02=0.The purpose of the GVF model is to make the following tube detection filter (TDF) [34, 35] invariant to contrast of airways [29]. Subsequently, the TDF is applied to detect tubular structures. Firstly, three eigenvectors of the Hessian matrix corresponding to the tube are estimated: one of three eigenvectors e⟶1 is related to the direction of the tube, and the other two eigenvectors e⟶2 and e⟶3 are associated with the cross-sectional plane of the tube. The main idea is to create a circle with a radius r in the cross-sectional plane, and N is the number of spaced points i on the circle. Therefore, the direction from the center to each point i is(5)d⟶i=e⟶2sin2πiN+e⟶3cos2πiN.The average dot product is(6)Tv⟶,r=1N∑i=0N−1v⟶+rd⟶i⋅−d⟶i.Whenr is increased, the average dot product is calculated again. In this paper, r ranges from 0.5 to 25. Lastly, the airway centerline is extracted by traversing a ridge in the TDF results. All values in TDF results above a threshold are extracted and treated as candidate starting points, starting from the candidate point with the positive tangent direction t0 and negative −t0. For a current centerline point xi and the corresponding tangent direction ti, all local neighbors xin are considered:(7)xixin⟶⋅ti>0.The neighbor point with the highest TDF response is chosen as the next centerline pointxi+1. Finally, the airway centerline is acquired and marked with S.Different subtrees in the airway tree belong to different lobes [15]. Based on this theory, a graph search method is used to search for major bifurcations for bronchial tree separation in appropriate orientations in Figure 3. In the first step, a centerline extraction method [29] is used to transfer the airway tree segmentation to its skeleton representation. Then, a directed graph is modeled from the airway tree with the trachea as the root. To extract the bifurcations, a simple but effective criterion [37, 38] is used:(8)N26p∩S∪p≥4,where p and N26 denote the pixel in the centerline of the airway tree S and its 26-neighborhood area, respectively. If there are 4 or more points within a 3 ∗ 3 ∗ 3 cube (p and its 26 neighborhood), we consider the pixel p as a bifurcation. Finally, a graph search method is used to remove major bifurcations labeled with yellow dots in Figure 3(b) and calculate for each subtree in Figure 3(c) along with the particular separation direction n⟶ in Figure 3(b), in which EndRMB and LMB denote the ends in the right main bronchus and the left main bronchus, respectively.Figure 3
Subtrees in an airway tree. (a) An airway tree. (b) Graphical representation of a typical skeletonized airway tree (figure taken from [36]). (c) Skeleton of an airway tree; subtrees belonging to different lobes are displayed in different colors.
(a)(b)(c)


### 2.3.2. Pulmonary Artery and Vein Separation

Airways are tubular structures and less adhered with each other. Furthermore, the GVF model has the advantage that it is feature preserving. Therefore, a method based on the GVF model is used to segment the airway tree, but the approach is not applicable for pulmonary artery and vein separation. As shown in Figure4, pulmonary arteries may adhere with pulmonary veins, and pulmonary veins may cross lobar fissures (marked with yellow arrows). To remedy the problem, a shape prior and graph cut method [30] is applied to separate and segment multiple interwoven vessel trees.Figure 4
Pulmonary arteries and veins. (a) Pulmonary arteries always adhered with pulmonary veins, and pulmonary veins may cross lobar fissures. (b) The local image corresponding to (a). As shown in the region labeled with yellow arrows, pulmonary arteries and veins are joined together. Furthermore, pulmonary veins may cross pulmonary fissures.Pulmonary artery and vein separation is a challenging task due to its complex structure. To solve the problem, the shape prior and graph cut method [30] is applied to segment the interwoven three-dimensional (3D) tubular trees. The first step is to identify all tubular objects by TDF. Three eigenvectors of the Hessian matrix corresponding to tubular objects are estimated: one of the three eigenvectors e⟶1 is related to the direction of the tube, and the other two eigenvectors e⟶2 and e⟶3 are associated with the cross-sectional plane of the tubular objects. The main idea is to create a circle with a radius r in the cross-sectional plane, and N is the number of spaced points i on the circle. Therefore, the direction from the center to each point i is(9)d⟶i=e⟶2sin2πiN+e⟶3cos2πiN.To segment the interwoven 3D tubular trees, a medialness method [30] is applied to separate pulmonary arteries and veins. The initial boundariness contribution and averaging boundariness contribution are(10)bi=σbV0v⟶+rd⟶id⟶i,b¯σb=1N∑i=0N−1σbV0v⟶+rd⟶id⟶i,where N represents the number of spaced points on the circle, r is the circle radius, di is the tube direction, and the gradient vector field V0 is with σb=rη0≤η≤1. The variance can be calculated as follows:(11)s2v⟶,r=1N∑i=1N−1bi−b¯2.Thus, the response of the final offset medialness is(12)R+v⟶,r=1N∑i=0N−1bi1−s2v⟶,rb¯2.To suppress the response away from the tube center, Bauer et al. [30] proposed a method based on an adaptive threshold:(13)R−v⟶,r=σbV0v⟶.The multiscale medialness response is(14)Rmultiv⟶=maxrmin≤r≤rmaxmaxR+v⟶,r−R−v⟶,r,0.Here,rmin=1 and rmax=5. After acquiring the multiscale medialness response, we set a threshold r=3 and take the medialness response larger than the threshold as seed points. Finally, a tree construction algorithm is used to separate pulmonary arteries and veins. The main idea is to extract the centerline of the multiscale medialness response by traversing a ridge map.Pulmonary arteries may adhere with pulmonary veins in Figure5(a). To separate the multiple interwoven vessel trees, we use the medialness method to extract the centerline of vessel trees, and the results are shown in Figure 5(b). As observed in the regions marked with yellow arrows in Figure 5, pulmonary arteries and veins have been separated.Figure 5
Pulmonary artery and vein separation. (a) Pulmonary arteries adhered with pulmonary veins, and the adhered regions are marked with yellow arrows. (b) Pulmonary arteries and veins are separated.
(a)(b)


### 2.3.3. Fissure ROI Identification

Based on the reality that pulmonary arteries are close to bronchial trees [39–41], we combine the bronchial trees (BT) with vessel subtrees (VST) to remove mounts of pulmonary veins. We dilate BT and use a connected-component analysis method (CCAM) to select 5 biggest objects for pulmonary artery segmentation. To give a mathematical expression,(15)PA=max1≤i≤5CCAMVST+DBT∩VST,where D denotes the morphological dilation.To illustrate the validation of the proposed scheme, we integrate with bronchial trees, pulmonary arteries, and lobar fissures to describe the distribution map from different angles of view. Bronchial trees and lobar fissures are labeled with yellow color and green color in Figure6, respectively. It can be seen that pulmonary arteries are close to bronchial trees [40, 41], and pulmonary arteries cannot cross lobar fissures [15, 39].Figure 6
The distribution map among bronchial trees, pulmonary arteries, and lobar fissures is displayed in different angles of view, in which bronchial trees, pulmonary arteries, and lobar fissures are marked with different colors.
(a)(b)(c)As a result, we consider different bronchial trees and the corresponding pulmonary arteries belonging to different pulmonary lobes as a whole in Figure7(a). Then, an alpha-shape algorithm [31, 32] is used to construct lobe shapes. The alpha-shape algorithm was presented by Edelsbrunner and Mücke [31] that constructs an entire family of shapes for a given set of points. We adopt for the implementation of alpha-shape algorithm to construct the geometric shape of different bronchial trees and their corresponding pulmonary arteries. Finally, surface fitting processing and region filling processing are used to approximately estimate the lobe distribution map for fissure ROI identification in Figure 7(b).Figure 7
Lobe shapes. (a) Bronchial trees and the corresponding pulmonary arteries belonging to different lobes are displayed in different colors. (b) Lobe shapes corresponding to (a). (c) The relationships among lobe shapes, fissures, and bronchial trees. (d, e) The relationships between lobe shapes and fissures.
(a)(b)(c)(d)(e)Fissure ROI identification is useful in lung lobe segmentation. As shown in Figure8(a), the estimated right upper, middle, and lower lobe distribution maps are marked with ST1, ST2, and ST3, respectively. We just regard ST2 and ST3 as an entirety and ST1 as another part; then, the alpha-shape algorithm is applied to identify the right oblique fissure ROI, and its corresponding results are shown in Figure 8(b).Figure 8
Lobe shapes. (a) Right fissure ROI, in which the estimated right upper, middle, and lower lobe distribution maps are marked withST1, ST2, and ST3, respectively. (b) Right oblique fissure ROI.
(a)(b)


## 2.3.1. Airway Segmentation

To accurately segment airways, a useful and valuable method [29] based on GPU is applied for tube extraction. Firstly, Gaussian smoothing with the parameter σ = 0.5 is used to blur the dataset I. Subsequently, a gradient vector field V⟶ is designed:(1)V⟶=−∇I∗σ.Normalize(2)V⟶nv⟶=V⟶v⟶V⟶v⟶ifV⟶nv⟶≥Fmax,V⟶v⟶Fmax,else.Here,Fmax controls the normalization. To prevent different tubular structures diffusing together, a gradient vector flow (GVF) method [33] is used:(3)EV⟶=∫μ∇V⟶v⟶2+V⟶0v⟶V⟶v⟶−V⟶0v⟶2dv⟶2,where V⟶0 and μ denote the initial gradient vector field and the weighting constant. Xu and Prince [33] proposed an iterative approach to solve the Euler equation:(4)μ∇2V⟶−V⟶−V⟶0V⟶02=0.The purpose of the GVF model is to make the following tube detection filter (TDF) [34, 35] invariant to contrast of airways [29]. Subsequently, the TDF is applied to detect tubular structures. Firstly, three eigenvectors of the Hessian matrix corresponding to the tube are estimated: one of three eigenvectors e⟶1 is related to the direction of the tube, and the other two eigenvectors e⟶2 and e⟶3 are associated with the cross-sectional plane of the tube. The main idea is to create a circle with a radius r in the cross-sectional plane, and N is the number of spaced points i on the circle. Therefore, the direction from the center to each point i is(5)d⟶i=e⟶2sin2πiN+e⟶3cos2πiN.The average dot product is(6)Tv⟶,r=1N∑i=0N−1v⟶+rd⟶i⋅−d⟶i.Whenr is increased, the average dot product is calculated again. In this paper, r ranges from 0.5 to 25. Lastly, the airway centerline is extracted by traversing a ridge in the TDF results. All values in TDF results above a threshold are extracted and treated as candidate starting points, starting from the candidate point with the positive tangent direction t0 and negative −t0. For a current centerline point xi and the corresponding tangent direction ti, all local neighbors xin are considered:(7)xixin⟶⋅ti>0.The neighbor point with the highest TDF response is chosen as the next centerline pointxi+1. Finally, the airway centerline is acquired and marked with S.Different subtrees in the airway tree belong to different lobes [15]. Based on this theory, a graph search method is used to search for major bifurcations for bronchial tree separation in appropriate orientations in Figure 3. In the first step, a centerline extraction method [29] is used to transfer the airway tree segmentation to its skeleton representation. Then, a directed graph is modeled from the airway tree with the trachea as the root. To extract the bifurcations, a simple but effective criterion [37, 38] is used:(8)N26p∩S∪p≥4,where p and N26 denote the pixel in the centerline of the airway tree S and its 26-neighborhood area, respectively. If there are 4 or more points within a 3 ∗ 3 ∗ 3 cube (p and its 26 neighborhood), we consider the pixel p as a bifurcation. Finally, a graph search method is used to remove major bifurcations labeled with yellow dots in Figure 3(b) and calculate for each subtree in Figure 3(c) along with the particular separation direction n⟶ in Figure 3(b), in which EndRMB and LMB denote the ends in the right main bronchus and the left main bronchus, respectively.Figure 3
Subtrees in an airway tree. (a) An airway tree. (b) Graphical representation of a typical skeletonized airway tree (figure taken from [36]). (c) Skeleton of an airway tree; subtrees belonging to different lobes are displayed in different colors.
(a)(b)(c)


## 2.3.2. Pulmonary Artery and Vein Separation

Airways are tubular structures and less adhered with each other. Furthermore, the GVF model has the advantage that it is feature preserving. Therefore, a method based on the GVF model is used to segment the airway tree, but the approach is not applicable for pulmonary artery and vein separation. As shown in Figure4, pulmonary arteries may adhere with pulmonary veins, and pulmonary veins may cross lobar fissures (marked with yellow arrows). To remedy the problem, a shape prior and graph cut method [30] is applied to separate and segment multiple interwoven vessel trees.Figure 4
Pulmonary arteries and veins. (a) Pulmonary arteries always adhered with pulmonary veins, and pulmonary veins may cross lobar fissures. (b) The local image corresponding to (a). As shown in the region labeled with yellow arrows, pulmonary arteries and veins are joined together. Furthermore, pulmonary veins may cross pulmonary fissures.Pulmonary artery and vein separation is a challenging task due to its complex structure. To solve the problem, the shape prior and graph cut method [30] is applied to segment the interwoven three-dimensional (3D) tubular trees. The first step is to identify all tubular objects by TDF. Three eigenvectors of the Hessian matrix corresponding to tubular objects are estimated: one of the three eigenvectors e⟶1 is related to the direction of the tube, and the other two eigenvectors e⟶2 and e⟶3 are associated with the cross-sectional plane of the tubular objects. The main idea is to create a circle with a radius r in the cross-sectional plane, and N is the number of spaced points i on the circle. Therefore, the direction from the center to each point i is(9)d⟶i=e⟶2sin2πiN+e⟶3cos2πiN.To segment the interwoven 3D tubular trees, a medialness method [30] is applied to separate pulmonary arteries and veins. The initial boundariness contribution and averaging boundariness contribution are(10)bi=σbV0v⟶+rd⟶id⟶i,b¯σb=1N∑i=0N−1σbV0v⟶+rd⟶id⟶i,where N represents the number of spaced points on the circle, r is the circle radius, di is the tube direction, and the gradient vector field V0 is with σb=rη0≤η≤1. The variance can be calculated as follows:(11)s2v⟶,r=1N∑i=1N−1bi−b¯2.Thus, the response of the final offset medialness is(12)R+v⟶,r=1N∑i=0N−1bi1−s2v⟶,rb¯2.To suppress the response away from the tube center, Bauer et al. [30] proposed a method based on an adaptive threshold:(13)R−v⟶,r=σbV0v⟶.The multiscale medialness response is(14)Rmultiv⟶=maxrmin≤r≤rmaxmaxR+v⟶,r−R−v⟶,r,0.Here,rmin=1 and rmax=5. After acquiring the multiscale medialness response, we set a threshold r=3 and take the medialness response larger than the threshold as seed points. Finally, a tree construction algorithm is used to separate pulmonary arteries and veins. The main idea is to extract the centerline of the multiscale medialness response by traversing a ridge map.Pulmonary arteries may adhere with pulmonary veins in Figure5(a). To separate the multiple interwoven vessel trees, we use the medialness method to extract the centerline of vessel trees, and the results are shown in Figure 5(b). As observed in the regions marked with yellow arrows in Figure 5, pulmonary arteries and veins have been separated.Figure 5
Pulmonary artery and vein separation. (a) Pulmonary arteries adhered with pulmonary veins, and the adhered regions are marked with yellow arrows. (b) Pulmonary arteries and veins are separated.
(a)(b)


## 2.3.3. Fissure ROI Identification

Based on the reality that pulmonary arteries are close to bronchial trees [39–41], we combine the bronchial trees (BT) with vessel subtrees (VST) to remove mounts of pulmonary veins. We dilate BT and use a connected-component analysis method (CCAM) to select 5 biggest objects for pulmonary artery segmentation. To give a mathematical expression,(15)PA=max1≤i≤5CCAMVST+DBT∩VST,where D denotes the morphological dilation.To illustrate the validation of the proposed scheme, we integrate with bronchial trees, pulmonary arteries, and lobar fissures to describe the distribution map from different angles of view. Bronchial trees and lobar fissures are labeled with yellow color and green color in Figure6, respectively. It can be seen that pulmonary arteries are close to bronchial trees [40, 41], and pulmonary arteries cannot cross lobar fissures [15, 39].Figure 6
The distribution map among bronchial trees, pulmonary arteries, and lobar fissures is displayed in different angles of view, in which bronchial trees, pulmonary arteries, and lobar fissures are marked with different colors.
(a)(b)(c)As a result, we consider different bronchial trees and the corresponding pulmonary arteries belonging to different pulmonary lobes as a whole in Figure7(a). Then, an alpha-shape algorithm [31, 32] is used to construct lobe shapes. The alpha-shape algorithm was presented by Edelsbrunner and Mücke [31] that constructs an entire family of shapes for a given set of points. We adopt for the implementation of alpha-shape algorithm to construct the geometric shape of different bronchial trees and their corresponding pulmonary arteries. Finally, surface fitting processing and region filling processing are used to approximately estimate the lobe distribution map for fissure ROI identification in Figure 7(b).Figure 7
Lobe shapes. (a) Bronchial trees and the corresponding pulmonary arteries belonging to different lobes are displayed in different colors. (b) Lobe shapes corresponding to (a). (c) The relationships among lobe shapes, fissures, and bronchial trees. (d, e) The relationships between lobe shapes and fissures.
(a)(b)(c)(d)(e)Fissure ROI identification is useful in lung lobe segmentation. As shown in Figure8(a), the estimated right upper, middle, and lower lobe distribution maps are marked with ST1, ST2, and ST3, respectively. We just regard ST2 and ST3 as an entirety and ST1 as another part; then, the alpha-shape algorithm is applied to identify the right oblique fissure ROI, and its corresponding results are shown in Figure 8(b).Figure 8
Lobe shapes. (a) Right fissure ROI, in which the estimated right upper, middle, and lower lobe distribution maps are marked withST1, ST2, and ST3, respectively. (b) Right oblique fissure ROI.
(a)(b)


## 2.4. Lobar Fissure Verification

To effectively verify pulmonary fissures, an oriented derivative of stick (ODoS) filter [26] is applied to enhance pulmonary fissures. Then, an orientation partition scheme is used to separate pulmonary fissures and adhering clutters. Next, a connected-component analysis approach is adopted to remove noise for plate-like structure extraction. Finally, the identified fissure ROI is employed to isolate lobar fissures from adhering clutters.


### 2.4.1. Retrospection of the DoS Filter and ODoS Filter

To efficiently detect the weak and thin pulmonary fissures, Xiao et al. [25] presented a derivative of stick (DoS) filter for pulmonary fissure enhancement. The filtering kernel consists of three parallel sticks (LS, MS, and RS), and θ is the stick orientation. Two 2D line strength measures lmaxθ and lminθ are(16)lmaxθx=max1L∑j=1LIjM−1L∑j=1LIjL,1L∑j=1LIjM−1L∑j=1LIjR−k⋅EIjM2−EIjM2,lminθx=min1L∑j=1LIjM−1L∑j=1LIjL,1L∑j=1LIjM−1L∑j=1LIjR−k⋅EIjM2−EIjM2,where IjL, IjM, and IjR denote the intensity of the jth pixel along the left, middle, and right stick, k=0.7 represents a positive coefficient, and E is an expected value operator. To suppress noise and other structures, Xiao et al. [25] integrated the line strength measures lmaxθ and lminθ in multidirections:(17)Fmaxx=maxmax0<θ≤πlmaxθ,0,Fminx=maxmax0<θ≤πlminθ,0.The DoS filtering response can be formulated as(18)Fο=FmaxοFmin.ο represents the cascading operator. To integrate the DoS filtering response in multidirections, Xiao et al. defined a shape-tuned fissureness function:(19)F3Dx=FοA+FοS+FοC⋅medianFοA,FοS,FοCmaxFοA,FοS,FοC.Here,FοA, FοS, and FοC denote the DoS filtering response Fο across the axial, sagittal, and coronal sections. However, the original DoS filter cannot efficiently suppress adhering clutters. To solve the problem, Peng and Xiao [26] improved and extended it by constructing its orientation response. The orientation response is(20)θmax=argmax0<θ≤πlmaxθiflmaxθ>0,NaN,else.A vector representation is(21)V⟶max=cosθmax,sinθmaxiflmaxθ>0,0,else.Therefore, the response of the ODoS filter can be described as follows:(22)V⟶F=F3D⋅V⟶max.


### 2.4.2. Plate-Like Structure Extraction

To suppress noise and other structures, Peng and Xiao [26] firstly used a threshold T to modulate the ODoS filtering response:(23)V⟶T=V⟶FifV⟶F>T,0,else,where the threshold T we set is equal to 1. Based on the fact that pulmonary fissures appear as plate-like structures and have similar directions, Peng and Xiao [26] presented an orientation partition scheme to separate lobar fissures and noise in Figure 8. They divided the enhanced 3D CT image into n overlapped subregions, and the subregions are marked with R1,R2,…Rn−1 and Rn. To give a mathematical expression,(24)FRi=1ifθmax∈RiandV⟶T≠0,0,else.In 3D partitioned imagesFRi, clutters such as airways and vessels usually appear as small objects. As a result, they used a connected-component analysis method to select the candidate fissures, and the results are displayed in Figure 9.Figure 9
Plate-like structure extraction in the right lung. After using an orientation partition scheme, the enhanced images are divided into a number of planar structures.


### 2.4.3. Oblique Fissure and Horizontal Fissure Verification

From the above processing, some plate-like structures are easily acquired, including lobar fissures and accessory fissures. As is shown in Figure10, accessory fissures have bad influence on fissure segmentation, but few methods have been proposed by using the priori knowledge of pulmonary structures to remove them. In this study, the identified fissure ROI is used to limit the search area and remove those unrelated accessory fissures for right fissure verification. Similarly, we can use the previously identified fissure ROI to verify right oblique fissures and left oblique fissures.Figure 10
Lobar fissure verification in the right lung, in which fissure ROI is used to limit the search area for fissure verification.After acquiring the right fissure (RF) and the right oblique fissure (ROF), we use a simple but effective criterion to verify the right horizontal fissure (RHF):(25)RHF=RF−ROF.As shown in Figure11, it can be seen that the proposed method has a good performance in lobar fissure verification.Figure 11
Right oblique fissure and horizontal fissure separation. (a) Right fissure (RF). (b) Right oblique fissure (ROF). (c) Right horizontal fissure (RHF).
(a)(b)(c)


## 2.4.1. Retrospection of the DoS Filter and ODoS Filter

To efficiently detect the weak and thin pulmonary fissures, Xiao et al. [25] presented a derivative of stick (DoS) filter for pulmonary fissure enhancement. The filtering kernel consists of three parallel sticks (LS, MS, and RS), and θ is the stick orientation. Two 2D line strength measures lmaxθ and lminθ are(16)lmaxθx=max1L∑j=1LIjM−1L∑j=1LIjL,1L∑j=1LIjM−1L∑j=1LIjR−k⋅EIjM2−EIjM2,lminθx=min1L∑j=1LIjM−1L∑j=1LIjL,1L∑j=1LIjM−1L∑j=1LIjR−k⋅EIjM2−EIjM2,where IjL, IjM, and IjR denote the intensity of the jth pixel along the left, middle, and right stick, k=0.7 represents a positive coefficient, and E is an expected value operator. To suppress noise and other structures, Xiao et al. [25] integrated the line strength measures lmaxθ and lminθ in multidirections:(17)Fmaxx=maxmax0<θ≤πlmaxθ,0,Fminx=maxmax0<θ≤πlminθ,0.The DoS filtering response can be formulated as(18)Fο=FmaxοFmin.ο represents the cascading operator. To integrate the DoS filtering response in multidirections, Xiao et al. defined a shape-tuned fissureness function:(19)F3Dx=FοA+FοS+FοC⋅medianFοA,FοS,FοCmaxFοA,FοS,FοC.Here,FοA, FοS, and FοC denote the DoS filtering response Fο across the axial, sagittal, and coronal sections. However, the original DoS filter cannot efficiently suppress adhering clutters. To solve the problem, Peng and Xiao [26] improved and extended it by constructing its orientation response. The orientation response is(20)θmax=argmax0<θ≤πlmaxθiflmaxθ>0,NaN,else.A vector representation is(21)V⟶max=cosθmax,sinθmaxiflmaxθ>0,0,else.Therefore, the response of the ODoS filter can be described as follows:(22)V⟶F=F3D⋅V⟶max.


## 2.4.2. Plate-Like Structure Extraction

To suppress noise and other structures, Peng and Xiao [26] firstly used a threshold T to modulate the ODoS filtering response:(23)V⟶T=V⟶FifV⟶F>T,0,else,where the threshold T we set is equal to 1. Based on the fact that pulmonary fissures appear as plate-like structures and have similar directions, Peng and Xiao [26] presented an orientation partition scheme to separate lobar fissures and noise in Figure 8. They divided the enhanced 3D CT image into n overlapped subregions, and the subregions are marked with R1,R2,…Rn−1 and Rn. To give a mathematical expression,(24)FRi=1ifθmax∈RiandV⟶T≠0,0,else.In 3D partitioned imagesFRi, clutters such as airways and vessels usually appear as small objects. As a result, they used a connected-component analysis method to select the candidate fissures, and the results are displayed in Figure 9.Figure 9
Plate-like structure extraction in the right lung. After using an orientation partition scheme, the enhanced images are divided into a number of planar structures.


## 2.4.3. Oblique Fissure and Horizontal Fissure Verification

From the above processing, some plate-like structures are easily acquired, including lobar fissures and accessory fissures. As is shown in Figure10, accessory fissures have bad influence on fissure segmentation, but few methods have been proposed by using the priori knowledge of pulmonary structures to remove them. In this study, the identified fissure ROI is used to limit the search area and remove those unrelated accessory fissures for right fissure verification. Similarly, we can use the previously identified fissure ROI to verify right oblique fissures and left oblique fissures.Figure 10
Lobar fissure verification in the right lung, in which fissure ROI is used to limit the search area for fissure verification.After acquiring the right fissure (RF) and the right oblique fissure (ROF), we use a simple but effective criterion to verify the right horizontal fissure (RHF):(25)RHF=RF−ROF.As shown in Figure11, it can be seen that the proposed method has a good performance in lobar fissure verification.Figure 11
Right oblique fissure and horizontal fissure separation. (a) Right fissure (RF). (b) Right oblique fissure (ROF). (c) Right horizontal fissure (RHF).
(a)(b)(c)


## 2.5. Lung Lobe Segmentation

Lung lobe segmentation is an arduous task due to incomplete, deformed, and disrupted fissures. To remedy the problem, a surface fitting model [6, 11, 12] is applied to construct a complete surface. Subsequently, the complete surface is employed to segment lung lobes.


### 2.5.1. Surface Fitting Model

Once oblique fissures and horizontal fissures are verified, we modeled a fissure surface in the form of(26)z=fx,y,where x, y, and z are the coordinate points in the 3D space. Imaging a thin and flexible plate for a given set of points by elastic bands is a nontrivial task. Each of these elastic bands draws the plate towards the given set of points, stretching only in the Z-direction. The given set of points was to fall on a simple planar surface [6, 11]:(27)z=w0+w1x1+w2y2+e.A matrix expression is(28)Z=AW+e.At every point of the surface, the surface fitting model attempts to force the partial derivatives of the surface in neighboring cells to be equal. A mathematical expression in the second set of linear equations is of the form(29)BW=0.To construct a smooth fissure surface, a Laplacian function is treated as the regularization term. It can be described as(30)B=α∂2fx,y∂2x∂2fx,y∂2y,where α is a scale factor that controls the effect of regularization. Then, we solve for the vector W by performing the following minimization:(31)Wmin=argminAW−Z2+BW2.Therefore, the surface fitting model is employed to construct a complete fissure surface, as shown in Figure12.Figure 12
Surface fitting model. (a) Right oblique fissure. (b) The surface fitting results corresponding to (a). (c) Right horizontal fissure. (d) The surface fitting results corresponding to (c).
(a)(b)(c)(d)


### 2.5.2. Cutting Model

To segment pulmonary lobes, we present a cutting model to subdivide human lungs into five lobes. As shown in Figure13, the right complete oblique fissure surface is used to divide the right lung into the special part and the lower lobe, and then the right complete horizontal fissure surface is employed to separate the special part into the upper and middle lobes. In the same way, we use the left complete fissure surface to divide the left lung into left upper and lower lobes. As a result, the human lung is subdivided into five lobes.Figure 13
Pulmonary lobe separation. (a) Right complete oblique fissure surface divides the right lung into two parts. (b) Right complete horizontal fissure surface divides the special part into two parts. (c) The flowchart of the right lung lobe segmentation.
(a)(b)(c)


## 2.5.1. Surface Fitting Model

Once oblique fissures and horizontal fissures are verified, we modeled a fissure surface in the form of(26)z=fx,y,where x, y, and z are the coordinate points in the 3D space. Imaging a thin and flexible plate for a given set of points by elastic bands is a nontrivial task. Each of these elastic bands draws the plate towards the given set of points, stretching only in the Z-direction. The given set of points was to fall on a simple planar surface [6, 11]:(27)z=w0+w1x1+w2y2+e.A matrix expression is(28)Z=AW+e.At every point of the surface, the surface fitting model attempts to force the partial derivatives of the surface in neighboring cells to be equal. A mathematical expression in the second set of linear equations is of the form(29)BW=0.To construct a smooth fissure surface, a Laplacian function is treated as the regularization term. It can be described as(30)B=α∂2fx,y∂2x∂2fx,y∂2y,where α is a scale factor that controls the effect of regularization. Then, we solve for the vector W by performing the following minimization:(31)Wmin=argminAW−Z2+BW2.Therefore, the surface fitting model is employed to construct a complete fissure surface, as shown in Figure12.Figure 12
Surface fitting model. (a) Right oblique fissure. (b) The surface fitting results corresponding to (a). (c) Right horizontal fissure. (d) The surface fitting results corresponding to (c).
(a)(b)(c)(d)


## 2.5.2. Cutting Model

To segment pulmonary lobes, we present a cutting model to subdivide human lungs into five lobes. As shown in Figure13, the right complete oblique fissure surface is used to divide the right lung into the special part and the lower lobe, and then the right complete horizontal fissure surface is employed to separate the special part into the upper and middle lobes. In the same way, we use the left complete fissure surface to divide the left lung into left upper and lower lobes. As a result, the human lung is subdivided into five lobes.Figure 13
Pulmonary lobe separation. (a) Right complete oblique fissure surface divides the right lung into two parts. (b) Right complete horizontal fissure surface divides the special part into two parts. (c) The flowchart of the right lung lobe segmentation.
(a)(b)(c)


## 3. Experimental Results

The proposed scheme is tested on 15 clinical CT images. All experiments are performed in 64 bit Windows 10 operating system, running on a computer with an Intel (R) Xeon (R) E5-1607 v2 3.00GHZ CPU and 20 GB RAM. The runtime of the proposed method for a 400 × 512 × 512-size 3D image is approximately 1580 s. In addition, the 3D visualization results are rendered by MeVisLab software [42].


### 3.1. Evaluation of Pulmonary Fissure Segmentation

#### 3.1.1. Evaluation Criteria

The manual reference and its adjacent 3 mm region are labeled withS1, and the segmented fissure and its adjacent 3 mm region are labeled with S2. The voxels of the segmented fissure inside S1 are considered true positive (TP1), and the rest are false positive (FP). Similarly, the voxels of the manual reference inside S2 are considered TP2, and the rest are false negative (FN). As a result, the false discovery rate (FDR), false negative rate (FNR), and F1-score (F1) are given as follows:(32)FDR=FPTP1+FP,FNR=FNTP2+FN,F1=21−FDR1−FNR2−FDR−FNR.


#### 3.1.2. Visual Inspection

As shown in Figure14, the segmentation results of the proposed method are shown in Figure 14(a). For comparison, the ODoS [26] and DoS [25] filtering segmentations are, respectively, displayed in Figures 14(b) and 14(c). The segmented results, the manual reference, and their overlapped regions are, respectively, rendered in green, yellow, and purple. In the areas labeled with red ellipses in Figure 14, the weak objects are treated as the undersegmentation error, which may lead to a high FNR index. As observed, the proposed method has a good performance in weak object detection.Figure 14
Fissure segmentation. (a) The proposed method. (b) ODoS method. (c) DoS method. The weak fissures are marked with red ellipses.
(a)(b)(c)


#### 3.1.3. Quantitative Evaluation

The proposed method is validated on 15 CT examinations in Figure15. The box plots of indices (F1, FDR, and FNR) corresponding to the proposed method (p), ODoS (od), and DoS (d) filtering scheme are drawn next to each other. In the left lung, the median values of box plots are 0.8865, 0.8828, 0.8374, 0.0755, 0.0489, 0.1971, 0.0978, 0.1368, and 0.1478. In the right lung, the median values of box plots are 0.9200, 0.8997, 0.8271, 0.0381, 0.0344, 0.1544, 0.0933, 0.1221, and 0.1906. Both visual inspection and quantitative evaluation exhibited that our scheme can outperform the compared methods [25, 26].Figure 15
Segmentation validation on 15 CT examinations. (a) Left lung. (b) Right lung.
(a)(b)


### 3.2. Evaluation of Lung Lobe Segmentation

#### 3.2.1. Evaluation Criteria

The volumes of the lung lobe ground truth and the segmented lung lobes are, respectively, marked withh1 and h2. The volumetric overlap results between h1 and h2 are marked with h3. We regard the percentage η=h3/h1 as the evaluation criteria.


#### 3.2.2. Visual Inspection

To provide a global impression, four representative segmentations are chosen from the LOLA11 dataset. As shown in Figure16, the upper, middle, and lower lobes are rendered in green, red, and yellow, respectively. The first, second, and third rows represent the visual lobar boundaries, the proposed method, and the constrained interpolation profile method [43]. As observed in the areas labeled with red arrows, the proposed method has a good performance in lung lobe segmentation.Figure 16
The segmented lung lobes are marked with different colors. The upper, middle, and lower lobes are rendered in green, red, and yellow, respectively. The first, second, and third rows represent the visual lobar boundaries, the proposed method, and the constrained interpolation profile method, in which the distinct areas of pulmonary lobes are marked with red arrows.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)(k)(l)


#### 3.2.3. Quantitative Evaluation

Table1 presents the evaluation results of lung lobes in 15 patients. The average percentages of the segmented LUL, LLL, RUL, RML, and RLL in the lung lobe ground truth are 0.960, 0.989, 0.973, 0.920, and 0.985, respectively. Lung lobe segmentation is useful in clinical diagnosis, especially in rapid diagnosis of new coronary pneumonia [44], in which we divide all patients into two parts: normal and abnormal groups. In abnormal groups, parts of bronchus and fissures are missed due to various pathological abnormalities. As shown in Figure 17(a), the disrupted bronchus has no power to guide and segment pulmonary arteries, and it may cause parts of clutters not to be removed. Another limitation of our scheme is drawn on incomplete fissures. As observed in Figure 17(b), incomplete fissures are used to construct complete lobe boundaries by using the surface fitting model; it may lead to a wrong lobar segmentation. Furthermore, the results can be seen in Table 1, and the abnormal cases have low volumetric overlap values with bold type.Table 1
Volumetric overlapped results. In this section, LUL, LLL, RUL, RML, and RLL are used to represent the left upper lobe, left lower lobe, right upper lobe, right middle lobe, and right lower lobe, respectively. We divide all cases into two parts: normal groups and abnormal groups.
LULLLLRULRMLRLLGroupCase 10.9940.9700.9850.9810.990NormalCase 20.9370.9990.9990.9980.977NormalCase 30.9960.9940.9980.9690.986NormalCase 40.9920.9900.9810.9540.960NormalCase 50.9860.9920.9840.9180.985NormalCase 60.9900.9720.9920.9280.987NormalCase 70.9330.9900.9940.9440.993NormalCase 80.9990.9940.9840.9690.991NormalCase 90.8340.9920.9900.9470.992AbnormalCase 100.9980.9860.9820.8640.979AbnormalCase 110.9890.9950.8910.9560.989AbnormalCase 120.9980.9920.9480.9150.990AbnormalCase 130.9880.9990.9990.8450.980AbnormalCase 140.9990.9890.8830.9420.995AbnormalCase 150.9570.9800.9850.6630.985AbnormalAverage0.9600.9890.9730.9200.985Maximum0.9990.9990.9990.9980.995Minimum0.8340.9700.8830.6630.960Median0.9900.9920.9850.9440.987Standard deviation0.0580.0090.0370.0820.009Figure 17
The limitations of our scheme. (a) Airways. (b) Fissures, in which the disrupted bronchus is marked with red arrows, and the incomplete fissure is marked with yellow arrows.
(a)(b)


## 3.1. Evaluation of Pulmonary Fissure Segmentation

### 3.1.1. Evaluation Criteria

The manual reference and its adjacent 3 mm region are labeled withS1, and the segmented fissure and its adjacent 3 mm region are labeled with S2. The voxels of the segmented fissure inside S1 are considered true positive (TP1), and the rest are false positive (FP). Similarly, the voxels of the manual reference inside S2 are considered TP2, and the rest are false negative (FN). As a result, the false discovery rate (FDR), false negative rate (FNR), and F1-score (F1) are given as follows:(32)FDR=FPTP1+FP,FNR=FNTP2+FN,F1=21−FDR1−FNR2−FDR−FNR.


### 3.1.2. Visual Inspection

As shown in Figure14, the segmentation results of the proposed method are shown in Figure 14(a). For comparison, the ODoS [26] and DoS [25] filtering segmentations are, respectively, displayed in Figures 14(b) and 14(c). The segmented results, the manual reference, and their overlapped regions are, respectively, rendered in green, yellow, and purple. In the areas labeled with red ellipses in Figure 14, the weak objects are treated as the undersegmentation error, which may lead to a high FNR index. As observed, the proposed method has a good performance in weak object detection.Figure 14
Fissure segmentation. (a) The proposed method. (b) ODoS method. (c) DoS method. The weak fissures are marked with red ellipses.
(a)(b)(c)


### 3.1.3. Quantitative Evaluation

The proposed method is validated on 15 CT examinations in Figure15. The box plots of indices (F1, FDR, and FNR) corresponding to the proposed method (p), ODoS (od), and DoS (d) filtering scheme are drawn next to each other. In the left lung, the median values of box plots are 0.8865, 0.8828, 0.8374, 0.0755, 0.0489, 0.1971, 0.0978, 0.1368, and 0.1478. In the right lung, the median values of box plots are 0.9200, 0.8997, 0.8271, 0.0381, 0.0344, 0.1544, 0.0933, 0.1221, and 0.1906. Both visual inspection and quantitative evaluation exhibited that our scheme can outperform the compared methods [25, 26].Figure 15
Segmentation validation on 15 CT examinations. (a) Left lung. (b) Right lung.
(a)(b)


## 3.1.1. Evaluation Criteria

The manual reference and its adjacent 3 mm region are labeled withS1, and the segmented fissure and its adjacent 3 mm region are labeled with S2. The voxels of the segmented fissure inside S1 are considered true positive (TP1), and the rest are false positive (FP). Similarly, the voxels of the manual reference inside S2 are considered TP2, and the rest are false negative (FN). As a result, the false discovery rate (FDR), false negative rate (FNR), and F1-score (F1) are given as follows:(32)FDR=FPTP1+FP,FNR=FNTP2+FN,F1=21−FDR1−FNR2−FDR−FNR.


## 3.1.2. Visual Inspection

As shown in Figure14, the segmentation results of the proposed method are shown in Figure 14(a). For comparison, the ODoS [26] and DoS [25] filtering segmentations are, respectively, displayed in Figures 14(b) and 14(c). The segmented results, the manual reference, and their overlapped regions are, respectively, rendered in green, yellow, and purple. In the areas labeled with red ellipses in Figure 14, the weak objects are treated as the undersegmentation error, which may lead to a high FNR index. As observed, the proposed method has a good performance in weak object detection.Figure 14
Fissure segmentation. (a) The proposed method. (b) ODoS method. (c) DoS method. The weak fissures are marked with red ellipses.
(a)(b)(c)


## 3.1.3. Quantitative Evaluation

The proposed method is validated on 15 CT examinations in Figure15. The box plots of indices (F1, FDR, and FNR) corresponding to the proposed method (p), ODoS (od), and DoS (d) filtering scheme are drawn next to each other. In the left lung, the median values of box plots are 0.8865, 0.8828, 0.8374, 0.0755, 0.0489, 0.1971, 0.0978, 0.1368, and 0.1478. In the right lung, the median values of box plots are 0.9200, 0.8997, 0.8271, 0.0381, 0.0344, 0.1544, 0.0933, 0.1221, and 0.1906. Both visual inspection and quantitative evaluation exhibited that our scheme can outperform the compared methods [25, 26].Figure 15
Segmentation validation on 15 CT examinations. (a) Left lung. (b) Right lung.
(a)(b)


## 3.2. Evaluation of Lung Lobe Segmentation

### 3.2.1. Evaluation Criteria

The volumes of the lung lobe ground truth and the segmented lung lobes are, respectively, marked withh1 and h2. The volumetric overlap results between h1 and h2 are marked with h3. We regard the percentage η=h3/h1 as the evaluation criteria.


### 3.2.2. Visual Inspection

To provide a global impression, four representative segmentations are chosen from the LOLA11 dataset. As shown in Figure16, the upper, middle, and lower lobes are rendered in green, red, and yellow, respectively. The first, second, and third rows represent the visual lobar boundaries, the proposed method, and the constrained interpolation profile method [43]. As observed in the areas labeled with red arrows, the proposed method has a good performance in lung lobe segmentation.Figure 16
The segmented lung lobes are marked with different colors. The upper, middle, and lower lobes are rendered in green, red, and yellow, respectively. The first, second, and third rows represent the visual lobar boundaries, the proposed method, and the constrained interpolation profile method, in which the distinct areas of pulmonary lobes are marked with red arrows.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)(k)(l)


### 3.2.3. Quantitative Evaluation

Table1 presents the evaluation results of lung lobes in 15 patients. The average percentages of the segmented LUL, LLL, RUL, RML, and RLL in the lung lobe ground truth are 0.960, 0.989, 0.973, 0.920, and 0.985, respectively. Lung lobe segmentation is useful in clinical diagnosis, especially in rapid diagnosis of new coronary pneumonia [44], in which we divide all patients into two parts: normal and abnormal groups. In abnormal groups, parts of bronchus and fissures are missed due to various pathological abnormalities. As shown in Figure 17(a), the disrupted bronchus has no power to guide and segment pulmonary arteries, and it may cause parts of clutters not to be removed. Another limitation of our scheme is drawn on incomplete fissures. As observed in Figure 17(b), incomplete fissures are used to construct complete lobe boundaries by using the surface fitting model; it may lead to a wrong lobar segmentation. Furthermore, the results can be seen in Table 1, and the abnormal cases have low volumetric overlap values with bold type.Table 1
Volumetric overlapped results. In this section, LUL, LLL, RUL, RML, and RLL are used to represent the left upper lobe, left lower lobe, right upper lobe, right middle lobe, and right lower lobe, respectively. We divide all cases into two parts: normal groups and abnormal groups.
LULLLLRULRMLRLLGroupCase 10.9940.9700.9850.9810.990NormalCase 20.9370.9990.9990.9980.977NormalCase 30.9960.9940.9980.9690.986NormalCase 40.9920.9900.9810.9540.960NormalCase 50.9860.9920.9840.9180.985NormalCase 60.9900.9720.9920.9280.987NormalCase 70.9330.9900.9940.9440.993NormalCase 80.9990.9940.9840.9690.991NormalCase 90.8340.9920.9900.9470.992AbnormalCase 100.9980.9860.9820.8640.979AbnormalCase 110.9890.9950.8910.9560.989AbnormalCase 120.9980.9920.9480.9150.990AbnormalCase 130.9880.9990.9990.8450.980AbnormalCase 140.9990.9890.8830.9420.995AbnormalCase 150.9570.9800.9850.6630.985AbnormalAverage0.9600.9890.9730.9200.985Maximum0.9990.9990.9990.9980.995Minimum0.8340.9700.8830.6630.960Median0.9900.9920.9850.9440.987Standard deviation0.0580.0090.0370.0820.009Figure 17
The limitations of our scheme. (a) Airways. (b) Fissures, in which the disrupted bronchus is marked with red arrows, and the incomplete fissure is marked with yellow arrows.
(a)(b)


## 3.2.1. Evaluation Criteria

The volumes of the lung lobe ground truth and the segmented lung lobes are, respectively, marked withh1 and h2. The volumetric overlap results between h1 and h2 are marked with h3. We regard the percentage η=h3/h1 as the evaluation criteria.


## 3.2.2. Visual Inspection

To provide a global impression, four representative segmentations are chosen from the LOLA11 dataset. As shown in Figure16, the upper, middle, and lower lobes are rendered in green, red, and yellow, respectively. The first, second, and third rows represent the visual lobar boundaries, the proposed method, and the constrained interpolation profile method [43]. As observed in the areas labeled with red arrows, the proposed method has a good performance in lung lobe segmentation.Figure 16
The segmented lung lobes are marked with different colors. The upper, middle, and lower lobes are rendered in green, red, and yellow, respectively. The first, second, and third rows represent the visual lobar boundaries, the proposed method, and the constrained interpolation profile method, in which the distinct areas of pulmonary lobes are marked with red arrows.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)(k)(l)


## 3.2.3. Quantitative Evaluation

Table1 presents the evaluation results of lung lobes in 15 patients. The average percentages of the segmented LUL, LLL, RUL, RML, and RLL in the lung lobe ground truth are 0.960, 0.989, 0.973, 0.920, and 0.985, respectively. Lung lobe segmentation is useful in clinical diagnosis, especially in rapid diagnosis of new coronary pneumonia [44], in which we divide all patients into two parts: normal and abnormal groups. In abnormal groups, parts of bronchus and fissures are missed due to various pathological abnormalities. As shown in Figure 17(a), the disrupted bronchus has no power to guide and segment pulmonary arteries, and it may cause parts of clutters not to be removed. Another limitation of our scheme is drawn on incomplete fissures. As observed in Figure 17(b), incomplete fissures are used to construct complete lobe boundaries by using the surface fitting model; it may lead to a wrong lobar segmentation. Furthermore, the results can be seen in Table 1, and the abnormal cases have low volumetric overlap values with bold type.Table 1
Volumetric overlapped results. In this section, LUL, LLL, RUL, RML, and RLL are used to represent the left upper lobe, left lower lobe, right upper lobe, right middle lobe, and right lower lobe, respectively. We divide all cases into two parts: normal groups and abnormal groups.
LULLLLRULRMLRLLGroupCase 10.9940.9700.9850.9810.990NormalCase 20.9370.9990.9990.9980.977NormalCase 30.9960.9940.9980.9690.986NormalCase 40.9920.9900.9810.9540.960NormalCase 50.9860.9920.9840.9180.985NormalCase 60.9900.9720.9920.9280.987NormalCase 70.9330.9900.9940.9440.993NormalCase 80.9990.9940.9840.9690.991NormalCase 90.8340.9920.9900.9470.992AbnormalCase 100.9980.9860.9820.8640.979AbnormalCase 110.9890.9950.8910.9560.989AbnormalCase 120.9980.9920.9480.9150.990AbnormalCase 130.9880.9990.9990.8450.980AbnormalCase 140.9990.9890.8830.9420.995AbnormalCase 150.9570.9800.9850.6630.985AbnormalAverage0.9600.9890.9730.9200.985Maximum0.9990.9990.9990.9980.995Minimum0.8340.9700.8830.6630.960Median0.9900.9920.9850.9440.987Standard deviation0.0580.0090.0370.0820.009Figure 17
The limitations of our scheme. (a) Airways. (b) Fissures, in which the disrupted bronchus is marked with red arrows, and the incomplete fissure is marked with yellow arrows.
(a)(b)


## 4. Discussion

In this paper, a unique scheme based on lung anatomy knowledge is presented for lung lobe segmentation. The described scheme has many benefits. First, the priori knowledge of airways and pulmonary arteries is used to identify fissure ROI. Many methods [13, 14, 22, 23] use only pulmonary airways and vessels to identify fissure ROI. Second, the identified fissure ROI and our previous works [25, 26] are combined to accurately verify oblique fissures and horizontal fissures. While parts of methods use the orientation information to verify oblique fissures and horizontal fissures, this will influence the results where parts of pulmonary fissures cannot be detected because the orientation can vary too much between patients [7]. Third, the proposed algorithm can effectively eliminate accessory fissures, which is a well-known challenge about pulmonary fissure detection [24–26, 45]. Fourth, the surface fitting model [46] is applied to complete the incomplete surface of pulmonary fissures. The approach ensures the local continuity of the estimated lobar boundaries. Last, the computerized scheme has a good performance in lobar fissure detection and lung lobe segmentation.Compared with manually segmented fissure references, our scheme achieved a medianF1-score of 0.8865 in the left lung and a median F1-score of 0.8920 in the right lung. The perfect performance of our computerized scheme can be exhibited from the lowest FNR index, which indicates the completeness of lobar fissure detection has been mended by our scheme. The major contribution is that the described scheme applied the priori knowledge of airways and pulmonary arteries to accurately identify fissure ROI and used a lower threshold to suppress clutters. While Xiao et al. [17] applied a sophisticated postprocessing pipeline and a volume sorting scheme for lobar fissure detection, a part of fissures was considered as clutters and removed. Using a different strategy, Peng et al. [18] used an appropriate threshold to suppress clutters and applied the volume sorting scheme to segment lobar fissures. However, the DoS and ODoS methods [25, 26] often faced a tradeoff between weak object detection and noise suppression.Our scheme is validated on 15 CT examinations. The average percentages of the segmented lung lobes in the lung lobe ground truth are 0.960, 0.989, 0.973, 0.920, and 0.985 for the left upper, left lower, right upper, right middle, and right lower lobes, respectively. The computerized scheme has a good performance in lung lobe segmentation. This is attributed to an ingeniously designed fusion of the surface fitting model, identified fissure ROI, and the oblique and horizontal fissure verification.In spite of the aforementioned benefits, the described scheme has many limitations. First, our scheme is heavily depending on airway segmentation, which is an arduous task and sensitive to the image quality. Second, the lung anatomy knowledge is employed to segment pulmonary fissures, which is time consuming. Third, our scheme is derived from the ODoS filter and inherits some drawbacks. For example, some deformed fissures cannot be detected by the computerized scheme due to their orientation distribution. Fourth, poor segmented airways and pulmonary arteries may cause parts of fissures to be undetected in some cases. Finally, the fitting accuracy of the estimated lobe boundaries is largely determined by verified fissures; a better segmented method of fissure detection and lung lobe segmentation is needed to pursue.


## 5. Conclusions

In this paper, we present a new lung lobe segmentation scheme. Motivated by the fact that pulmonary veins may cross pulmonary fissures, the priori knowledge of airways and pulmonary arteries is used to identify pulmonary fissure ROI. The second contribution of our scheme emphasizes on pulmonary fissure segmentation. Using an ingeniously designed framework, adhering clutters are reasonably removed. Finally, the surface fitting model, fissure line strength, and orientation are fused to complete the incomplete fissure surface for lobar boundary estimation. The proposed method is validated on 15 CT examinations. Experimental results showed that the described method has a good performance in pulmonary fissure detection and lung lobe segmentation.



---

*Source: 5588629-2021-04-20.xml*