Cystic Fibrosis (CF) is a chronic disease that mainly affects the respiratory and the digestive system. Its pathology in lungs tends to increase in severity over time . Computed Tomography (CT) Scans have shown to be an excellent tool for assessing different stages of the CF lung pathologies. They can be used in conjunction with standard lung function tests that provide functional information corresponding to the morphological information obtained from CT Scans. Quantifying CF features manually is time consuming and requires operation by trained individuals . In addition, the evaluations of the pathological stages often show a considerable variance among different individuals. Hence, an automated pathological feature detection of CF-CT Scans is desired as a decision support tool for doctors . In this paper, we focus on developing an automatic scheme to detect and evaluate Bronchiectasis in CF patients belonging to different stages of the disease.
Given chest CT slices of a CF patient, the proposed airway automatic detection methods consist of several sequential steps: preprocessing and lung segmentation as well as airway detection and their evaluation.
Grayscale reconstruction and the intensity thresholding techniques have been combined to boost the performance of airway detection in CT slice images of pathological lungs. This prototype scheme has been implemented with MATLAB R2015b.
2.1. Preprocessing and Lung Segmentation
Firstly, CT slices were cropped in order to focus more on the lung regions. Following this, preprocessing was started by applying low-pass filtering on the cropped CT image to eliminate noise. Contrast adjustment was carried out on the cropped CT image for mapping the pixel intensities of the airway walls to a higher value in the grayscale. Thresholding was done on this contrast adjusted image to produce a binary version of the preprocessed CT image, that was used for the intensity thresholding technique. For lung segmentation, the original CT slice image in Hounsfield Units (HU) was used to generate the binary lung mask. A histogram of the pixel intensities in HU was created to obtain the air threshold value (ATV), which was computed as an average of the HU values corresponding to the two largest peaks from the histogram . All pixels below the ATV were selected and subjected to grayscale morphological processing to produce the final lung mask. The resultant binary lung mask was multiplied by the preprocessed image to present the lung regions for detecting the airways.
2.2. Airways Detection Using Grayscale Reconstruction
In CT slices the enlargement of small or medium sized airways was used to denote the presence of Bronchiectasis in CF patients, which is analyzed as a vital feature for diagnosis. The airways were perceived to be the local minima of the histogram . The grayscale reconstruction technique was used to detect the small airways in the CT slice. Steps of the grayscale reconstruction technique have been listed below:
Step 1: Grayscale closing. This step marks airways of different sizes thereby, producing a marker image (denoted by J). The morphological closing step was applied on the preprocessed CT image. Explicitly, we state this step as:
J = I∙B0 (1)
where I is the preprocessed CT image, B0 is a 4-connected binary Structuring Element (SE) and the symbol ∙ denotes morphological closing operation. This step was performed for the initial identification of potential airways by increasing their gray values.
Step 2: Erosion. The marker image was eroded with the same sized SE B0. We represent this step by:
E = JɵB0 (2)
Eroded image E was obtained by performing the morphological erosion operation which is denoted by the symbol ɵ. This helped in removing small parts of the lung fields situated close to high attenuation pixels such as blood vessels, and that could be detected as non-airway regions, as shown in Figure 2, because of their intensity values being similar to that of the airways.
Step 3: Grayscale reconstruction, which is represented by:
G = max(E, I) (3)
where G is the grayscale reconstructed image, which compared each pixel value between E and I and displays the greater one. In the final grayscale reconstructed image, regions of local minima (airways) were replaced with an intensity value proportional to the difference between the maximum and minimum intensity values within the B-sized neighborhood of the airway lumen .
Step 4: Final reconstruction. This step calculates the difference image D, which was computed to highlight the regions of local minima that were detected in the previous steps. This step is represented by:
D = G − I (4)
The final grayscale reconstructed image appeared to have similar intensity values as that of I, except for the airways and other regions of local minima whose intensities were dramatically increased in the grayscale. The difference image D displayed an image with only the airways and small regions of local minima in the lung fields as clusters of pixels with higher gray values. Thresholding was applied on the difference image to retain only the airways.
Steps 1) to 4) are repeated with different sized SEs for detecting airways of different sizes. In our study, 4 SEs of different sizes were used. The first size of B0 enabled the detection of the smallest airway captured in the image. However, for getting the eroded image, the size of the SE was limited to having a total of 13 neighbors only.
2.3. Airways Detection Using Intensity Thresholding
Before we move onto explaining intensity thresholding, it is important to understand the terms parent object and children object. Parent objects are objects in a binary image, where the object pixels have a gray value of 1 and the background pixels have a gray value of 0. Children objects are objects present in the holes (consisting of pixels with gray value 0) of parent objects. Contrast adjustment done in section 2.1 was especially useful for carrying out intensity thresholding. The airways (children objects) which appeared to occupy holes in the lung regions (parent objects) were easily detected by this technique. Airways in the range of [20, 250] pixels were retained by a morphological opening step, that was implemented for improving the performance of intensity thresholding.
Finally, the airways detected by grayscale reconstruction and intensity thresholding were combined to generate a single image that produced a binary image of the detected airway pixels. This binary image was projected onto the lung regions extracted in Section 2.1. As an example, the overlapped image for one CT slice is shown in Figure 1 (right). The proposed airway detection scheme is summarized by flowcharts in Figure 2. The percentage of the airway pixels with respect to lung regions in each CT slices was calculated.
Three-dimensional CT Scans of four CF patients were employed for presenting the first results of the proposed automatic airway detection scheme. These CT images belong to patients with different stages of the disease. Detection of Bronchiectasis was done on two-dimensional CT slices between the first and the fifth intercostal spaces of each patient. Slices below the fifth intercostal space were excluded as they did not contain sufficient lung area for processing. Algorithm parameters were fixed for all patients, with a fully automatic generation of the air threshold value (ATV) which was used for lung segmentation (Section 2.1). While creating the lung mask for segmentation of the lung regions from original CT slices, the parameters for area opening and grayscale morphological closing were adjusted differently for slices above and below the carina. For slices above the carina, area opening was done to retain image objects containing 900 or
Figure 1. Original CT slice (left). Detected airways marked in pink in lung regions segmented out from the same CT slice (right).
Figure 2. Flowchart of the proposed automatic system for detection of enlarged airways.
more pixels. Whereas, for slices below the carina, area opening was done to retain image objects containing 600 or more pixels. Flat disk shaped structuring element having 5 neighbors was used for carrying out grayscale morphological closing on slices above the carina. On the contrary, for slices below the carina, a flat disk shaped structuring element having a total of 25 neighbors was used. All the parameters were kept unaltered for slices between the first and the fifth intercostal spaces during the detection of airways. For each slices, the area percentage of airways with respect to lung regions were calculated. Within each intercostal region, the mean of the airway percentages of each slice in the region has been presented in Figure 3.
This study was conducted to test the feasibility of the automatic system created for obtaining quantitative measurements of Bronchiectasis in CF patients. The Bronchiectasis stage was expressed as the percentage of visible airway counts in the corresponding lung regions of each CT slice.
The final plot fairly accurately depicts the severity of the disease in each patient. Patient 1, being at a very initial stage of the disease, shows very low percentages as compared to its counterparts. Patients 2 and 3 appear to be in similar stages of the disease as per the final plot. On the contrary, patient 4 displays the highest level of severity and therefore, belongs to a much higher stage of the disease.
However, small areas of non-airway regions were detected as airways as a result of a closing operation performed in the grayscale reconstruction technique for creating the marker image. Figure 4 shows these regions that could be described as small areas situated very close to high attenuation pixels (blood vessels) thereby, giving it a similar appearance to that of the airways which are low attenuation regions enclosed within
Figure 3. The mean airway percentage of each CT slice within each intercostals. This plot corresponds the 1st to the 5th intercostal regions of each patient.
Figure 4. Non-airway region (circled in red), having a similar appearance to that of an actual airway, falsely detected by the system.
airway walls, having intensities alike the blood vessels. The size of the structuring element had to be reduced for performing erosion; in order to remove these non-airway regions, while retaining the complete boundaries of the airways.
This paper describes the methods applied for creating a fully automatic system for the detection of Bronchiectasis using CT scans of CF patients. The area of pathological enlarged airway regions is presented as percentage of lung areas. The results obtained fairly match the severity of the disease stage of each patient.
This work is partially supported by the Federal Ministry of Education and Research (BMBF) under grant no. 03FH038I3 (MOSES).