Quantitative analysis of computed tomography of the lungs in patients with lymphangioleiomyomatosis treated with sirolimus.

OBJECTIVES: We aimed to study sirolimus-related lung parenchymal changes by quantitative analysis of computed tomography (CT) of the lungs in patients with lymphangioleiomyomatosis (LAM).
METHODS: We studied 20 participants from the Multicenter Lymphangioleiomyomatosis Sirolimus Trial for Safety study, who had undergone both thin-section CT scans and pulmonary function tests at baseline, 12, and 24 months. Quantitative CT parameters such as CT-derived total lung capacity, percentage of low attenuation area (LAA%), lung density histogram, fractal property of low attenuation area, and airway dimensions were analyzed, and correlations were conducted between the longitudinal change in each quantitative CT measurement and changes in pulmonary function were examined. Among 20 participants, pre-trial (n = 8) and post-trial (n = 16) CT data were also analyzed to deduce pathophysiologic implications of the serial changes in CT parameters during trial periods.
RESULTS: FEV1 significantly increased from baseline to 24 months (slope 3.71 ± 1.50 ml/month) whereas FVC didn't during sirolimus therapy. Strikingly, LAA%, and skewness and kurtosis of density histogram significantly increased from baseline to 24 months, while mean and mode CT values significantly decreased from baseline to 24 months. Statistically significant positive correlations were found between ΔFEV1 and Δskewness (r = 0.465, p = 0.045). Taking the changes in lung density during pre-trial period into consideration, sirolimus decreases the area of -800 to -750 Housefield unit (HU) density and inhibits the decrease of -950 to -800 HU area during treatment, then producing the increased LAA% during the trial and post-trial periods. Given few sirolimus-related changes in airway dimensions, possible changes in lung mechanics may have contributed to increased FEV1.
CONCLUSION: Our study suggests that the lung density histogram parameters, kurtosis, and skewness, may be useful as indicators of the efficacy of sirolimus.

Introduction

Lymphangioleiomyomatosis (LAM) is a rare, progressive neoplastic lung disease that primarily affects women. LAM is characterized by proliferation of abnormal smooth muscle-like cells (LAM cells) in the lungs and along lymphatic vessels, including the lymph nodes and thoracic ducts. LAM presents with cystic lung destruction, lymphatic involvement (e.g., chylous pleural effusions, lymphangioleiomyomas), and renal angiomyolipomas (AMLs) [1, 2]. Recent clinical trials, the Multicenter International Lymphangioleiomyomatosis Efficacy and Safety of Sirolimus (MILES) trial and the Multicenter Lymphangioleiomyomatosis Sirolimus Trial for Safety (MLSTS) successfully demonstrated that the mTOR inhibitor, sirolimus, stabilized lung function decline and improved the quality of life in adult females with LAM [3, 4]. Furthermore, serum vascular endothelial growth factor D (VEGF-D) was markedly reduced by sirolimus [5], but tended to increase again when the drug was withdrawn [5].

Previous studies have shown the utility of the quantitative analyses of chest computed tomography (CT) images for evaluating disease severity in LAM patients [6, 7]. Using these methods, Yao et al. reported that sirolimus therapy slowed down the increase in cystic lesions in LAM [8]. Argula et al. also demonstrated that sirolimus therapy attenuated progressive gas trapping in the disease [9]. Although both studies evaluated cysts in the condition by CT lung imaging, no other quantitative CT measurements were assessed, such as CT-derived total lung capacity, low attenuation area (LAA), lung density histogram, fractal property, and airway dimensions, all of which are suitable for evaluating the efficacy of sirolimus therapy because they reflect pathophysiological aspects of LAM.

The purpose of this study is to examine changes in quantitative CT measurements during sirolimus treatment in a subset of MLSTS participants to investigate whether the CT variables more closely reflect the effects of sirolimus in the patients with LAM.

Materials and methods

Study population

We used the data of 25 LAM patients who participated in the MLSTS [4] at our hospital. The MLSTS was a phase II, multicenter, open trial of sirolimus of 2 years duration in order to confirm the safety profile in Japanese patients with LAM, approved by our institutional review board (JIRB2012-003 and 18–193), and listed at Japan Medical Association, Center for Clinical Trials (www.jmacct.med.or.jp/en/; JMA-IIA00096). All participants received sirolimus for 2 years at the doses adjusted to maintain a trough blood level of 5–15 ng/ml. The participants had pulmonary function tests at baseline and every 6 months thereafter. They underwent thin-section CT scans at baseline, 12, and 24 months. LAM was diagnosed from tissue biopsies or characteristic clinical pictures (recurrent pneumothorax and/or chylous pleural effusion) and CT findings (diffusely scattered thin-walled pulmonary cysts) [10]. We excluded 4 patients with pleural effusion and 1 patient who had LAM-associated protein-losing enteropathy but very mild lung disease from the analysis, and 20 patients were finally included.

Thin-section CT techniques

All patients underwent chest thin-section CT with a 64-detector row CT scanner (Aquilion 64 scanner; Toshiba Medical, Tokyo, Japan) using 0.5-mm and 2-mm slice thickness. Scanning parameters were 120 kVp, 150 mAs and a field of view of 320 mm or 340 mm. No contrast media were used. CT images were obtained in the supine position to assess total lung capacity. During the scan, the patients held their breath after a deep inspiration in the supine position. Each CT image was composed of a 512 × 512 matrix of numeric data (CT numbers) in Hounsfield units (HU) reconstructed using a kernel of FC03 and FC85 [10].

Quantitative analysis of CT images

All images were transferred to a personal computer for postprocessing and image analysis. First, for the correction of calibration due to X-ray tube ageing, we measured air density sampled from the intrathoracic tracheal air in the raw CT images of each subject using and semiautomatic quantitative image-processing software (Image J version 1.48, available at http://rsb.info.nih.gov/ij/) [11, 12]. We confirmed that these air densities showed no statictically significant changes during the observations (Table S1).

CT-derived total lung capacity. The analysis for CT-derived total lung capacity (CT-TLC) was carried out with CT images of 2-mm slice thickness (FC85) using free open-source software (Airway Inspector, Brigham and Women's Hospital, Boston, MA, USA) [www.airwayinspector.org], as conducted previously [13]. The software automatically segmented the lung parenchyma from the chest wall and the hilum, and measured CT-TLC.

Analysis of low attenuation area, lung density histogram, and fractal property. The analyses for LAA, lung density histogram, and fractal property were carried out with CT images of 2-mm slice thickness (FC 85) using ImageJ. Three slices from each patient were analyzed (the upper slice taken 1 cm above the upper margin of the aortic arch; the middle slice taken 1 cm below the carina; and the lower slice taken 1 cm above the top of the diaphragm), and a mean score of all images was considered as a representative value for each patient. We defined lung fields as areas with CT numbers less than -200 HU, whereas the cut-off level between LAA was set at -960HU [14]. The percentage of low attenuation area (LAA%) was determined as the percentage of LAA per total lung area.

A density histogram of lung area in each CT image was generated, and mean and mode CT values, kurtosis, and skewness were measured from the histogram (Figure 1). Kurtosis and skewness represent the distortion and the disparity deviation of a histogram when compared with a normal distribution, and these indicators have been reported to correlate with changes in lung structure [15]. To investigate longitudinal changes in lung density in detail, the number of pixels in the lung field of -700 HU or less were calculated at every 50 HU interval and the percentage of each 50HU density area occupying in the lung field of -700 HU or less (pixel%) was evaluated for each image.

Figure 1

The representative lung density histogram of the middle lung field in a 43-year-old patient. The lung density histogram showed the number of pixels at every 50 HU interval. Black and gray bars represent the values of baseline and 48 months, respectively. LAA%, kurtosis, and skewness increased from baseline to 48 months. Mean and mode CT values decreased from baseline to 48 months.

Fractals were self-similar structures characterized by power-law functions and noninteger dimensions (fractal dimension) [16]. The details regarding the fractal properties of the distributions of LAA sizes are described in the online supplement.

Airway dimensions. The analyses for airway dimensions were carried out with CT images of 0.5-mm slice thickness (FC03) using Airway Inspector [13, 17]. The inner and outer airway walls were identified by using both the full width at half maximum and phase congruency edge-detection methods [18]. We analyzed the cross-sectional airway parameters; luminal area (Ai), wall area (WA), and the percentage of WA to the total area of the airway (WA%). The details are described in the online supplement.

Pulmonary function tests

Forced vital capacity (FVC) and forced expiratory volume in 1 s (FEV1) after inhalation of short-acting β2-agonist were measured using an electronic spirometer (Autospiro; Minato Medical Science, Osaka, Japan) at baseline, 6, 12, and 24 months after sirolimus treatment.

Statistical analysis

All statistical analyses were performed using IBM SPSS version 19 software (IBM JAPAN, Tokyo, Japan). A linear mixed-effects model was used to estimate longitudinal changes in each parameter. The Wilcoxon signed-rank test for paired samples was used for within-group comparisons. Relationships between data obtained from pulmonary function tests and measured CT parameters were assessed by the Spearman's rank correlation test. P values <0.05 were considered statistically significant.

Results

Characteristics of the study population

Results are presented as mean ± SEM (standard error of mean) unless otherwise specified. The baseline characteristics of 20 participants are shown in Tables 1 and 2. Mean age was 41.9 years old (range = 28–61), FVC was 2.82 L (92.1% of the predicted value), and FEV1 was 1.51 L (57.5% of the predicted value). At baseline, 11 patients were treated with gonadotropin-releasing hormone analogue (GnRH) therapy, and 1 patient received supplemental oxygen therapy. 2 patients withdrew from the trial after completion of examinations at 12 months, and leaving 18 patients for analysis at 24 months (Table 2). After completion of the MSLST study, 16 patients continued to visit and received chest CTs at 48 months, using the study protocol and scanner (Table 2).

Table 1

Baseline clinical characteristics of the study patients (n = 20).

Characteristics		Value (%)
Age	(yr)	41.9 ± 1.61
FVC	(L)	2.82 ± 0.10
FVC %pred	(%)	92.1 ± 2.5
FEV1	(L)	1.51 ± 0.12
FEV1 %pred	(%)	57.5 ± 4.0
FEV1/FVC x100	(%)	52.5 ± 3.0
Supplemental oxygen therapy	n	1 (5)
GnRH therapy	n	11 (55)
Lymphangioleiomyomas	n	7 (35)
Renal angiomyolipoma	n	3 (15)
Liver angiomyolipoma	n	2 (10)

FEV1 = forced expiratory volume in 1 s; FVC = forced vital capacity; GnRH = gonadotropin releasing hormone analogue; %pred = percentage of the predicted value.

Table 2

Longitudinal changes of pulmonary function and CT parameters in all patients.

		Baseline	6 months	12 months	24 months	48 months
No. of patients		20	20	20	18	16
GnRH therapy		11	10	10	10	9
Lung Function
FVC	(L)	2.82 ± 0.10	2.86 ± 0.11	2.86 ± 0.11	2.79 ± 0.12	NA
FEV1	(L)	1.51 ± 0.12	1.61 ± 0.13	1.60 ± 0.13	1.57 ± 0.14†	NA
CT parameters
CT-TLC	(mL)	4584 ± 126	4555 ± 102	4536 ± 125	4525 ± 109	4478 ± 124
LAA%	(%)	29.8 ± 2.1	31.2 ± 1.7	31.5 ± 2.0	33.4 ± 2.0‡	32.9 ± 1.7‡
Histogram analysis
Mean CT value	(HU)	-880 ± 4.2	-885 ± 3.1	-884 ± 4.5	-890 ± 3.7‡	-887 ± 3.6‡
Mode CT value	(HU)	-937 ± 8.2	-941 ± 5.5	-941 ± 7.7	-947 ± 7.0	-946 ± 6.8†
Kurtosis		5.78 ± 0.25	6.24 ± 0.28	6.16 ± 0.31	6.44 ± 0.31†	6.35 ± 0.31*
Skewness		2.06 ± 0.04	2.16 ± 0.05	2.14 ± 0.05	2.20 ± 0.05‡	2.21 ± 0.06‡
Fractal analysis
r2 value		0.950 ± 0.006	0.953 ± 0.005	0.958 ± 0.004	0.966 ± 0.003†	0.968 ± 0.004‡
No. of patients with r2 ≥ 0.95		13/20 (65.5%)	10/20 (50.0%)	13/20 (65.5%)	15/18 (83.3%)	15/16 (93.8%)
D-value		1.030 ± 0.045	0.994 ± 0.038	0.996 ± 0.039	0.943 ± 0.034†	0.885 ± 0.025‡
D-value in patients with r2 ≥ 0.95		0.965 ± 0.051	0.965 ± 0.045	0.912 ± 0.041	0.928 ± 0.038†	0.876 ± 0.025‡
Airway dimensions
RB1 Ai	(mm2)	18.7 ± 1.2	18.6 ± 1.2	18.5 ± 1.2	18.4 ± 1.3	17.4 ± 1.2
WA%	(%)	51.9 ± 1.1	52.1 ± 0.9	52.6 ± 0.8	52.8 ± 0.89	53.7 ± 1.1†
RB10 Ai	(mm2)	23.6 ± 1.6	24.7 ± 1.9	23.9 ± 1.7	23.2 ± 1.8	22.2 ± 2.2
WA%	(%)	49.0 ± 1.0	48.3 ± 1.0	48.6 ± 0.9	49.2 ± 1.0	50.4 ± 1.9

Ai = luminal area; CT-TLC = computed tomography-derived total lung capacity; FEV1 = forced expiratory volume in 1 s; FVC = forced vital capacity; GnRH = gonadotropin releasing hormone analogue; LAA% = percentage of low-attenuation area; No = number; RB1 = the right upper lobe apical bronchus; RB10 = the right lower lobe posterior basilar bronchus; WA% = percentage of wall area to total airway area.

*p < 0.05, †p < 0.01 and ‡p < 0.001 were calculated by comparing the linear mixed-effects model with a slope of zero.

Longitudinal changes of pulmonary function

FEV1 significantly increased from baseline to 24 months whereas FVC showed no significant change (Table 2). The slope of longitudinal change of FEV1 from baseline to 24 months was 3.71 ± 1.50 ml per month (Table S2). TLC was evaluated as CT-TLC, but no significant change was found from baseline to 24 months.

Quantitative analysis of CT images and their longitudinal changes

We found significant changes in various parameters determined by quantitative measurement of CT images (Table 2 and their slopes of longitudinal changes in Table S2). LAA% significantly increased at 24 months as well as after the MLSTS (i.e., at 48 months).

Histogram analysis of lung parenchymal CT density, using 4 different parameters (mean and mode CT values, kurtosis, and skewness) showed the most remarkable changes during the MLSTS and post-MLSTS periods. The decrease of both mean and mode CT values together with the increase of both kurtosis and skewness indicates that the lung density histogram had shifted toward air density with a sharper peak and more asymmetrical density distribution, as illustrated in Figure 1. To identify what lung density category showed the most prominent changes during sirolimus therapy, we investigated the pixel% at every 50 HU in the lung field of -700 HU or less from baseline to 48 months. As shown in Figure 2A, the pixel% of -1000 to -950 HU significantly increased from baseline to 24 and 48 months (p = 0.0001 and 0.00001, respectively), whereas the pixel% of -850 to -800 HU significantly decreased from baseline to 24 and 48 months (p = 0.0008 and 0.0002, respectively) (details in the changes of pixel% provided in the online supplement Table S3).

Figure 2

The changes in the lung density profile prior to and during sirolimus therapy. (A) The differences of pixel% from baseline during sirolimus therapy in all patients. Pixel% of -1000 to -950 HU significantly increased from baseline to 24 and 48 months, whereas pixel% of -850 to -800 HU significantly decreased from baseline to 24 and 48 months. *p < 0.05 and ∗∗p < 0.01 were calculated by comparing the linear mixed-effects model with a slope of zero. pp = percent point. (B) The comparison of changes in between from pre-treatment to baseline (closed column) and from baseline to 24 months (open column) in the 8 patients whose pre-trial data were available for analysis. Pixel% of -900 to -850 HU decreased from pre-treatment to baseline, but was stable from baseline to 24 months; the difference in the trends were statistically significant. Additionally, pixel% of -800 to -750 HU significantly increased from pre-treatment to baseline whereas decreased during sirolimus therapy. †p < 0.05 generated by the Wilcoxon signed-rank test. pp = percent point.

Fractal analysis of LAA (i.e., pulmonary cysts) revealed that the r2 value significantly increased from baseline to 24 and 48 months together with the increased number of patients whose r2 value was ≥0.95. On the other hand, D-values of the entire cohort as well as that of patients with r2 ≥ 0.95 significantly increased from baseline to 24 and 48 months. Serial changes of these fractal parameters suggest that terminal airspace geometry in LAM lungs became simpler during sirolimus therapy.

When airway dimensions were analyzed, only WA% of the right apical bronchus (RB1) showed a statistically significant increase at 48 months, but the other measurements were not significant.

Comparisons of longitudinal changes of pulmonary function and CT quantitative parameters between pre-trial and trial periods

8 of 20 patients had undergone both pulmonary function tests and CT scans before participating in the MLSTS (pre-trial data), using the same protocol and equipment that was employed for the MLSTS. Having pre-trial data for only a limited number of patients, we compared the longitudinal changes of various parameters between pre-trial and trial periods (Table 3 and Figures 2B and 3). Both pulmonary function tests and CT scan at pre-treatment were performed 24 ± 2.4 months (mean ± SD) before the initiation of sirolimus. As shown in Table 3 and Figure 3, FEV1 tended to decrease in the pre-trial period whereas increased during the trial, and the difference in trends between pre-trial and trial periods were statistically significant. LAA% consistently increased in the both periods. However, both kurtosis and skewness of the lung density histogram showed a similar trend to that of FEV1, a decrease in the pre-trial period but an increase during the trial. Consistent with the trends in both kurtosis and skewness, we found significant differences in the trend of pixel% between pre-trial and trial periods (Figure 2B). In the pre-trial period, the pixel% of both -1050 to -950 HU and -800 to -700 HU area increased together with the decrease in pixel% of -950 to -800 HU area, suggesting the increase of LAA% due to parenchymal destruction at -950 to -800 HU area. In contrast, the decrease in pixel% of -900 to -850 HU area markedly suppressed together with the decrease in the pixel% of -800 to -750 HU area 24 months after sirolimus therapy. Thus, the pathobiological mechanism(s) contributing to the increase in LAA% appear to be different between pre-trial and trial periods. LAA% is likely to increase due to decrease of -950 to -800 HU area (probably destruction of the alveoli) caused by proliferating LAM cells, whereas the complete recovery of the decrease of -950 to -800 HU area and the decrease of -800 to 750 HU area result in the increased LAA% after sirolimus therapy. These changes of histogram distribution during the trial period might reflect the inhibition of LAM cell proliferation and LAM-associated lymphatic congestion in the lung parenchyma.

Table 3

Longitudinal changes of pulmonary function and CT parameters in patients with pre-trial data (n = 8).

		Pre-treatment	Baseline	24 months	p value
Lung function
FVC	(L)	2.87 ± 0.15	2.87 ± 0.16	2.92 ± 0.13	NA
Difference from Baseline	(mL)	0 ± 55.2	NA	53.8 ± 89.9	0.484
FEV1	(L)	1.46 ± 0.17	1.36 ± 0.16	1.57 ± 0.22	NA
Difference from Baseline	(mL)	92.5 ± 55.3	NA	203.8 ± 85.3	0.0361
CT parameters
CT-TLC	(mL)	4586 ± 178	4795 ± 236	4776 ± 191	NA
Difference from Baseline	(mL)	-209 ± 117	NA	-19.4 ± 89.2	0.208
LAA%	(%)	27.9 ± 3.7	31.2 ± 3.1	34.1 ± 2.8	NA
Difference from Baseline	(%)	-3.29 ± 1.27	NA	2.91 ± 1.30	0.889
Histogram analysis
Mean CT value	(HU)	-879 ± 7	-884 ± 6	-892 ± 5.5	NA
Difference from Baseline	(HU)	5.05 ± 3.66	NA	-8.65 ± 3.67	0.327
Mode CT value	(HU)	-939 ± 16	-945 ± 12	-952 ± 11	NA
Difference from Baseline	(HU)	5.46 ± 5.65	NA	-7.38 ± 4.05	0.326
Kurtosis		6.53 ± 0.57	6.01 ± 0.50	6.73 ± 0.63	NA
Difference from Baseline		0.52 ± 0.23	NA	0.72 ± 0.26	0.0361
Skewness		2.23 ± 0.10	2.10 ± 0.08	2.24 ± 0.10	NA
Difference from Baseline		0.13 ± 0.04	NA	0.14 ± 0.05	0.0171
Fractal analysis
r2 value		0.937 ± 0.006	0.959 ± 0.007	0.969 ± 0.005	NA
Difference from Baseline		-0.022 ± 0.005	NA	0.010 ± 0.004	0.327
D-value		1.049 ± 0.071	1.019 ± 0.072	0.937 ± 0.050	NA
Difference from Baseline		0.030 ± 0.019	NA	-0.082 ± 0.028	0.161

Ai = luminal area; CT-TLC = computed-tomography-derived total lung capacity; FEV1 = forced expiratory volume in 1 s; FVC = forced vital capacity; GnRH = gonadotropin releasing hormone analogue; LAA% = percentage of low-attenuation area; NA = not applicable.

p < 0.05 generated by the Wilcoxon signed-rank test.

Figure 3

Longitudinal changes of FEV1 (A), LAA% (B), kurtosis (C) and skewness (D) in the 8 patients with pre-trial data. Gray and black lines represent the individual data and mean values, respectively. FEV1, kurtosis, and skewness decreased from pre-treatment to baseline, but increased from baseline to 24 months. The change in LAA% from pre-treatment to 24 months was not statistically significant. p values were generated by the Wilcoxon signed-rank test.

Correlation between pulmonary function and CT measurements

We examined the correlations between the changes (Δ) from baseline to 24 months in pulmonary function and CT-measured parameters (Figure 4). The correlation between ΔFEV1 and Δkurtosis (r = 0.431, p = 0.066) was borderline and the correlation between ΔFEV1 and Δskewness (r = 0.465, p = 0.045) was marginally significant.

Figure 4

Correlations between the changes of FEV1 and kurtosis (A), and FEV1 and skewness (B) (Baseline to 24 months). From baseline to 24 months, no statistically significant correlation was found between ΔFEV1 and Δkurtosis (Fig. 4A, r = 0.431, p = 0.066), but a statistically positive correlation was demonstrated between ΔFEV1 and Δskewness (Fig. 4B, r = 0.465, p = 0.045). P values were calculated using the Spearman's rank-correlation test.

Discussion

Our study is the first to demonstrate that the parameters reflecting lung density histogram, kurtosis and skewness, exhibited considerable changes during sirolimus therapy and had a statistically significant correlation with improved pulmonary function. Additionally, we found changes in lung density histograms between pre-trial and trial periods, clearly suggesting that the lung density histogram reflects pathobiological changes driven by sirolimus therapy. Preceding studies reported the relationship between pulmonary function and kurtosis or skewness in patients with chronic obstructive pulmonary disease (COPD) and interstitial lung diseases (ILD) [15, 19]. Lower kurtosis and skewness were observed in more progressed disease and severely compromised lung function in ILD [11, 19].

Proliferating LAM cells in lung parenchyma and distal airways lead to connective tissue matrix degradation due to protease-antiprotease imbalance in local milieus and air trapping due to airway narrowing, respectively, and ultimately generate pulmonary cysts [10, 20, 21]. Moreover, LAM cells induce lymphatic vessel hyperplasia and lymphatic obstruction resulting in pulmonary lymphedema, and also involve pulmonary vessels leading to pulmonary hemorrhage [10, 21, 22, 23]. These pathological changes can explain the characteristic alteration of pulmonary function in LAM patients: airflow obstruction; air-trapping; hyperinflation; and decrease in diffusing capacity. Sirolimus cannot eliminate LAM cells and exert just cytostatic effect [24], but it will cause cell size reduction since the mTOR signaling pathway is well-known to regulate cell size [25]. Hyperactivated mTOR pathways result in plumpy and large cells. On the other hand, sirolimus impedes lymphangiogenesis and inhibits function of lymphatic endothelial cells [26], and also has been reported to decrease circulating LAM cells in blood and serum VEGF-D concentrations [5, 27]. These actions of sirolimus seem to clear lymphatic congestion in the lungs [28], reduce or eliminate chylous effusion [4, 29], and shrink lymphangioleiomyomas [29]. In the light of this LAM pathophysiology and the action of sirolimus as well earlier study reporting that the pixel% of -900 to -750 HU represents lung tissue in the denser normal lung range [30], the reduction of lung density of -850 to -800 HU during the sirolimus therapy may reflect the reduction of cell size, lymphatic flow, parenchymal hemorrhage, etc. and convert these density areas to areas of increased air density (i.e., -1050 to -950 HU).

The analysis of the fractal property of LAA cluster (represented by r ≥ 0.95) is an established method to evaluate terminal airspace geometry in COPD [16, 31]. It maintains that regardless of the disease severity, the D-value sensitively detects alveolar tissue destruction [16, 31]. We previously reported that the percentage of CT images with fractal property tended to decrease when the disease became more severe in LAM patients [10], since medium-sized LAA clusters increased as LAA% increased, suggesting that LAA clusters had low tendency to fuse in LAM as compared with COPD. Therefore, pathological changes other than pulmonary cysts (e.g., LAM cell proliferation, lymphatic vessel hyperplasia, and lymphatic obstruction) being greatly different from COPD, may play a key role in preventing expansion and fusion of cysts. In the present study, the percentage of patients with fractal property (r ≥ 0.95) tended to increase during sirolimus treatment, suggesting the modifying effect of sirolimus on LAM progression. However, the trends of both r and D-values were similar in both pre-trial and treatment periods; the mean value of r tended to increase and the mean D-value tended to decrease. The pathophysiological mechanisms behind the changes of these parameters would not be identical in pre-trial and treatment periods, but r and D-values might not be appropriate indicators of the efficacy of sirolimus.

The quantitative analysis of airway dimensions was reported to be useful for evaluating disease progression and the effect of treatment in patients with COPD and asthma [32, 33]. Since our Juntendo cohort showed preserved or improved pulmonary function (i.e., FVC and FEV1), we measured airway dimension parameters (WA% and Ai) using cross-sectional images of airways at RB1 and basal bronchus (RB10). Although involvement of LAM cells into the bronchus and bronchioles had already been reported by histopathologic studies [34, 35], these parameters showed no significant longitudinal changes during sirolimus treatment. We measured AWT-Pi10, a parameter reflecting the dimension of peripheral airway [36] as well as the irregularity along the longitudinal lower lobe bronchus [33] to find no significant correlation with pulmonary function during sirolimus therapy (data not shown).

Argula et al. determined the baseline to 12-month change in CT-derived lung volumes and the volume of lung occupied by cysts as well as a Markov transition chain analysis of respiratory cycle cyst volume changes using data of inspiratory and expiratory CT images in the 31 MILES participants (17 in sirolimus group; 14 in placebo group) [9]. They found greater dynamic variation of respiratory cycle cyst volume changes in the sirolimus group than in the placebo group, indicating improvement of airway obstruction facilitating cyst emptying, which probably derived from a reduction in LAM cell burden around small airways and cyst walls. In the MILES cohort, they did not analyze airway dimensions as was done in our study, but suggested that sirolimus changed dynamic properties of airways rather than static ones, which lead to greater cyst volume emptying. Earlier quantitative CT studies consistently found a correlation between quantitative CT measures and severity of physiologic dysfunction; the more cyst area (LAA%) [6, 7, 37], cyst volume [38], or cyst score by texture analysis [8] increased, the more FEV1 decreased. However, it is apparent that increased LAA% resulted in the increase of FEV1 during sirolimus therapy. On inspiratory CT, it has been reported that there is a significant relationship between lung density and lung weight (eg. pulmonary oedema). In LAM patients, lymphatic vessel hyperplasia, lymphatic obstruction resulting in pulmonary lymphedema, and pulmonary hemorrhage may increase lung density, and these changes are reported to improve by sirolimus therapy. Therefore, the changes of LAA% and lung density on inspiratory CT images may represent these effects of sirolimus in our Juntendo cohort. Moreover, considering the remarkable lung density changes during sirolimus therapy that was revealed in this study and the results of dynamic evaluation of cyst volume changes by Argula et al., we speculate that the changes of LAA% and lung density seems to be one of the factors leading to the improvement of airflow limitation in our cohort. It is difficult to distinguish whether an increase of LAA% is caused by (I) the reduction of interstitial edema/LAM cell burden around small airways and cyst walls or by (II) increase of air volume accompanied by improvement of dynamic function of the pulmonary cysts. Due to the limited spatial resolution of CT scanners, the reduction of interstitial edema/LAM cell burden around small airways and cyst walls are not directly visualized on CT images. Therefore, dynamic study such as the combinatory analysis of inspiratory and expiratory CT prior and after sirolimus treatment like the study of Argula et al. [9] is needed to confirm. However, it is difficult to perform expiratory CT in clinical settings because of its radiation exposure and cost. Therefore, we believe our results obtained by inspiratory CT only are useful for clinical practice. Clinicians should be aware that LAA% increases both when sirolimus works and when LAM progresses, and need to evaluate its change carefully.

There are several study limitations. First, this is a subgroup analysis of the MLSTS study and the number of LAM patients analyzed was limited. Moreover, we combined 3 different data sets for the analysis: pre-trial data from 8 patients; trial data from 20 patients; and follow-up data from 16 patients. However, since the scanner and the data acquisition methods are identical, incorporating pre-trial trends and analysis over long time periods (approximately 2 years in pre-trial and 4 years in sirolimus treatment) may be meritorious. Second, we failed to conduct histopathological examination to explain the sirolimus-induced lung density change. We believe, however, that our quantitative CT analysis compensates for the lack of histopathological evaluation. Third, we used the kernels for bone (FC85), which can lead an overestimation of LAA%. In this study, we could not obtain CT images reconstructed with more adequate kernels (i.e. FC03). However, we consider that our method is appropriate because this method is often used in past reports and also the investigation of change over time is the main purpose of this research. Finally, the 3-slice measurement is old-fashioned and voxel analysis (analyzing volumetric data) is now popular. However, because pulmonary cysts distribute diffusely in the lungs of LAM patients, we do not think that the whole lung analysis is necessary. Most of the previous reports on fractal analysis of emphysema or pulmonary cysts used 2D analysis, and they reported that “three-slice methods” were sufficient for the analysis of emphysema or pulmonary cysts [10, 16].

In conclusion, our study suggests that lung density histogram parameters, kurtosis, and skewness, may be useful as indicators of the efficacy of sirolimus. Stabilization and improvement of airflow limitation obtained by sirolimus therapy may be in part due to the change in lung mechanics derived from the change in lung density histogram.

Declarations

Author contribution statement

Y. Ko: Conceived and the designed the experiments; Performed the experiments; Analyzed and interpreted the data; Wrote the paper.

K. Seyama, K. Tobino: Conceived and the designed the experiments; Performed the experiments; Analyzed and interpreted the data; Wrote the paper.

K. Asakawa, T. Oguma: Performed the experiments; Contributed reagents, materials, analysis tools or data.

T. Hirai, T. Takada, K. Takahashi: Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data.

Funding statement

This work was supported by a grant from the Ministry of Health, Labour, and Welfare (grant number H24-Nanchi-Ippan-010).

Competing interest statement

The authors declare no conflict of interest.

Additional information

The clinical trial described in this paper was registered at Japan Medical Association, Center for Clinical Trials under the registration number JMA-IIA00096.

Supplementary content related to this article has been published online at https://doi.org/10.1016/j.heliyon.2020.e03345.