Sputum Exosomal microRNAs Profiling Reveals Critical Pathways Modulated By Pseudomonas aeruginosa Colonization In Bronchiectasis.

Background: Pseudomonas aeruginosa (PA) colonization confers poor prognosis in bronchiectasis. However, the biomarkers and biological pathways underlying these associations are unclear. Objective: To identify the roles of PA colonization in bronchiectasis by exploring for sputum exosomal microRNA profiles.
Methods: We enrolled 98 patients with clinically stable bronchiectasis and 17 healthy subjects. Sputum was split for bacterial culture and exosomal microRNA sequencing, followed by validation with quantitative polymerase chain reaction. Bronchiectasis patients were stratified into PA and non-PA colonization groups based on sputum culture findings. We applied Gene Ontology and Kyoto Encyclopedia of Genes and Genome pathway enrichment analysis to explore biological pathways corresponding to the differentially expressed microRNAs (DEMs) associated with PA colonization.
Results: Eighty-two bronchiectasis patients and 9 healthy subjects yielded sufficient sputum that passed quality control. We identified 10 overlap DEMs for the comparison between bronchiectasis patients and healthy subjects, and between PA and non-PA colonization group. Both miR-92b-5p and miR-223-3p could discriminate PA colonization (C-statistic >0.60) and independently correlated with PA colonization in multiple linear regression analysis. The differential expression of miR-92b-5p was validated by quantitative polymerase chain reaction (P<0.05), whereas the differential expression of miR-223 trended towards statistical significance (P=0.06). These DEMs, whose expression levels correlated significantly with sputum inflammatory biomarkers (interleukin-1β and interleukin-8) level, were implicated in the modulation of the nuclear factor-κB, phosphatidylinositol and longevity regulation pathways.
Conclusion: Sputum exosomal microRNAs are implicated in PA colonization in bronchiectasis, highlighting candidate targets for therapeutic interventions to mitigate the adverse impacts conferred by PA colonization.

Introduction

Chronic airway infection and inflammation are implicated in the pathogenesis and progression of bronchiectasis.1 The dilated bronchi provide the niche for the expansion of opportunistic pathogens that elicit airway inflammation.2 Compared with other potentially pathogenic bacteria, Pseudomonas aeruginosa (PA) conferred significantly greater adverse impacts on bronchiectasis,3 including the more rapid lung function decline and frequent exacerbations.4 In fact, colonization of PA correlated with a substantial morbidity and mortality.5

Micro-ribonucleic acids (miRNA), the 18–25nt non-coding ribonucleic acids, have been implicated in the airway inflammatory responses.6 Up-regulated expression of miR-301b correlated with the heightened inflammation associated with PA infection.7 Exosomes are extracellular vesicles laden mainly with non-coding RNAs such as miRNAs,8 which are transferred to the target cells for modulating the host immunity.9 For instance, exosomal miR-155 enhanced whereas miR-146a suppressed inflammatory gene expression,10 suggesting that exosomal miRNA might be capable of modulating airway inflammation.

Currently, the mechanisms underlying the adverse impacts associated with PA colonization in bronchiectasis are unclear. Here, we sought to identify the differentially expressed miRNAs (DEMs) in sputum exosomes that were associated with PA colonization. This might help clarify whether and which miRNAs could be the biomarkers associated with PA colonization, and through which pathways PA modulates chronic airway inflammation in bronchiectasis.

Methods

Participants

Between November 2016 and March 2018, we consecutively recruited 98 bronchiectasis patients and 17 healthy subjects. Bronchiectasis was diagnosed with chest high-resolution computed tomography in the out-patient clinics of The First Affiliated Hospital of Guangzhou Medical University. Eligible patients remained clinically stable (respiratory symptoms not significantly exceeding normal daily variations) for at least 4 weeks. Healthy subjects had normal chest X-ray findings, and had no history of chronic respiratory diseases. None of the participants had malignancy, acute upper airway infection or antibiotic use within 4 weeks before recruitment. Exacerbations were defined based on the criteria stated in international guidelines.11 The study protocol was approved by The Ethics Committee of The First Affiliated Hospital of Guangzhou Medical University (Medical Ethics Year 2012, the 33rd). Written informed consent was obtained.

Study Design

Patients underwent history inquiry, spirometry (COSMED Inc., Italy) and collection of spontaneous sputum when clinically stable. Patients were followed-up at 3-to-6-month intervals to determine the bacterial colonization status, and record exacerbations and hospitalizations. We excluded patients followed-up for less than 6 months. Modified Reiff score was used to evaluate the radiologic severity of bronchiectasis.12,13 Etiologies, Bronchiectasis Severity Index (BSI), and E-FACED score were assessed5,14,15 at the initial stable visit. We validated the DEMs with quantitative polymerase chain reaction (qPCR) in 34 bronchiectasis patients who had sufficient paired sputum.

Procedures

Patients thoroughly gargled their mouth with distilled water, followed by forceful coughing. Sputum was expectorated in a sterile plastic pot. We selected the most purulent portion and performed cytology assessment under the light microscope. Sputum that met quality control criteria was split for bacterial culture (within 2 hrs), sputum sol preparation for inflammatory biomarker [interleukin (IL)-8, IL-1β, and tumor necrosis factor-α (TNF-α)] assays, and storage in −80°C freezers for subsequent extraction of the exosomes and miRNA detection. Sputum exosomes were identified with the H-7650 transmission electron microscope (HITACHI, Japan), and the DEMs associated with PA colonization were validated with qPCR by using ABI 7500 Real-Time System (Thermo Scientific, USA).

Definition Of PA Colonization

PA colonization was defined as sputum culture positive of PA for at least two occasions (at least 3 months apart) within 1 year, when clinically stable during the longitudinal follow-up. Bronchiectasis patients were stratified into PA and non-PA group (unmet the criteria of PA colonization).

Sputum miRNA Expression Profiling

After thawing, sputum was mixed with phosphate buffer solution. Exosomes were isolated with exoEasy Maxi Kit, and the total RNA was extracted with exoRNeasy Serum/Plasma Maxi Kit (Qiagen, Germany). Agilent Bioanalyzer 2100 system (Agilent Technologies, USA) was used to assess the purity and integrity of RNA, followed by sequencing with Hiseq 2500 platform (Illumina, USA). After denoising the raw reads, we calculated the false discovery rate (FDR), P value and fold-change for high-quality reads. MiRNA expression was corrected by transcripts per million.16,17 Mature human miRNAs were annotated with miRBase 22.0. See regarding the flowchart of sequencing schemes. See further details in .

Statistical Analysis

Statistical analyses were conducted with R software version 3.5.0 and Graphpad Prism 5.0 (Graphpad Inc., USA). Kolmogorov–Smirnov test was applied to determine the normality of continuous variables, which were presented as mean ± standard deviation or median (interquartile range) as appropriate. Count (percentage) was demonstrated for categorical data.

Principal component analysis (PCA) was applied for quality control and assessment of the global gene expression patterns. Identification of the DEMs was based on the FDR (cut-off: 0.10).18–20 The DEMs were initially screened via comparison between all bronchiectasis patients and healthy subjects, followed by the comparison between PA and non-PA group in bronchiectasis patients. Venn diagram was applied to demonstrate the unique or overlap DEMs. Multiple-group comparisons were made using analysis-of-variance or non-parametric analysis and adjusted with Benjamini–Hochberg correction. The diagnostic value of miRNA in differentiating PA group from non-PA group was analyzed with the C-statistic and 95% confidence interval (95% CI). To select the minimum set of DEMs that best discriminated PA colonization, machine learning was applied to perform least absolute shrinkage and selection operator linear regression analysis (LASSO)21 to avoid overfitting. Multiple linear regression was performed to analyze the correlation between DEMs and PA colonization. We finally included the miRNA for validation by performing qPCR assays, provided that the C-statistic was greater than 0.6, and that the expression level correlated significantly with PA colonization in the above-mentioned multiple regression models. The biological function and pathways of DEMs were analyzed with Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes Pathway Analysis (KEGG). Significant enrichment was defined as the unadjusted P value <0.01 and Padj <0.05 (corrected with the Benjamini–Hochberg algorithm).

Results

Subject Recruitment And Clinical Characteristics

Of the 115 participants who underwent screening, 82 bronchiectasis patients and 9 healthy controls were included in the final analysis. Twenty-four participants were excluded because of incomplete bacterial culture data or low-quality sequencing results (Figure 1).

Figure 1

Study participants recruitment flow chart. We screened for 115 study participants, 22 of whom yielded sputa with a poor quality and were therefore excluded for sequencing. Two patients whose data for sequencing with low quality were also excluded from further analysis. Hence, the final cohort consisted of 91 study participants (82 with stable bronchiectasis and 9 healthy subjects).

Table 1 shows the demographic and clinical characteristics of the study participants. Most patients had post-infectious and idiopathic bronchiectasis. Thirty-four bronchiectasis patients (41.5%) had PA colonization, who had worse lung function and greater disease severity than those without. The demographic and clinical characteristics of PA and non-PA group are shown in , and those in RT-qPCR cohort are shown in .

Table 1

Clinical Characteristics Of Patients With Clinically Stable Bronchiectasis And Healthy Subjects

Parameter	Bronchiectasis Patients	Healthy Subjects	P value
No. of subjects	82	9	–
PA colonization (No., %)	34 (41.5%)	NA	NA
Anthropometry
Age (years)	45.4±14.0	45.8±16.5	0.352
Females (No., %)	49 (59.8%)	4 (44.4%)	0.483
Body-mass index (kg/m2)	20.2±2.8	22.6±2.0	0.432
Never-smokers (No., %)	6 (7.3%)	1 (11.1%)	0.530
FEV1 predicted%	57.7±22.3	95.4±8.8	<0.001
Exacerbation history
No. of exacerbations within 1 yr	1.4±1.1	NA	NA
Severity
Bronchiectasis Severity Index	6.6±4.0	NA	NA
E-FACED score	3.7±1.9	NA	NA
Modified Reiff score	9.8±4.1	NA	NA
Aetiology
Idiopathic (No., %)	25 (30.5%)	NA	NA
Post-infectious (No., %)	33 (40.2%)	NA	NA
Kartagener syndrome (No., %)	8 (9.8%)	NA	NA
Asthma (No., %)	2 (2.4%)	NA	NA
Immunodeficiency (No., %)	5 (6.1%)	NA	NA
GORD-aspiration (No., %)	4 (4.8%)	NA	NA
Others (No., %)a	7 (8.5%)	NA	NA
Bacteria isolated during the initial visit
Pseudomonas aeruginosa	36 (43.9%)	NA	NA
Haemophilus influenzae	9 (11.0%)	NA	NA
Escherichia coli	4 (4.9%)	NA	NA
Other pathogensb	7 (8.5%)	NA	NA
Commensal	26 (31.7%)	NA	NA
Sputum physical characteristic
Viscosity score	2.5±0.7	2.3±0.9	0.747
Purulence score	6.2±0.9	4.0±0.5	<0.001

Notes: P-value for comparisons between patients with clinically stable bronchiectasis and healthy subjects. Othersa included congenital airway defects, diffuse panbronchiolitis and connective tissue disease. bOther pathogens included Klebsiella pneumiae, Proteus mirabilis, Streptococcus pneumoniae, Haemophilus parainfluenzae, Haemophilus parahaemolyticus, Sphingomonas sphingomonas, Enterobacter aerogenes.

Abbreviations: E-FACED, exacerbation-FACED; GORD, gastroesophageal reflux disease-aspiration.

Main Sequencing Findings

We first applied the transmission electron microscope to confirm the successful extraction of extracellular vesicles (50–200 nm, ). On average, 10 million raw reads per sample were sequenced. We identified 1254 miRNAs in PA group, 1403 miRNAs in non-PA group, and 1132 miRNAs in healthy subjects. Next, we applied PCA and dendrogram to demonstrate the global miRNA expression patterns (). There was some separation between PA group and healthy subjects, but less so between PA and non-PA group. The flowchart of the analysis is shown in .

DEMs Between PA And Non-PA Group

After quality control of sequencing data, we identified 190 (129 up-regulated and 61 down-regulated) DEMs between all bronchiectasis patients and healthy subjects, and 22 DEMs (FDR<0.10) between PA and non-PA group in bronchiectasis patients. Of these, 10 (5 up-regulated and 5 down-regulated) DEMs in PA group discriminated from both non-PA group and healthy subjects (Figure 2A–C). Figure 2D demonstrates the expression levels of these 10 DEMs in PA and non-PA group. We noted significantly higher levels of miR-628-5p in PA group than in non-PA group (fold-change: 1.72, FDR: 0.03). The expression levels of the overlap DEMs are shown in .

Figure 2

miRNA expression differentiating PA colonization status. (A) 190 DEMs were differentially expressed between PA group and healthy subjects, and 22 DEMs between PA group and non-PA group. There were 10 overlapped DEMs for the comparison of PA and non-PA group. Five of the DEMs were down-regulated, whereas five were up-regulated. (B) The volcano plot showing the log2(fold-change) and log2FDR of the 190 differentially expressed miRNAs between PA group and healthy subjects. The vertical line was set at the Log2(fold-change) of +1 and −1, respectively. (C) The volcano plot showing the log2(fold-change) and log2FDR of these 22 differentially expressed miRNAs between bronchiectasis patients with PA colonization and those without (non-PA colonization). The vertical line was set at the Log2(fold-change) of +1 and −1, respectively. (D) A heat map demonstrating 10 DEMs (5 up- and 5 down-regulated) after correction for multiple testing. Each column represents an individual patient sample, and each row represents an miRNA. The color code indicates the expression level, with red indicating higher expression levels and blue indicating lower expression levels. (E, F) Receiver operation characteristic curve demonstrating the diagnostic value of the 10 DEMs and the LASSO regression model to differentiate bronchiectasis patients with PA colonization from those without PA colonization.

Abbreviations: PA, P. aeruginosa; DEM, differential expressed miRNA; FDR, false-discovery rate; AUC, area under the curve; Down, down-regulated (lower expression levels) in PA group compared with non-PA; Up, up-regulated (higher expression levels) in PA group compared with non-PA; LASSO, least absolute shrinkage and selection operator linear regression analysis.

The Diagnostic Value Of DEMs

Next, we analyzed how well the 10 DEMs discriminated against the PA colonization status (Figure 2E and F). Overall, most miRNAs had the modest diagnostic value in discriminating the PA group from non-PA group. miR-628-5p (C-statistic: 0.71, 95% CI: 0.60, 0.83), miR-223-3p (C-statistic: 0.71, 95% CI: 0.59, 0.83), miR-23a-3p (C-statistic: 0.71, 95% CI: 0.59, 0.83) yielded the highest diagnostic value. Other miRNAs yielded the C-statistics ranging from 0.40 to 0.69. To improve the discriminative capacity, we also performed the LASSO regression model in which we included 74 bronchiectasis patients (>90%) in the training set and eight in the validation set. The accuracy of the LASSO model (see details in ) to discriminate PA colonization in the validation set was 87.5%. The LASSO model yielded higher diagnostic value (C-statistic: 0.87, 95% CI: 0.79, 0.95) than any single miRNA in discriminating PA colonization.

DEMs Associated With PA Colonization

We performed multiple linear regression analysis to further validate the correlation between PA colonization and miRNA expression. miR-223-3p, miR-23a-3p, miR-628-5p, miR-92b-5p and miR-548f-3p expression significantly and independently correlated with PA colonization (all P<0.05, Table 2), but the remaining five DEMs did not (all P>0.05, ). We validated the former five DEMs with qPCR assay among 34 bronchiectasis patients who had sufficient paired sputum samples. The correlation coefficient of the readings in duplicate wells was 0.943. miR-92b-5p was significantly down-regulated in PA group (fold-change: −2.94, P<0.05), while miR-223-3p trended towards an increased expression in PA group (fold-change: 2.06, P=0.063) compared with non-PA group (, ). In light of the discriminative capacity, the significant correlation with PA colonization, we have included miR-223-3p in pathway analysis despite that it narrowly missed the statistical significance in qPCR assays.

Table 2

Multiple Linear Regression Between The DEMs Associated With PA Colonization And Other Clinical Parameters

miRNA ID	Adjusted R2	P-value Of The Model	Predictor Variables	Β	P value Of β
miR-23a-3p	0.170	0.005	Age	21.46	0.099
			Sex	220.09	0.550
			BMI	−95.53	0.150
			FEV1 pred%	−17.69	0.107
			No. of lobes affected	32.10	0.877
			HRCT Reiff scores	−15.19	0.849
			PA colonization	982.60	0.012
			Other bacterial colonization	−11.50	0.983
miR-92b-5p	0.127	0.020	Age	4.88	0.490
			Sex	57.74	0.775
			BMI	69.43	0.058
			FEV1 pred%	−5.19	0.386
			No. of lobes affected	119.62	0.296
			HRCT Reiff scores	−76.27	0.086
			PA colonization	−541.08	0.011
			Other bacterial colonization	−491.09	0.102
miR-548f-3p	0.116	0.028	Age	0.68	0.135
			Sex	−0.30	0.981
			BMI	−0.76	0.745
			FEV1 pred%	−0.394	0.306
			No. of lobes affected	−3.83	0.601
			HRCT Reiff scores	3.11	0.272
			PA colonization	31.03	0.023
			Other bacterial colonization	26.80	0.163
miR-223-3p	0.197	0.002	Age	30.51	0.106
			Sex	−336.37	0.529
			BMI	−102.57	0.286
			FEV1 pred%	−23.80	0.135
			No. of lobes affected	49.27	0.870
			HRCT Reiff scores	26.99	0.816
			PA colonization	1748.27	0.002
			Other bacterial colonization	930.22	0.239
miR-628-5p	0.164	0.006	Age	0.22	0.478
			Sex	5.12	0.555
			BMI	0.52	0.740
			FEV1 pred%	−0.27	0.286
			No. of lobes affected	0.726	0.882
			HRCT Reiff scores	0.72	0.704
			PA colonization	27.37	0.003
			Other bacterial colonization	−2.93	0.818

Note: P-value that corresponded to the β value of less than 0.05 is presented in bold.

Abbreviations: BMI, body mass index; FEV1, forced expiratory volume in 1 s; PA, P. aeruginosa.

miR-92b-5p And miR-223-3p Correlated With The Clinical Parameters In Bronchiectasis

We next explored the correlation between miR-92b-5p and miR-223-3p and the clinical parameters in bronchiectasis (corrected with Benjamini–Hochberg algorithm, Figure 3). The expression level of miR-223-3p correlated negatively with that of miR-92b-5p (r=−0.65 in PA group and −0.61 in non-PA group, both Padj<0.05). The expression level of miR-223-3p positively correlated with that of sputum IL-8 (r=0.61, Padj<0.05), and the expression level of miR-92b-5p negatively correlated with that of sputum IL-1β (r=−0.63, Padj<0.05) and the BSI (r=−0.58, Padj<0.05) in PA group. Moreover, the expression level of miR-92b-5p negatively correlated with that of sputum IL-1β (r=−0.49, Padj<0.05). Neither the expression level of miR-223-3p nor miR-92b-5p correlated with the number of bronchiectatic lobes, lung function or the E-FACED scores. We noted a less prominent correlation between the expression level of miR-223-3p and miR-92b-5p and that of sputum IL-1β in non-PA group than in PA group. Nonetheless, the BSI correlated with the expression level of miR-223-3p (r=0.18, P<0.001) but not miR-92b-5p (). However, there was no marked difference in the expression levels of miR-223-3p and miR-92b-5p between severe and mild-to-moderate bronchiectasis (both P>0.05, ).

Figure 3

miR-92b-5p and miR-223-3p correlated with the clinical parameters of bronchiectasis. (A) Correlation matrix in PA group. (B) Correlation matrix in non-PA group. In both matrices, the lower left quadrant of the matrix represents the crude P values, and the upper right quadrant of the matrix presents the P values corrected with the Benjamini–Hochberg algorithm. The intersection within the circle represents the adjusted P value >0.05. Blue circle indicate positive correlation, while red circles denote negative correlation. Greater color intensity indicates a greater magnitude of correlation (correlation coefficient).

Abbreviations: IL-1β, interleukin-1β; IL-8, interleukin-8; TNF-α, tumor necrosis factor-α; lobe, The number of affected lobes; FEV1%, forced expiratory volume in 1 s; BSI: bronchiectasis severity index; E-FACED, exacerbation-FACED score.

Biological Pathways Of miR-92b-5p And miR-223-3p

Finally, we analyzed the canonical pathways of miR-223-3p and miR-92b-5p by using clusterProfiler package,22 based on the 826 target mRNAs of these two DEMs that were predicted by the DIANA-microT web server version 5.0 and validated with the DIANA-TarBase version 8.0 database.23,24

Gene Ontology Analysis

We identified the biological processes associated with inflammation, such as leukocyte differentiation, cellular response to oxidative stress. See further details in .

KEGG And Genomic Analysis

KEGG analysis revealed 12 canonical pathways, five of which were not further analyzed because they were related to cancer pathways. Table 3 shows the most significant pathways of the annotated mRNAs. The interaction network of mRNAs was classified into four domains: inflammatory processes, hormonal responses, lung aging and stem cell (in which the nuclear factor-κB [NF-κB] signaling pathway, phosphatidylinositol signaling system were the top enriched pathways in the inflammatory processes). See regarding all enriched signaling pathways.

Table 3

The Significant Biological Pathways Of The DEMs Associated With PA Colonization In Bronchiectasis Patients

Description And ID	Ontology	Count	P value	FDR	Representative Genes
Part A: Inflammatory processes
Leukocyte differentiation (GO: 0002521)	BP	38	8.56E-05	2.69E-02	mTOR, PIK3CD, MMP9, STAT3
Response to oxygen levels (GO: 0070482)	BP	31	1.82E-04	4.11E-02	FoxO1, ATP7A
Cellular response to oxidative stress (GO: 0034599)	BP	25	2.89E-04	4.38E-02	CHUK, FoxO1, FoxO3
Myeloid leukocyte differentiation (GO: 0002573)	BP	22	1.12E-05	5.46E-03	PARP1, VEGFA
NF-kappa B signaling pathway (KEGG: hsa04064)	KEGG	11	1.23E-03	3.36E-02	PRKCB, PARP1, ICAM1, BCL2, CHUK
Phosphatidylinositol signaling system (KEGG: hsa04070)	KEGG	11	1.72E-03	3.54E-02	PIK3CD, PRKCB, PI3K-C2α
Part B: Hormone response
Insulin resistance (KEGG: hsa04931)	KEGG	13	2.81E-04	2.30E-02	STAT3, FoxO1
Endocrine resistance (KEGG: hsa01522)	KEGG	11	1.59E-03	3.54E-02	IGF1R, MMP9, PIK3CD, mTOR
Part C: Lung aging
Longevity regulating pathway (KEGG: hsa04211)	KEGG	11	7.08E-04	3.36E-02	KL, IGF1R, FoxO1, FoxO3
Longevity regulating pathway – multiple species (KEGG: hsa04213)	KEGG	11	2.54E-05	6.27E-03	FoxO1, mTOR
Part D: Stem cell
Signaling pathways regulating pluripotency of stem cells (KEGG: 04550)	KEGG	14	1.13E-03	3.36E-02	IL6ST, IGF1R, PIK3CD, LIF, SMAD5

Notes: A total of 826 targeted genes (DEGs) were predicted by DIANA-microT and TarBase databases for miR-92b-5p and miR-223-3p. Hence, we identified 59 gene ontologies by GO analysis, which included 32 biological processes related to these DEGs, 20 related to cellular component and 7 related to molecular function. Meanwhile, a total of 12 biological pathway were enriched by KEGG analysis. We selected diseases-relevant ontology items as shown above.

Abbreviations: BP, biological processes; P-adj, P-value corrected with Benjamini–Hochberg algorithm; mTOR, mammalian target of rapamycin; PIK3CD, phosphatidylinositol-4,5-bisphosphate 3-kinase catalytic subunit delta; MMP9, matrix metalloprotein-9; STAT3, signal transducers and activators of transcription-3; FoxO1, forkhead box protein O-1; ATP7A, Menkes’ protein; CHUK, conserved helix-loop-helix ubiquitous kinase; FoxO3, forkhead box protein O-3; PARP1, poly(ADP-ribose) polymerase 1; VEGFA, vascular endothelial growth factor A; PRKCB, protein kinase C-β; ICAM1, intercellular adhesion molecule 1; BCL2, BCL2, apoptosis regulator; PI3K-C2α, class II phosphatidylinositol 3-kinase; KL, klotho; IL6ST, interleukin 6 signal transducer; IGF1R, insulin-like growth factor-1; LIF, leukemia inhibitory factor; SMAD5, SMAD family member 5.

Discussion

PA is a common opportunistic pathogen which correlates with poor prognosis in bronchiectasis.25 The known mechanisms included airway epithelial injury26,27 and toxin secretions (such as cyanide and pyocyanin).28,29 However, little is known regarding other candidate biological pathways underlying the adverse impacts of PA in bronchiectasis. We have detected sputum exosomal DEMs and analyzed their ability to discriminate PA colonization status and the correlation with clinical parameters in bronchiectasis. miR-223-3p and miR-92b-5p, which correlated significantly with sputum inflammatory biomarkers, were implicated in nuclear factor-κB, phosphatidylinositol and longevity regulation pathways. However, neither miR-223-3p nor miR-92b-5p correlated significantly with bronchiectasis severity (the BSI).

Previous studies demonstrated that sputum exosomes could be successfully extracted for biomarker assays,30–34 which was supported by our findings. Exosomal miRNAs play vital roles in mediating infection and inflammation in many chronic airway inflammatory diseases. For instance, miRNA expression in airway epithelial cells differed significantly between patients with cystic fibrosis and healthy subjects, and between patients with and without PA infection.35 Specifically, miR-223 was up-regulated along with increased IL-8 levels in the cells from patients with ΔF508 mutation35 in those with PA infection. miR-223 reportedly modulated the physiologic function of cystic fibrosis transmembrane conductance regulator and airway IL-8 expression, which might be the target of PA colonization in bronchiectasis. In neutrophilic asthma, miR-223 was up-regulated and positively correlated with sputum IL-1β levels.36 These observations hinted on how PA colonization might interact with the airway inflammation via miRNA modulation in bronchiectasis.

The inflammatory signaling pathways and the accelerated lung senescence might be associated with PA colonization in bronchiectasis. Upon PA infection, airway epithelial-derived exosomes37 stimulated toll-like receptor 4/MyD88/NF-κB signaling pathways.7 Increased levels of IL-8 of epithelial cells have also been detected in cystic fibrosis patients with PA infection.38 Congruent with these findings, our findings hinted on the potential link among PA colonization, altered miRNA expression, and NF-κB signaling pathway in bronchiectasis. Our findings suggested the significant correlation between miR-223-3p expression and sputum IL-8 levels, and between miR-92b-5p expression and sputum IL-1β levels. Phosphatidylinositol signaling pathway might also play a role in modulating inflammation associated with PA colonization in bronchiectasis. In fact, PA infection reportedly triggered the formation of airway epithelial cell protrusions leading to activation of phosphatidylinositol and NF-κB signaling pathways.39 Meanwhile, the toxin pyocyanin released by PA could mediate oxidative stress of airway epithelial cells by inducing the expression of antioxidative gene via the epidermal growth factor receptor-PI3K signaling pathway,40 which implied that phosphatidylinositol signaling pathway might play a role in bronchiectasis patients with PA colonization.

Apart from the inflammation-related pathways, longevity regulation pathway (ageing) might have also affected the prognosis of bronchiectasis associated with PA colonization. Ageing of the immune system predisposed to infection of opportunistic pathogens such as PA.41 In turn, PA infection accelerated cell senescence via NF-κB signaling pathway, which further hampered the eradication of PA.42 Collectively, our findings provide some hints on the complex network underlying the pathogenesis of bronchiectasis associated with PA colonization.

Our findings have clinical implications. Identification of sputum exosomal miRNA helps explore the important biological pathways underlying the pathogenesis of bronchiectasis, thus miRNAs and their pathways might serve as the complementary tool to explore for the candidate targets for future interventions of PA colonization. Apart from canonical inflammatory pathways, our results have implied the modulation of longevity regulation pathway which opens a new avenue for exploring how PA colonization interacts with the airway epithelium. The significant correlation between sputum inflammatory biomarkers and miR-92b-5p and miR-223-3p provided further evidence on the unresolved inflammation in the PA-colonized microenvironment.

There are notable limitations of our study. First, we compared miRNA expression in spontaneous sputum in bronchiectasis and in induced sputum in healthy subjects. Because all bronchiectasis patients could spontaneously produce sputum, we did not collect induced sputum in bronchiectasis patients. Salivary contamination has been minimized based on our experimental collection procedures. Second, the differential expression of miR-223-3p between PA and non-PA group did not reach to statistical significance (P=0.063) according to the qPCR assays given the insufficient sputum volume in the remaining patients. However, based on the discriminative capacity and the significant correlation with PA colonization and the fold-change of differential expression in both sequencing and qPCR assays, our findings of miR-223-3p expression patterns remained valid. We did not investigate the interactions between miRNA and mRNA. Besides, the association among airway neutrophilic inflammation, miRNA expression and PA colonization in bronchiectasis merits further investigation. Despite the efforts of identifying patients with PA colonization during longitudinal follow-up, we might have still underestimated the risk of misclassification (false-negative cases). Finally, all biological pathways were explored via bioinformatic analysis, further confirmatory experiments are warranted.

In summary, we have identified the differentially expressed sputum exosomal miRNA (miR-223-3p and miR-92b-5p) in bronchiectasis patients with PA colonization. The activation of inflammation and aging regulation pathways might provide clues for the intervention of PA colonization in bronchiectasis.