Inspiratory Lung Expansion in Patients with Interstitial Lung Disease: CT Histogram Analyses

This study aimed to evaluate inspiratory lung expansion in patients with interstitial lung disease (ILD) using histogram analyses based on advanced image registration between inspiratory and expiratory CT scans. We included 16 female ILD patients and eight age- and sex-matched normal controls who underwent full-inspiratory and expiratory CT scans. The CT scans were sequentially aligned based on the surface, landmarks, and attenuation of the lung parenchyma. Histogram analyses were performed on the degree of lung expansion (DLE) of each pixel between the aligned images in x-, y-, z-axes, and 3-dimensionally (3D). The overall mean registration error was 1.9 mm between the CT scans. The DLE3D in ILD patients was smaller than in the controls (mean, 17.6 mm vs. 26.9 mm; p = 0.023), and less heterogeneous in terms of standard deviation, entropy, and uniformity (p < 0.05). These results were mainly due to similar results in the DLEZ of the lower lungs. A forced vital capacity tended to be weakly correlated with mean (r2 = 0.210; p = 0.074), and histogram parameters (r2 = 0.194–0.251; p = 0.048–0.100) of the DLE3D in the lower lung in ILD patients. Our findings indicate that reduced and less heterogeneous inspiratory lung expansion in ILD patients can be identified by using advanced accurate image registration.

Idiopathic pulmonary fibrosis (IPF), the most common form of interstitial lung disease (ILD), is a distinct type of chronic, progressive, fibrosing interstitial pneumonia of unknown cause, predominantly occurring in older adults [1][2][3] . The overall prognosis of IPF remained poor until the last few years, with a median survival ranging from 3 to 5 years and a 5-year survival rate ranging from 30% to 50% 4,5 . There had been no proved treatment other than lung transplantation for IPF. However, recent studies proved that nintedanib and pirfenidone could delay the progression of the fibrotic process in IPF, subsequently improving the overall prognosis. The use of these anti-fibrotic drugs is currently recommended for the treatment of IPF 2,6 .
Chest computed tomography (CT) is a key component of the management of ILD, as it is useful for the diagnosis, stratification of severity, assessment of treatment response, and prognosis prediction. Indeed, subjective evaluations of the degree of fibrosis on chest CT scans have been found to be useful for assessing severity and disease progression in IPF patients 7,8 . However, substantial variability existed in the subjective assessment across readers 9 . The objective quantification of severity and disease progression in ILD patients is becoming increasingly important.
Three-dimensional respiratory lung motion -in other words, the degree of respiratory lung inflation -is quantifiable by image registration between inspiratory and expiratory CT scans 10 . It was reported to be a useful option for assessing pulmonary functional loss in cases of chronic obstructive lung disease 11 . However, it has not yet been reported in patients with ILD, presumably due to the difficulty of accurate image registration of fibrotic lungs between inspiratory and expiratory CT scans. Recently introduced advanced image registration algorithms using attenuation-and landmark-based registration have shown to be allow more accurate registration in both animal and human study 10 . This algorithm may give a chance to overcome the difficulty of accurate image registration of fibrotic lungs during the respiration in patients with ILD. Thus, the purpose of our study was to evaluate inspiratory lung expansion in patients with ILD using histogram analyses based on advanced image registration of inspiratory and expiratory CT scans.

Results
Patient population. A total of 24 patients (mean age, 59.6 ± 10.2 years; age range, 41-76 years), encompassing 8 normal female controls (mean age, 58.8 ± 11.1 years; age range, 42-73 years) with normal CT scans, and 16 age-and sex-matched patients with ILD (mean age, 60.0 ± 10.1 years; age range, 41-76 years) who underwent paired full inspiratory-end expiratory thin-section chest CT scans were included. Male patients were excluded since there were no normal male subjects. Among the 8 controls, 5 had connective tissue disease and 3 were healthy subjects. Among the 16 patients, 10 were diagnosed as having connective tissue disease-related ILD (4 cases of systemic sclerosis, 3 cases of dermatomyositis, 2 cases of mixed connective tissue disease, and 1 case of rheumatoid arthritis), 4 as having IPF, and 2 as having nonspecific interstitial pneumonia (NSIP). Three patients had a pathologically proven diagnosis by means of surgical lung biopsy, and other diagnoses were made through a multi-disciplinary discussion between an experienced pulmonologist and radiologist.
Accuracy of image segmentation and registration. The Dice similarity coefficient (DSC) showed a good intra-observer agreement between the repeated lung segmentation (mean DSC, 0.972; range, 0.814-1.000). The detailed data of intra-observer agreement was shown in Supplementary Table 1.
The mean registration errors without image registration and with surface-based affine registration were 13.1 ± 6.7 mm and 15.6 ± 8.7 mm, and 8.2 ± 6.0 mm and 7.5 ± 6.5 mm in the right and left lung, respectively. When attenuation-based deformation was additionally applied without landmark-based registration, the error was 5.7 ± 4.7 mm and 7.7 ± 7.4 mm, respectively. The errors further decreased to 4.1 ± 2.9 mm and 4.8 ± 2.9 mm, respectively, by applying a bronchial landmark, and the application of both bronchial and peripheral vascular landmarks finally resulted in errors of 1.9 ± 1.2 mm and 1.9 ± 0.9 mm in the right and left lung, respectively (Supplementary Table 2). In terms of lung volume, the mean volume difference between inspiratory and expiratory lung volume was 484.2 ± 206.9 cm 3 and 421.4 ± 229.5 cm 3 in the right and left lung before registration, respectively. The volume difference was sequentially decreased with applying affine, landmark-based, and deformable registration, and the volume difference between inspiratory and registered expiratory lung volume became 91.6 ± 81.3 cm 3 and 99.9 ± 93.6 cm 3 in the right and left lung after all registrations were applied, respectively (Supplementary Table 3).
Inspiratory lung expansion analysis. The degree of fibrosis of the whole lung on visual analysis was 17.2% ± 14.8% in ILD patients. The degrees of fibrosis in the upper, middle, and lower lungs were 8.1% ± 13.9%, 10.3% ± 13.1%, and 33.1% ± 21.6%, respectively. The degree of fibrosis in the lower lung was significantly higher than in the upper and middle lungs in ILD patients (p = 0.001).
Degree of lung expansion (DLE) of each pixel between the full-inspiratory and expiratory CT scans was measured in horizontal axis (x-axis; DLE x ), ventrodorsal axis (y-axis; DLE y ), craniocaudal axis (z-axis; DLE z ), and 3-dimensionally (3D; DLE 3D ) (Reply to the grammar issue 4). Regarding whole lung expansion, the mean DLE 3D of ILD patients was significantly smaller than those of normal controls (17.6 mm vs. 26.9 mm, p = 0.023) ( Table 1 and Fig. 1). The histogram of ILD patients had significantly higher kurtosis (p = 0.016) than that of normal controls and tended to be more positively skewed (p = 0.235) (Fig. 2). Statistical significance was found in the 5th to 95th percentiles of DLE 3D (all p < 0.05) ( Fig. 3 and Supplementary Table 4). The DLE 3D values were more homogeneous in ILD patients than in normal controls, as shown by a smaller standard deviation (p = 0.015), less entropy (p = 0.018), and greater uniformity (p = 0.015) ( Table 1). Among the three axes, similar results were found in the histogram parameters of DLE z , with statistically significant differences between the two groups for mean (p = 0.011), standard deviation (p = 0.021), skewness (p = 0.012), kurtosis (p = 0.003), entropy (p = 0.013), uniformity (p = 0.005), and 20th to 95th percentiles. A similar trend was also observed in the y-axis, but none of the histogram parameters showed statistically significant differences (Table 1).
When analyzed in the upper and lower lungs, the aforementioned results of the whole lung were similar to those of the lower lung except for kurtosis ( Table 2). The mean (p = 0.007 and 0.004) and 5th to 95th percentiles of the DLE 3D and DLE z in the lower lung were significantly lower in ILD patients (Supplementary Table 5 and Supplementary Figure 1). The DLE 3D and DLE z of the lower lung were more homogeneous in ILD patients, as shown by a smaller standard deviation (p = 0.028 and 0.048), less entropy (p = 0.020 and 0.019), and greater uniformity (p = 0.016 and 0.023). Skewness of the DLE z was also significantly higher in ILD patients (p = 0.003) ( Supplementary Figures 2 and 3). In the upper lung, however, the parameters were mostly similar between the two groups.
The results of mean and standard deviation of DLE were similar after the standardization of DLE to the reference value of vital capacity (VC). In the analysis of the whole lung, the standardized mean of DLE 3D was significantly smaller in the ILD patients (6.44 mm/L vs. 9.90 mm/L, p = 0.033) as well as DLE z (p = 0.021). The standardized standard deviation of DLE 3D and DLE z was also significantly smaller in the ILD patients (p = 0.021 and 0.047) ( Table 3). The results of the lower lung was similar to those of the whole lung, with smaller standardized mean of DLE 3D , DLE y , and DLE z (p = 0.010, 0.047 and 0.008) and standardized standard deviation of DLE 3D and DLE x (p = 0.050 and 0.031). However, there was no significant difference in the standardized parameters of DLE between the two groups in the upper lung (Table 3).

Correlation between inspiratory lung expansion and pulmonary function test (PFT).
Among the histogram parameters, the mean (r 2 = 0.210; p = 0.074), entropy (r 2 = 0.182; p = 0.100) and 5th to 95th percentiles (r 2 = 0.194 to 0.251; p-value, 0.048 to 0.093) of the DLE 3D in the lower lung tended to show a weak positive correlation with forced vital capacity (FVC), and the uniformity (r 2 = 0.210; p = 0.074) of DLE 3D in the lower lung showed a weak negative correlation (Table 4).

Discussion
Our quantitative 3-dimensional CT analysis revealed reduced and less heterogeneous inspiratory lung expansion in ILD patients, mainly attributable to reduced craniocaudal expansion of the lower lung. These results were obtained by advanced accurate registration of paired inspiratory and expiratory chest CT images based on the surface, landmarks, and attenuation of lung parenchyma 12,13 , resulting in a mean registration error of less than 2 mm.
Little is known about the pixel-level distribution of inspiratory lung expansion in normal subjects and ILD patients. In normal subjects, the distribution of DLE 3D mimics an inverted U-shape with relatively preserved symmetry, and the distribution of DLE z seems to be the main contributor to the shape of DLE 3D (Fig. 2). Based on a separate analysis of the upper and lower lungs, the distribution of DLE 3D and DLE z in the whole lung was mainly determined by the DLE 3D and DLE z in the lower lungs. In ILD patients, although the distribution of DLE in all 3-axes tended to be positively skewed and asymmetric, the decrease of DLE z was the most pronounced (in comparison to DLE x and DLE y ), especially in the lower lungs (Fig. 2). Those differences between normal subjects and ILD patients were supported by the statistically significant decrease of the DLE 3D and DLE z in ILD patients.
Our results also revealed less heterogeneous 3D lung expansions in ILD patients, which were especially pronounced for the lower lungs and in the craniocaudal direction. At inspiration, the lung parenchyma typically expands more in the peripheral and basal areas than in the central and upper areas 14,15 . Accordingly, the degree of peripheral and basal lung expansion is the main determinant of the degree of heterogeneity of overall lung expansion. In our ILD patients, fibrotic changes predominantly involved the basal and peripheral lungs, as in most IPF and NSIP patients 16,17 . The less heterogeneous lung expansion of ILD patients probably resulted from decreased expansion of the basal and peripheral lungs due to fibrotic changes (Fig. 1).
The reduced lung expansion in ILD patients presumably originates from reduced lung compliance due to fibrosis 18,19 . The elastic compliance of the lung parenchyma decreased in ILD patients, and the degree to which it was reduced was correlated with the degree of fibrosis due to increased resistance 20 . The increased resistance was mainly due to the abnormality of the lung parenchyma, rather than altered chest wall mechanics 19 . Although we were not able to correlate the degree of lung expansion with the degree of elastic resistance, our analysis of inspiratory lung expansion may be a noninvasive indirect imaging biomarker to assess lung compliance in ILD patients, since the elastic resistance limits the inspiratory inflation of the lung parenchyma.
The histogram parameters of lung expansion in the lower lungs were marginally correlated with FVC, which is an important predictive marker of survival in IPF patients 21,22 . Furthermore, there is no established tool to quantify inspiratory inflation of the basal lungs, where ILD is often involved. Although diaphragmatic motion can be evaluated as a proxy for basal lung movement, the compensatory lengthening and curving of the diaphragm hindered its use as an appropriate evaluation in IPF patients 23 ; in contrast, our histogram analysis showed a distinct differentiation of lower lung movements between normal subjects and ILD patients. This direct quantification of lower lung movements may also be used for monitoring respiratory muscle strengthening in pulmonary rehabilitation, the importance of which is increasing in ILD patients 24 .
Our study has several limitations. First, our patient population was relatively small and exclusively consisted of women due to the paucity of normal male controls who underwent paired inspiratory and expiratory chest CT scans. Nevertheless, we expect that similar observations would be obtained in male ILD patients if there were a lower lobar predilection of fibrotic changes, as was the case in our study population. Second, it is well recognized that pulmonary function can vary with body habitus, as well as age and sex 25 . Thus, the degree of lung expansion may also be influenced by age, sex, and body habitus. Although we matched ILD patients with normal controls based on age and sex and we standardized DLE values with the reference lung volume, this may not be sufficient to precisely reflect the differences lung expansion due to age, sex, and body habitus. Therefore, further study may be required to validate our results with a sufficient number of patients, along with adjustment for age, sex, and body habitus. Third, full inspiration and end expiration might not have been sufficiently achieved during the CT scanning due to the patients' non-compliance, although we gave detailed instructions before scanning. This might have led to an underestimation of the lung expansion between full inspiration and end expiration. Fourth, a minor portion of the destroyed peripheral lung parenchyma in ILD patients might not have been segmented due to severe fibrosis, although we adjusted the attenuation threshold values for lung segmentation in ILD patients to maximize the inclusion of fibrotic lung areas. Fifth, although we assessed regional lung expansions, our results were not based on the pulmonary lobar anatomy. However, lobar segmentation was not practically achievable in some cases due to an incomplete fissure or severe destructive fibrosis hampering a proper identification of the fissure.
In conclusion, reduced and less heterogeneous 3D inspiratory lung expansions were found in ILD patients using histogram analyses of DLE based on the advanced accurate image registration of paired inspiratory and expiratory CT scans. These results were mainly attributable to reduced craniocaudal lung expansion in the lower lung. The quantitative analysis of DLE at inspiration on chest CT scans may be a unique potential noninvasive imaging biomarker in ILD patients for indirectly evaluating elastic lung compliance and directly assessing fibrotic basal lung movements and respiratory muscle function.

Materials and Methods
Our institutional review board (Seoul national university hospital institutional review board) approved this retrospective study, and informed consent was waived.  Patient selection. A total of 392 subjects who underwent paired full-inspiratory and expiratory chest CT scans from January 2013 to February 2016 were selected from the electronic medical records of our institution. Among the subjects, 8 female subjects had normal chest CT scan, and 384 had abnormal chest CT scan. All female subjects with normal chest CT scan were selected as normal female controls. As there were no normal male subjects who underwent the paired chest CT scans, we excluded male patients. Of 114 female subjects with abnormal chest CT scan, 51 subjects had interstitial lung abnormalities. 16 age interval-matched females were randomly selected within the subjects with interstitial lung abnormalities. Finally, 24 subjects (8 controls, 16 patients) were included in this study (Fig. 4).
CT acquisition. All CT examinations were performed using multidetector CT machines with two 64-channel scanners (n = 4; Brilliance 64, Philips Medical Systems, the Netherlands) (n = 20; Ingenuity, Philips Medical Systems, the Netherlands). The scanning parameters were as follows: detector configurations, 64 × 0.625 mm; tube voltage, 120 kVp; tube current, reference mAs of 200 with automatic exposure control; pitch, 1.015; reconstruction kernel; YC0; slice thickness, 1 mm; and reconstruction interval, 1 mm. We used the standardized CT protocol for these noncontrast chest CT scans. All patients were trained in breathing prior to CT scan. Patients were instructed to perform maximal inflation so that the lung volume could be near total lung capacity (TLC) at full inspiration, and relaxed exhalation to allow the lung volume to become functional residual capacity (FRC) at end expiration as in the previous study by Kauczor H.U. et al. 26 . After providing detailed instructions on breathing, thin-section CT images were obtained in the supine position at full-inspiration and end expiration.   Image segmentation and registration. The paired inspiratory and expiratory CT images were aligned using in-house software (Fig. 5). To fully extract the regions of the lung containing fibrosis caused by ILD, lung segmentation was performed by one researcher (J.J.) using an appropriate threshold between −400 and −200 Hounsfield units. To check the intra-observer variability of the segmentation, the researcher repeated the segmentation after several months. The segmented lung in the expiratory CT image (source image) was then globally aligned with the segmented lung in the inspiratory CT image (target image) by inflating the expiratory lung using surface-based affine registration 12 . To improve the alignment of the internal structures of the lung parenchyma, landmark-based registration using a thin-plate spline warping was performed: one radiologist chose 31-34 landmarks of the bronchial tree and pulmonary vessel per lung 13 . The landmarks were marked carefully by placing small dots in the side-by-side display of 1-mm-thick inspiratory and expiratory CT images to ensure that identical locations were annotated in the CT images. The branching points of the lobar and segmental bronchi served as the bronchial landmarks, while pulmonary vessels located within 2 cm from the pleura were used for the vascular landmarks. The landmarks were evenly distributed in the upper to lower lungs (Fig. 6). Finally, the registered lungs were locally aligned using attenuation-based demons deformable registration 12 . The main idea of the demons deformable registration is to position 'demons' at certain voxels of the target image and to make the source image diffuse to the target image using the estimated force based on the difference in attenuation between the target and source images.

Accuracy evaluation of image registration.
To quantitatively evaluate the accuracy of image registration, we calculated the registration error for each patient as the mean distance (mm) of the corresponding bronchial and vascular landmarks between the inspiratory CT and registered expiratory CT images. The registration error was assessed at each step of applying three registration methods (surface-based affine registration, landmark-based registration, and attenuation-based deformable registration) in order to identify the effect of combining the registration approaches, especially for landmark-based registration.
Inspiratory lung expansion analysis. After checking the accuracy of image registration between the inspiratory and expiratory CT scans, we measured DLE of each pixel between the two CT scans in the following axes: horizontal axis (x-axis; DLE x ), ventrodorsal axis (y-axis; DLE y ), craniocaudal axis (z-axis; DLE z ), and 3-dimensionally (3D; DLE 3D ). Histogram analyses were then performed for the right and left lungs in all axes. Each lung was divided into the upper and lower lungs at the branching point of the lower lobar bronchus from either the right bronchus intermedius or the left main bronchus. We performed the analyses for the whole lung and for the upper and lower lungs separately.
The following histogram parameters were calculated using the probability distribution of the DLEs in each axis: mean, standard deviation, skewness, kurtosis, and the 5th, 10th, 20th, …, 90th and 95th percentiles, entropy and uniformity.
The mean (M) is the average value, and is calculated as equation (1): where N is the number of discrete bins of the histogram and P(x i ) is the probability of bin x i . The standard deviation (SD), which is a measure of contrast, is calculated as equation (2); it describes the spread in the data, so a high contrast in DLE will have high variance and a low contrast in DLE will have low variance.
Skewness (SK), which is a measure of asymmetry about the mean in the DLE distribution, is calculated as equation (3); it is positive if the tail of the histogram spreads to the right, and negative if the tail of the histogram spreads to the left. Kurtosis (K), which is a measure of the relative flatness of the histogram, is calculated as equation (4); a high value of kurtosis indicates that the peak of the distribution is sharp and the tail is longer and fat. A low value of kurtosis indicates that the peak of the distribution is rounded and the tail is shorter and thinner. Entropy (E), which is a measure of the variability of the DLE distribution, is calculated as equation (5), and increases as the variability of P(x i ) increases.
i 1 N i i Uniformity (U), which is a measure of the homogeneity of the DLE distribution, is defined as equation (6), and is maximized when all P(x i ) are equal.
Since the DLE increases with increasing lung volume, we calculated standardized mean and standard deviation of DLE in each axis by dividing the values by the reference value of VC to adjust the differences of DLE according to the reference lung volume. The reference lung volume was calculated according to the previous study 27 , and the reference value of VC was calculated by subtracting the reference value of RV from TLC. The standardized parameters of DLE were separately compared between the ILD patients and the normal controls.
CT visual analysis. Visual analysis was performed by a resident (J.H.P.) from the Department of Radiology under the supervision of an experienced thoracic radiologist (S.H.Y., 12 years of thoracic CT experience) to assess the degree of pulmonary fibrosis. The mean extent of fibrotic changes, including reticular opacity and honeycombing, was scored at 5% intervals in the upper, middle, and lower zones of each lung using the inspiratory CT images. Each zone was defined according to the previous study by Best et al. 8 . The degree of pulmonary fibrosis was summarized as the mean value of the fibrosis score of the three zones.
Correlation with pulmonary function tests. Sixteen of the 24 subjects underwent a PFT within four weeks of chest CT scanning. The mean interval between the PFT and CT scanning was 4.4 ± 7.9 days (range, 0-26 days). We obtained FVC, diffusing capacity of carbon monoxide (DLCO) and alveolar volume (VA) values in the PFT. Results were expressed as both absolute value (FVC and DLCO) and percentage of predicted performance (%FVC and %DLCO).

Figure 6.
A representative map of landmarks of the bronchial tree and pulmonary vessels. In this patient, a total of 67 landmarks (33 in right lung, 33 in left lung, and on at carina) were placed in a side-by-side display of 1-mm-thick inspiratory and expiratory CT images. The landmarks in the inspiratory (red dots) and expiratory CT (blue dots) are paired to represent the same location.
SCIeNtIfIC REPORTs | (2018) 8:15265 | DOI:10.1038/s41598-018-33638-x Statistical analysis. The Dice similarity coefficient was used to calculate the intra-observer agreement of the segmentation process 28 . Nonparametric repeated-measures analysis of variance was used to compare the histogram parameters of DLE in each axis between the ILD patients and the normal controls. The Pearson correlation coefficient or Spearman correlation coefficient by rank was used to assess the relationship between PFT results and the histogram features of DLE 3D in the lower lung depending on whether or not they showed a normal distribution. A two-tailed p-value of less than 0.05 was considered to indicate a significant difference. Statistical analyses were carried out using SAS software (version 9.4; SAS Institute Inc., Cary, NC, USA) and SPSS software (version 23; IBM Corp., Armonk, NY, USA).