http://orcid.org/0000-0001-6936-103X
Figures
Abstract
Background
Differences in asthma severity may be related to inflammation in the airways. The lower airway microbiota has been associated with clinical features such as airway obstruction, symptom control, and response to corticosteroids.
Objective
To assess the relationship between local airway inflammation, severity of disease, and the lower airway microbiota in atopic asthmatics.
Methods
A cohort of young adult, atopic asthmatics with intermittent or mild/moderate persistent symptoms (n = 13) were assessed via bronchoscopy, lavage, and spirometry. These individuals were compared to age matched non-asthmatic controls (n = 6) and to themselves after six weeks of treatment with fluticasone propionate (FP). Inflammation of the airways was assessed via a cytokine and chemokine panel. Lower airway microbiota composition was determined by metagenomic shotgun sequencing.
Results
Unsupervised clustering of cytokines and chemokines prior to treatment with FP identified two asthmatic phenotypes (AP), termed AP1 and AP2, with distinct bronchoalveolar lavage inflammatory profiles. AP2 was associated with more obstruction, compared to AP1. After treatment with FP reduced MIP-1β and TNF-α and increased IL-2 was observed. A module of highly correlated cytokines that include MIP-1β and TNF-α was identified that negatively correlated with pulmonary function. Independently, IL-2 was positively correlated with pulmonary function. The airway microbiome composition correlated with asthmatic phenotypes. AP2, prior to FP treatment, was enriched with Streptococcus pneumoniae. Unique associations between IL-2 or the cytokine module and the microbiota composition of the airways were observed in asthmatics subjects prior to treatment but not after or in controls.
Citation: Turturice BA, McGee HS, Oliver B, Baraket M, Nguyen BT, Ascoli C, et al. (2017) Atopic asthmatic immune phenotypes associated with airway microbiota and airway obstruction. PLoS ONE 12(10): e0184566. https://doi.org/10.1371/journal.pone.0184566
Editor: Aran Singanayagam, Imperial College London, UNITED KINGDOM
Received: July 12, 2017; Accepted: August 25, 2017; Published: October 20, 2017
Copyright: © 2017 Turturice et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: Data for the samples that were sequenced have been uploaded onto MG-RAST and are available under the Project ID: 12564.
Funding: This work was supported by the National Institutes of Health and the National Institute of Allergy and Infectious Diseases, 5R01AI053878.
Competing interests: The authors have declared that no competing interests exist.
Introduction
Asthma is a heterogeneous disease. Asthma phenotypes have been established on the basis of underlying inflammation that may be predictive of treatment response and severity of disease. For example, atopic asthma has been associated with T-helper (Th) 2 responses and development of symptoms early in life, whereas other phenotypes have been associated with non-Th2 inflammation and development of symptoms later in life [1–3]. Recently, additional T-cell subsets have been determined to have a role in the pathogenesis of atopic asthma, such as Th17 and T-regulatory (Treg) cells [4–9]. These subsets are responsible for immunity to specific microbes or regulation of inflammation.
The composition of the airway microbiome in asthmatics has been associated with clinical features, such as symptom control, severity of airway obstruction, and response to corticosteroid therapy [10–13]. Pathogenic bacteria (e.g., Haemophilus, Streptococcus, Moraxella, Klebsiella, and Rickettsia) in the upper and lower airways have been associated with the diagnoses of asthma and potentially hold prognostic value for future diagnosis of asthma. [12, 14–16]. Bacteria of the upper airway and oral cavity can also be identified within the lower airways. Notably, a high abundance of supraglottic taxa in the lower airways of healthy individuals has been associated with Th17-like cytokine profiles in bronchoalveolar lavage (BAL) along with increased lymphocyte and neutrophil counts in BAL [17, 18]. Thus, there is a strong link with the airway microbiome and local airway inflammation, specifically with respect to T-cell phenotypes, in healthy individuals.
Recently, the relationship between the microbiome of the airways and immune response in the airway of asthmatics has been explored [12, 19, 20]. Prior studies focused on changes in the airway microbiome related to the presence of Th2 or Th17 responses in bronchial epithelium and the cellularity of eosinophils and neutrophils in the airways. We hypothesized that there are unique associations between inflammation, airway microbiome, and pulmonary function in asthmatics. We analyzed cytokine concentrations in bronchoalveolar lavage of young adult, atopic asthmatics (n = 13) and their relationship to the airway microbiome and lung function. In addition, we tested the impact of a short term, commonly used asthma medication, fluticasone propionate (FP) on airway inflammation and microbiome. Asthmatics were compared to age matched non-asthmatic individuals (n = 6) and to themselves after treatment with fluticasone propionate. The subjects were free of smoking, signs of respiratory infection, and corticosteroid therapy in the preceding months to control for potential confounders that can affect the microbiome and inflammation [11, 21, 22]. We identified two unique cytokine profiles that were associated with airway obstruction and compositional changes in the microbiome. These inflammatory changes were reduced after FP therapy, as well as association between the cytokines and the microbiome composition. This indicates that the relationship between inflammation and the microbiome is associated with pulmonary function in asthmatics.
Materials and methods
Subjects
The study was approved by the Ethics Review Committee of Royal Prince Alfred Hospital in the Sydney South West Area Health Service, protocol number X02-0137. Written informed consent was obtained from all subjects, and no minors were included. BAL samples were transferred and processed under our University of Illinois at Chicago approved IRB protocol number 2015–0085.
Method and subject recruitment for bronchoalveolar lavage, allergen testing, serum IgE, exhaled nitric oxide (NO), and spirometry has been described previously [23]. Briefly, bronchoscopy took place following fasting for the previous 12 hours. Subjects were sedated with doses of midazolam and alfentanil, based upon the subject’s initial response to treatment. Topical anaesthesia (lidocaine) was introduced through the fiberoptic bronchoscope, with the maximum dose determined by the subject’s weight. Following wedging of the bronchoscope in the sub segmental airways, 4 x 60 ml (i.e. 240ml) aliquots of warmed sterile sodium chloride (0.9% W/V) were instilled, with suctioning after each, into a single, sterile, ice-cold glass reservoir. Pooled bronchoalveolar (BAL) fluid was placed on ice and transported to the laboratory for processing. BAL was centrifuged at 500 x g for 10, and BAL fluid aspirated and stored at -80 until use. Adult volunteers between the ages of 19–30 years were recruited to participate in the study. Asthmatic volunteers with intermittent or mild persistent or moderate persistent atopic asthma (according to Global Initiative for Asthma, Global Strategy for Asthma Management and Prevention, NIH Publication, Updated 2005) were included if they met the following inclusion criteria: asthma symptoms in the preceding 12 months and positive mannitol bronchial provocation test. Individuals were excluded from the study if they had baseline FEV1 less than 60% predicted, inhaled or systemic corticosteroid treatment in the preceding 2 months, smoking in the preceding 12 months, greater than 5 pack year smoking history and symptoms of upper respiratory tract infection in the past 6 weeks.
Study design
Subjects were initially screened at a baseline visit for their medical history, physical examination, exhaled NO, spirometry, mannitol challenge, allergen skin test, and Juniper Asthma Control Questionaire Score. One to four weeks post-intial visit subjects underwent bronchodilator challenge, blood draw for serum IgE, and bronchoscopy with lavage and biopsy. As part of the original Baraket et al. study [23], asthmatic subjects were then randomized into two dosing strategies of inhaled fluticasone propionate (100 μg/BID and 500 μg/BID) delivered by metered dose aerosol and large volume spacer. Subjects were provided with salbutamol metered dose aerosol (100 μg) to be used as needed for asthmatic symptoms. Subjects were on no other medications during the study period. Asthmatic subjects were assessed six weeks after initiating of inhaled fluticasone propionate for exhaled NO, spirometry, mannitol challenge, and Juniper Asthma Control Questionaire Score. Seven weeks after initiating inhaled fluticasone propionate subjects underwent bronchodilator challenge, and bronchoscopy with lavage and biopsy.
Spirometry and lung function reference values
Spirometric measurements were performed with a Vmax spirometer (SensorMedics, Yorba Linda, California) using the closed circuit technique. Subjects were asked to withhold short acting bronchodilators for at least 8 hours prior to the visit. No subjects were taking long acting beta-2 agonists. Subjects were asked to not to consume caffeine or engage in vigorous exercise on the morning of the visit. Predicted values were determined based on race, height, age, and sex using the reference equations derived from NHANES III. Spirometry values were interpreted using the ATS/ERS guidelines [24, 25].
Multiplex cytokine and chemokine assay
100 μL of BAL was centrifuged at 10,000 x g for 10 min to clear any cellular debris. Acellular BAL was assessed for concentrations of sCD40L, MCP-1, MIP-1β, TNF-α, G-CSF, GM-CSF, IFN-γ, IL-1β, -2, -4, -5, -6, -7, -8, -10, -12(p70), -13, -17A, -17F, -21, -22, -23, -25, -31, and -33 using Bio-Plex Pro™ Human Th17 Cytokine Panel 15-plex and Human Cytokine 17-plex assays (Bio-Rad). Concentrations were determined using the Bio-Plex® 200 system and software. IL-5 was excluded from analysis as no samples had concentrations within the detectable range. Samples outside detectable limits were given the concentration of zero.
DNA extraction, metagenomic DNA library construction, and sequencing
BAL was removed of eukaryotic cells by centrifugation at 1000 x g for 10 min. To extract microbial DNA from BAL, 1 mL of BAL was centrifuged at 22,000 x g for 2 hours at 4°C. Supernatants were removed and the pellet was treated with Lysozyme (Fisher Bio) in AL Buffer (Qiagen) to lyse microbial cells. DNA was extracted from lysed pellets using QIAamp Virus Spin Minelute kit (Qiagen) per manufacturer’s protocol. DNA libraries for sequencing were constructed using Nextera XT DNA Library Prep Kit (Illumina). The quality and quantity of all the DNA libraries were analyzed with DNA High Sensitivity Kit (Agilent) on the 2100 Bioanalyzer Instrument (Agilent) and Qubit (Invitrogen). The samples were sequenced as described previously [26]. DNA libraries were pooled in equimolar concentration and were sequenced following manufacturer’s protocol by multiplexing on Illumina MiSeq using the v3-600 kit for 301 paired-end read length and included an additional 6 cycles for the index. The average sequencing depth was 5,065,879 reads per sample.
De novo genomic assembly and annotation of metagenomic sequences
FASTQ files for each sample were mapped against the human genome at a fraction size of 0.500 bp and 95 percent identity to remove contaminating human genomic DNA. There was an average of 116,200 unmapped reads (i.e. non-human or microbial reads) and 18,900,000 base pairs per sample. Unmapped reads were assembled into contiguous sequences at minimum 300 base pair using CLC Genomics Workbench (CLC Bio). The assembly files were uploaded on to the public database server MG-RAST [27] and taxonomy was assigned with the minimum parameters of 60% identity, 10−5 e-value, and 60 amino acids alignment length using MG-RAST’s M5NR database. Data for the samples that were sequenced have been uploaded onto MG-RAST and are available under the Project ID: 12564.
Statistical analysis
All statistical analyses were performed in R (Version 3.3) unless otherwise specified. Cytokines were standardized using z-score and principal component analysis (PCA) was performed using Euclidean distances. Differential centroid testing for PCA was performed by permutation ANOVA (PERMANOVA) using Adonis. For determination of asthmatic phenotypes, unsupervised clustering of the first six principal components was used for hierarchical clustering using Pearson’s correlations as distances and Unweighted Pair Group Method with Arithmetic Mean (UPGMA) as the clustering algorithm. Optimal cluster number was assessed using average silhouette width and cluster stability was determined using consensus clustering [28]. Cytokine and chemokine concentration differences were tested using limma and pair-wise contrasts [29]. Linear regression and ordinal logistic regression for cytokines or modules was performed using univariate or multiple predictors. For multiple adjustment, modules were adjusted for treatment with FP and the interaction of the module with treatment. Modules of cytokines were identified using significant Pearson’s correlations controlled for multiple hypothesis testing (P < 0.05, Q < 0.05) and clustered using complete-linkage algorithm. Clusters were said to be modules if five or more cytokines had interactions with each other to control for clusters composed of spurious interactions. Detrended correspondence analysis and microbiome correlation with phenotypes and clinical variables was performed using vegan. Differential abundance of specific species was tested by pair-wise contrasts using DESeq2 [30]. Slopes of univariate linear regressions with DCA1 as a linear predictor of cytokines or modules were compared using ANCOVA. All subject’s clinical data was tested for significant differences via Student’s T-test for continuous data or Chi-squared analysis with bootstrap approximation for distribution for categorical data. P < 0.05 was considered statistically significant. For multiple testing procedures, Q values were determined to control for Type I error [31]. P < 0.05 and Q < 0.1 was considered significant for multiple testing procedures.
Results
Unsupervised clustering of cytokines and chemokines identifies asthmatic phenotypes
We first assessed markers of inflammation between asthmatics and controls. We analyzed bronchoalveolar lavage (BAL) cytokines and chemokines associated with T-cell responses and found no significant differences between asthmatics and controls using principal coordinate analysis (PCoA) and permutational MANOVA (PERMANOVA) (Fig 1A). We next queried whether sub-classifications of asthmatics were present utilizing unsupervised hierarchical clustering (S1 Fig). Two clusters of asthmatics were identified. We termed these clusters asthma phenotype (AP) 1 and AP2. AP1 and AP2 were found to have high within cluster consensus (AP1: 0.99 [95% CI 0.99–1.00], AP2: 0.97 [95% CI 0.96–0.98] Consensus Index) and low out of cluster consensus (0.04 [95% CI 0.02–0.05] Consensus Index) (Fig 1B). Using the new classification, there was significant separation between controls, AP1, and AP2 utilizing PCoA and PERMANOVA (Pseudo-F test, P = 0.005) (Fig 1C).
(A) Principal coordinate analysis of 24 cytokines and chemokines measured in bronchoalveolar lavage (BAL) between control subjects (red) and asthmatics (blue). (B) Optimal clustering solution (k = 2) shown with heatmap of consensus indices (Legend: 1 = dark blue, 0 = white). (C) Asthmatics were reclassified based on their phenotype, as determined by unsupervised clustering: control (red), asthmatic phenotype (AP)1 (blue), and AP2 (orange). Significance of grouping was determined by PERMANOVA: P < 0.05 was considered significant and R2 is goodness of fit. (D) Clinical markers, exhaled NO (left), bronchial basement membrane thickness (middle), and total serum IgE (right), of asthma severity between controls, AP1, and AP2. Boxplot represent the median and quartiles. Exhaled NO and total serum IgE were log10 transformed for normalization. Significance determined by one-way ANOVA and Tukey’s post-hoc t-test, P < 0.05 was considered significant. BAL cytokines and chemokines were compared between controls and asthmatic phenotypes by pair-wise comparisons using limma. Cytokines and chemokines are shown that were found to be significantly (P < 0.05, Q < 0.1) enriched or depleted in (E) AP1 and (F) AP2. Each group is represented as mean and error bars are standard error of mean.
The immunologic phenotypes, AP1 and AP2, differed in pulmonary function. AP2 was significantly more obstructed with a lower FEV1 (% predicted). More individuals had mild or moderate obstruction compared to AP1 (Table 1). Both phenotypes had positive bronchial provocation tests with AP2 individuals requiring less mannitol to induce airway hyper-reactivity indicating more airway hyper-responsiveness in AP2 individuals. Though there were no significant increases in AP2 compared to both control and AP1, there was a significant linear trend to have higher fractional exhaled nitric oxide (geometric mean (GM) [95% CI]: 7[2–18] vs. 13[7–25] vs. 24[12–47] (ppb), control vs AP1 vs. AP2, P = 0.007), basement membrane thickness (mean [SD]: 8.1[1.0] vs. 9.1[1.2] vs. 9.5[0.7] (μM), control vs AP1 vs. AP2, P = 0.038), and total serum IgE (GM [95% CI]: 51[11–229] vs. 204[60–698] vs. 303[88–1043] (kU/L), control vs AP1 vs. AP2, P = 0.027), suggesting a more severe asthmatic phenotype (Fig 1D). With regards to age, sex, or body mass index (BMI), there were no significant differences between controls, AP1, or AP2 (Table 1).
Significant differences in concentrations of BAL cytokines and chemokines were determined between controls, AP1, and AP2. AP1 had significantly decreased concentrations of G-CSF, IFN-γ, and IL-17F compared to AP2 and controls (Fig 1E, Table A in S1 File). Conversely, AP2 had significantly decreased IL-2 and -10 compared to both controls and AP1. AP2 also displayed increased concentrations of MIP-1β (CCL4) and TNF-α (Fig 1F, Table A in S1 File).
Differentially expressed cytokine associate with airway obstruction
We assessed whether the cytokines and chemokines that were differentially abundant between controls, AP1, and AP2 were also associated with the enhanced airway obstruction noted in AP2. To test for associations with airway obstruction, we constructed models where cytokine or chemokine concentrations were used as linear predictors of FEV1 (% predicted). IL-2 (β [95% CI] = 9.922[1.243–18.600], R2 = 0.255, P = 0.027) and MIP-1β (β [95% CI] = -17.220[-33.480–-0.921], R2 = 0.226, P = 0.039) significantly associated with FEV1 (% predicted) using univariate linear regression (Fig 2A and 2B). To confirm that these were not only correlated with low FEV1 (% predicted) in non-obstructed individuals, ordinal probabilistic models were constructed for no, mild, and moderate obstruction. Both inflammatory mediators, IL-2 (β [95% CI] = -1.841[-3.911–-0.160], P = 0.047) and MIP-1β (β [95% CI] = 3.861[0.897–7.471], P = 0.018) were significantly associated with the type of obstruction (Fig 2C and 2D) and suggest a role for IL-2 and MIP-1β in defining obstructive airway status.
Linear regression of (A) MIP-1β and (B) IL-2 associations with FEV1 (% predicted). Ordinal logistic regression of (C) MIP-1β and (D) IL-2 predicting type of obstruction: none or no obstruction (solid), mild obstruction (FEV1 (% predicted) > 70%) (dashed), moderate obstruction (FEV1 (% predicted) < 70%) (dotted). Solid and dashed lines for all graphs depict trend lines and 95% confidence intervals, respectively. All models were univariate. P < 0.05 was considered significant for the whole model and linear predictor.
Altered cytokine concentrations and reduced obstruction is observed after inhaled corticosteroids
After six weeks of therapy with fluticasone propionate (FP) individuals in both phenotypes had improved lung function with all but one individual in each group demonstrating airflow obstruction (Table B in S1 File). There was no significant difference between AP1 and AP2 with regards to pulmonary function or symptom control after FP treatment. This improvement in lung function corresponded to significantly (Tukey Post-hoc T-test, P < 0.05) increased concentrations of IL-2 and decreased concentrations of MIP-1β, TNF-α, and fractional exhaled nitric oxide in AP2 but not AP1 (Fig 3). There was no significant increase or reduction in IL-10, G-CSF, IFN-γ, or IL-17F (data not shown).
Comparison of IL-2, MIP-1β, TNF-α, and exhaled NO between AP1 and AP2 pre- and post-six weeks of fluticasone propionate (FP) treatment. Comparison was performed using two-way ANOVA and Tukey’s post-hoc t-test, P < 0.05 was considered significant. Cytokines shown as mean and error bars represent standard error of mean. Boxplot represent the median and quartiles. Exhaled NO was log10 transformed for normalization.
Cytokine modules associate with asthmatic phenotypes and pulmonary function
To determine if the cytokines that were differentially abundant between controls and asthmatic phenotypes have effects that are independent or are correlated with other cytokines, we constructed a correlation network of cytokines across all subject samples (n = 30) using significant associations (Pearson’s correlations P < 0.05, Q < 0.05) (S2 Fig). We identified two modules of cytokines, Module 1 and Module 2 (Fig 4A and 4B). Module 1 was comprised of IL-1β, -8, -31, -33, TNF-α, G-CSF, and MIP-1β. Module 2 was comprised of IL-17F, -21, -22, -25, sCD40L, and IFN-γ. Module 1 abundance significantly varied with controls and asthmatic phenotypes pre- and post-FP (ANOVA P = 0.004), whereas Module 2’s variance was not significantly associated with controls and asthmatic phenotypes pre- and post-FP (ANOVA P = 0.349) (Fig 4C and 4D). In addition to Module 1 having increased abundance in AP2 pre-FP, Module 1 had a significant negative association with both FEV1 (% predicted) and FEV1/FVC (Table 2). Module 2 had a significant positive association with exhaled NO. Neither module had associations with basement membrane thickness. When assessing for hub cytokines within each module, we identified both MIP-1β and IFN-γ to have the highest correlation (eigennode significance) with the abundance of their modules and influence over their module (closeness centrality) (Fig 4E and 4F).
Two modules, (A) Module 1 and (B) Module 2, were identified after clustering of cytokines with significant associations. Eigenvalue, linear transformation of cytokine abundances in each module, for (C) Module 1 and (D) Module 2. Hub cytokines were identified by eigennode significance and closeness centrality for (E) Module 1 and (F) Module 2.
The airway microbiome varies between asthmatic phenotypes, non-asthmatics, and after inhaled corticosteroid treatment
We hypothesized that the difference in cytokine profile in BAL was related to microbial composition that could be detected in the BAL. To test this we utilized shotgun sequencing of BAL to determine the composition of the microbiome. We detected 656 taxa at the species level with mean read count of one across all samples. We assessed for differences in α-diversity across samples (S3A Fig). We found no difference between phenotypes and controls. There was a significant increase in diversity after corticosteroid treatment in both phenotypes. To assess between group diversity (β-diversity), taxa counts were normalized using median sum scaling (see Methods) and fit onto detrended correspondence analysis to determine if there was separation of samples along species gradients (S3B Fig). There was significant separation of the phenotypes along the first correspondence axis (DC1) (Table C in S1 File).
To test which taxa were differentially abundant between controls, AP1, and AP2, normalized taxa counts were fit to generalized linear models using the negative binomial distribution (S4 Fig) and sample classification (control, AP1 pre-FP, AP2 pre-FP, AP1 post-FP, and AP2 post-FP) as linear predictors (Table D in S1 File). A trend toward higher error was noted with taxa with low fitted counts, likely due to insufficient modeling of excess zeroes (S4C Fig). We therefore choose to focus our analysis on the top 25 most abundant species, which all had a mean abundance greater that 30 normalized counts. We detected seven species that were differentially abundant (FDR adjusted P < 0.05) between controls, AP1 pre-FP, and AP2 pre-FP (Fig 5A).
(A) β coefficients and 95% CI from negative binomial generalized linear model fits for each differentially abundant species in controls and asthmatic phenotypes pre- and post-FP. (B) LOEWSS regressions of normalized counts of differentially abundant species and DCA1. Univariate linear regressions of (C) IL-2 and (D) Module 1 using DCA1 as a linear predictor. Line type and point color represents controls, asthmatics pre-, or post-FP.
An increased abundance of Streptococcus pneumoniae was found in AP2 pre-FP compared to AP1 pre-FP (Log2 Fold Change (L2FC) [95% CI]: 3.504[1.815–5.192], adjusted P = 0.002, AP2/AP1) and controls (L2FC [95% CI]: 2.909[1.219–4.598], adjusted P = 0.002, AP2/control). There was also a trend toward increased abundance in Neisseria meningitidis in AP2 pre-FP compared to AP1 pre-FP (L2FC [95% CI]: 1.998[0.586–3.409], adjusted P = 0.055, AP2/AP1) and controls (L2FC [95% CI]: 3.216[1.774–4.659], adjusted P = 0.002, AP2/control). The increase in these pathogenic species was not observed in AP1 pre-FP (S. pneumoniae adjusted P = 0.816 and N. meningitidis adjusted P = 0.654). There was decreased abundance of Rothia mucilaginosa, Actinomyces odontolyticus, and Actinomyces sp. oral taxon 180 in AP1 compared to controls and AP2. Whereas, AP2 had decreased abundance of Enterococcus faecalis and faecium compared to AP1 and controls (Tables E-G in S1 File).
After six weeks of FP treatment there were significant reductions in the abundance of S. pneumoniae, N. meningitis, E. faecium, and E. faecalis (Fig 4). AP1 showed a significant reduction in E. faecalis (L2FC [95% CI]: -1.576[-0.854–2.298], P < 0.001, AP1 pre-FP/post-FP) and E. faecalis (L2FC [95% CI]: -1.802[-2.987–4.277], P < 0.001, AP1 pre-FP/post-FP). AP2 had a significant decrease in abundance of S. pneumoniae (L2FC [95% CI]: -3.139[-1.390–4.889], P < 0.001, AP2 pre-FP/post-FP) and N. meningitis (L2FC [95% CI]: -1.708[-0.241–3.174], P = 0.022, AP1 pre-FP/post-FP). There was no significant change in R. mucilaginosa, A. odontolyticus, and A. sp. oral taxon 180 in either phenotype after fluticasone propionate treatment (Tables H-J in S1 File).
The airway microbiome is associated with differential cytokine production observed in asthmatics
To further assess if the differences in the microbiome were related to differences in cytokines we performed multiple regressions with DC1 as a linear predictor for cytokines in controls, asthmatics pre-FP, and asthmatics post-FP. DC1 represents a gradient for the differentially abundant species identified using generalized linear models (Fig 5B). Notably, E. faecalis and E. faecium decrease along the gradient, whereas A. odontolyticus, A. sp. oral taxon 180, R. mucilaginosa, and S. pneumoniae increased along the gradient. The relationship between cytokines and DC1 was significantly different for IL-2, IFN-γ, TNF-α, and Module 1 (ANCOVA P = 0.006, 0.030, 0.002, 0.007, respectively) when comparing controls, asthmatics pre- and post-FP (Table 3). There were also trends in G-CSF, MIP-1β, and Module 2 (ANCOVA P = 0.096, 0.058, respectively). IL-2 had a significant negative association with DC1 in asthmatics pre-FP not seen in controls or asthmatics post-FP. Conversely, Module 1 had significant positive associations with DC1 in asthmatics pre-FP not seen in controls or asthmatics post-FP (Fig 5C and 5D). Differences in microbiome composition, pulmonary function, and cytokine concentrations were not attributable to differences in FP dosing (S5 Fig).
Discussion
Within a homogenous group of atopic asthmatics (i.e. positive skin prick tests, elevated total serum IgE, and elevated exhaled NO), we identified two distinct asthmatic phenotypes (AP) termed AP1 and AP2. These are associated with differences in the local inflammatory milieu, lower airway microbiota composition, and pulmonary function. We identified MIP-1β and IL-2 as cytokines associated with the obstructive status in asthmatics. The abundance of these cytokines was associated with composition of the microbiome. IL-2 dominant profiles were associated with increased Enterococcus species, whereas MIP-1β dominant profiles were associated with Streptococcus pneumoniae and Oral taxa. We observed a reduction in MIP-1β and an increase in IL-2 after inhaled corticosteroid therapy. Additionally, the associations between the microbiome composition and cytokine profiles were diminished after inhaled corticosteroid therapy.
Interestingly, with regards to baseline pulmonary function, AP1 and AP2 resemble cluster 1 and 2 from Moore WC et al. (2010). In this large analysis, asthma phenotypes were identified using unsupervised clustering of asthmatics based on clinical parameters [3]. All individuals in our study had onset of asthma in childhood, positive skin prick tests, and elevated IgE, which is consistent with atopic asthma and the individuals identified in cluster 1 and 2. Similar to cluster 1, AP1 had preserved baseline pulmonary function parameters (FEV1 (% predicted): 102 vs. 90, and FEV1/FVC: 0.78 vs. 0.79, cluster 1 vs. AP1), whereas AP2 moderately reduced pulmonary function baseline more comparable to cluster 2 (FEV1 (% predicted): 82 vs. 76, and FEV1/FVC: 0.74 vs. 0.73, cluster 2 vs. AP2). Importantly, this suggests a link between pulmonary function and the immunology of the airway as we had defined these clusters based on immunologic profiles. In asthmatics prior to treatment the immunologic profiles were highly associated with the composition of the microbiota in airway. These data indicate that in atopic asthmatics differences in pulmonary function may be mediated by unique immune responses to the airway microbiome.
AP1 exhibited significantly decreased concentrations of cytokines (IFN-γ, IL-17F, and G-CSF) associated with microbial pathogens [12, 18]. In contrast, AP2 had decreased concentrations of cytokines associated with reduced airway inflammation (IL-2 and IL-10) [9, 32, 33]. Additionally, AP2 exhibited higher concentrations of TNF-α and MIP-1β. Notably, the concentrations of IL-2 and IL-10 were independent of other cytokines. In contrast, TNF-α and MIP-1β were highly correlated with IL-1β, -8, -31, and -33, which we identified as Module 1. MIP-1β was identified as a hub in Module 1. The pro-inflammatory Module 1 was associated with a greater likelihood of airway obstruction, negatively correlating FEV1 (% predicted) and FEV1/FVC. Importantly, IL-8, -31, and -33 have been shown to negatively correlate in independent associations with FEV1 (% predicted) [34–36]. The abundance of Module 1 was reduced by treatment with the inhaled corticosteroid, fluticasone propionate. Future investigations may identify specific cell populations secreting cytokines related to changes in the microbiome.
The composition of the airway microbiome also significantly differed between asthmatic phenotypes. AP2 was enriched with species from the upper airways, Actinomyces ondontolyticus, Actinomyces oral taxon 180, Neisseria menigiditis, and Streptococcus pneumoniae, whereas AP1 was enriched with Enterococcus species. These species varied along the constructed ecological gradient DCA1. Similar to previous reports, these data suggest that the lung microbiome can vary with regards to abundance of upper airway taxon [17, 18]. Interestingly, asthmatic individuals may have impaired mucocilliary clearance upon allergen exposure and therefore might have different rates of elimination of upper-airway species [37]. Unlike what has been previously reported in non-asthmatic individuals, in asthmatics prior to FP treatment the increasing abundance of these upper airway species was not associated with Th17 responses but was highly associated with decreasing IL-2 and increasing Module 1.
Interestingly, the abundance of Streptococcus pneumoniae was reduced after FP treatment which correlated with the reduction in concentration of Module 1 (IL-1β, -8, -31, -33, TNF-α and MIP-1β). Previously, the abundance of Streptococcus in sputum has been shown to positively correlate with IL-8 and neutrophils [35]. Actinomyces and Neisseria species have been shown to be enriched in asthmatics but also increased in abundance with lower BAL eosinophils [20]. In our study, DCA1 represented an increasing gradient of Actinomyces ondontolyticus, Actinomyces oral taxon 180, Neisseria menigiditis, and Streptococcus pneumoniae. DCA1 had unique associations with IL-2 and Module 1 in asthmatics prior to FP treatment. In murine infection models, Streptococcus pneumoniae induces increased production of many of the cytokines in Module 1 [38].
Our study, though limited in sample size and racial diversity, highlights the heterogeneous nature of asthma. Even within asthmatics that fit the picture of atopic asthma, there appears to significant variation within local inflammation and the lower airway microbiome, thus allowing us to discern two distinct phenotypes. Our observation may have been strengthened by the strong exclusion criteria: no history of smoking in the past year, less than a five pack year history of smoking, and no symptoms of respiratory infection in the proceeding 6 weeks. Although the inclusion and exclusion criteria control for several potential confounders, the small sample size of this study limits the data to be less generalizable to all asthmatics and effects covariates associated with disease severity might be missed, such as BMI which made up 50% of AP1 subjects. In addition to larger sample sizes, future investigations of smokers and those after an infection or asthma exacerbation would be of interest to test for subtypes similar to AP1 and AP2.
In addition to small sample size, our analysis did not include assessment of BAL cellular profiles (e.g. % of eosinophils in BAL) nor did it include a placebo arm for treatment comparison. It should be noted that BAL cellular profile are associated with differences in the airway microbiome [39, 40]. Notably, in these studies the taxa Streptoccocus, Neisseria, Actinomyces, and Rothia were identified as differentially abundant between asthmatics with different granulocytic profiles. The lack of placebo arm reduces conclusions that can be drawn after steroid treatment, as it is not known whether the changes in the microbiome or inflammation are due to random temporal exposures to bacteria. Due to clinical limitations of obtaining serial samples from the same individuals, this area of literature is lacking. The association between the exposure to different microbiota and the detected inflammation was not significant after treatment (Table 3). These data suggest that even if the bacterial variation is due random temporal variation, ICS therapy may reduce the inflammation associated with exposure to specific bacterial species, such as those from the oral cavity.
Currently, more personalized approaches are being developed for treatment of the most severe forms of asthma. Our study indicates that even in mild to moderate asthmatics with intermittent or persistent symptom there are distinct sub-phenotypes (AP1 and AP2) that differ in their immune response and microbiome. These sub-phenotypes are associated with important aspects of airway physiology. Additionally, monitoring for microbiological changes, such as an abundance of Streptococcus pneumoniae and oral-pharyngeal taxa, may allow for personalized approaches for patient monitoring.
Supporting information
S1 Fig. Hierarchical clustering of asthmatic subjects pre-fluticasone propionate identifies two clusters.
Bronchoalveolar lavage (BAL) cytokines and chemokines from asthmatic subjects were standardized using z-score and sample distances were determined using Euclidean distances. Sample wise distance matrix was used for metric multidimensional scaling and (A) a Scree plot of the components was used to determine how many components to use in the clustering algorithm. The first six components were selected (containing > 90% of all sample variance) and underwent hierarchical clustering using Pearson’s correlation’s as distances and Unweighted Pair Group Method with Arithmetic Mean (UPGMA) as the clustering algorithm. (B) Average silhouette width for each cluster solution shown with cross hairs showing optimal cluster solution. (C) Consensus clustering was used to determine stability of various clustering solution k = 2–8 (colors) and cumulative distribution functions of consensus indices for each solution are shown.
https://doi.org/10.1371/journal.pone.0184566.s001
(TIF)
S2 Fig. Correlation network of cytokines.
(A) Heatmap showing significant (P < 0.05 and Q < 0.05) absolute Pearson’s correlation between cytokines. Cytokines were clustered using the complete linkage algorithm. (B) Correlation network of cytokines with significant associations with each other. Modules are highlighted by color, Module 1 (blue) and Module 2 (red). Size of edge is proportional to strength of association. Straight lines signify positive correlations and sine wave lines signify negative correlations.
https://doi.org/10.1371/journal.pone.0184566.s002
(TIF)
S3 Fig. α- and β-Diversity of broncho-alveolar lavage microbiome.
(A) Rarefaction curves and (B) Detrendend correspondence analysis of 656 microbial species showing the first two scaled ordination axes, DCA1 and DCA2, which account for 50.8% and 25.9% of variance between all samples. Colors represent groups: control (light blue), AP1 pre-FP (black), AP2 pre-FP (green), AP1 post-FP (red), and AP2 post-FP (blue). Ellipses are standard deviation from centroids of each group. Significant factors and vectors plotted as text in red and determined centroids. Significance determined by R2 and degrees of freedom, P < 0.05 was considered significant.
https://doi.org/10.1371/journal.pone.0184566.s003
(TIF)
S4 Fig. Negative Binomial for modeling taxa counts.
(A) Dispersion parameter for individual taxa (blue) is determined by initial estimation (black) and Bayesian shrinkage using the prior knowledge of local average regression of dispersion with mean of normalized count as the predictor (trend line, red). (B) Fitted counts plotted against observed counts with a pseudo-counts of one added to each to visualize zeros. (C) Residuals plot showing error vs. fitted counts with counts greater than 30 displaying a more normally distributed error. (D) Q-Q plot of observed quantiles versus the theoretical quantiles of random data generated using the negative binomial distribution and parameters for each taxa, NB(μ = mean, α = dispersion).
https://doi.org/10.1371/journal.pone.0184566.s004
(TIF)
S5 Fig. No difference in pulmonary function, cytokines, or microbiome between low and high dose fluticasone propionate.
Boxplots of FEV1 (% predicted), FEV1/FVC, Module 1, IL-2, and DCA1 between individuals randomized to low (100 μg/BID) or high dose (500 μg/BID). P < 0.05 was considered significant, Welch’s T-test.
https://doi.org/10.1371/journal.pone.0184566.s005
(TIF)
S1 File.
Supplemental tables A-J. Supplemental tables A-J.
https://doi.org/10.1371/journal.pone.0184566.s006
(DOCX)
Acknowledgments
The authors would like to thank Drs. Dean Schraufnagel and Thomas Ferkol for their helpful comments regarding this manuscript.
References
- 1. Miranda C, Busacker A, Balzar S, Trudeau J, Wenzel SE. Distinguishing severe asthma phenotypes: role of age at onset and eosinophilic inflammation. J Allergy Clin Immunol. 2004 January;113(1):101–8. pmid:14713914
- 2. Haldar P, Pavord ID, Shaw DE, Berry MA, Thomas M, Brightling CE, et al. Cluster Analysis and Clinical Asthma Phenotypes. Am J Respir Crit Care Med. 2008 August 1,;178(3):218–24. pmid:18480428
- 3. Moore WC, Meyers DA, Wenzel SE, Teague WG, Li H, Li X, et al. Identification of Asthma Phenotypes Using Cluster Analysis in the Severe Asthma Research Program. Am J Respir Crit Care Med. 2010 February 15,;181(4):315–23. pmid:19892860
- 4. Thorburn AN, McKenzie CI, Shen S, Stanley D, Macia L, Mason LJ, et al. Evidence that asthma is a developmental origin disease influenced by maternal diet and bacterial metabolites. Nat Commun. 2015;6:7320. pmid:26102221
- 5. Schaub B, Campo M, He H, Perkins D, Gillman MW, Gold DR, et al. Neonatal immune responses to TLR2 stimulation: Influence of maternal atopy on Foxp3 and IL-10 expression. Respir Res. 2006;7(1):40.
- 6. Irvin C, Zafar I, Good J, Rollins D, Christianson C, Gorska MM, et al. Increased frequency of dual-positive TH2/TH17 cells in bronchoalveolar lavage fluid characterizes a population of patients with severe asthma. J Allergy Clin Immunol. 2014 November;134(5):1186.e7.
- 7. Singh A, Yamamoto M, Ruan J, Choi JY, Gauvreau GM, Olek S, et al. Th17/Treg ratio derived using DNA methylation analysis is associated with the late phase asthmatic response. Allergy Asthma Clin Immunol. 2014;10(1):32. pmid:24991220
- 8. Choy DF, Hart KM, Borthwick LA, Shikotra A, Nagarkar DR, Siddiqui S, et al. TH2 and TH17 inflammatory pathways are reciprocally regulated in asthma. Sci Transl Med. 2015 August 19,;7(301):301ra129. pmid:26290411
- 9. Bohm L, Maxeiner J, Meyer-Martin H, Reuter S, Finotto S, Klein M, et al. IL-10 and Regulatory T Cells Cooperate in Allergen-Specific Immunotherapy To Ameliorate Allergic Asthma. Journal of Immunology. 2015 Feb 1,;194(3):887–97.
- 10. Huang YJ, Nelson CE, Brodie EL, DeSantis TZ, Baek MS, Liu J, et al. Airway microbiota and bronchial hyperresponsiveness in patients with suboptimally controlled asthma. Journal of Allergy and Clinical Immunology. 2011 February;127(2):381.e3.
- 11. Denner DR, Sangwan N, Becker JB, Hogarth DK, Oldham J, Castillo J, et al. Corticosteroid therapy and airflow obstruction influence the bronchial microbiome, which is distinct from that of bronchoalveolar lavage in asthmatic airways. J Allergy Clin Immunol. 2015 November 25,.
- 12. Huang YJ, Nariya S, Harris JM, Lynch SV, Choy DF, Arron JR, et al. The airway microbiome in patients with severe asthma: Associations with disease features and severity. Journal of Allergy and Clinical Immunology. 2015 October;136(4):874–84. pmid:26220531
- 13. Goleva E, Jackson LP, Harris JK, Robertson CE, Sutherland ER, Hall CF, et al. The effects of airway microbiome on corticosteroid responsiveness in asthma. Am J Respir Crit Care Med. 2013 November 15,;188(10):1193–201. pmid:24024497
- 14. Kloepfer KM, Lee WM, Pappas TE, Kang TJ, Vrtis RF, Evans MD, et al. Detection of pathogenic bacteria during rhinovirus infection is associated with increased respiratory symptoms and asthma exacerbations. J Allergy Clin Immunol. 2014 May;133(5):3.
- 15. Teo S, Mok D, Pham K, Kusel M, Serralha M, Troy N, et al. The Infant Nasopharyngeal Microbiome Impacts Severity of Lower Respiratory Infection and Risk of Asthma Development. Cell Host & Microbe. 2015 May 13,;17(5):704–15.
- 16. Bisgaard H, Hermansen MN, Buchvald F, Loland L, Halkjaer LB, Bønnelykke K, et al. Childhood Asthma after Bacterial Colonization of the Airway in Neonates. New England Journal of Medicine. 2007 October 11,;357(15):1487–95. pmid:17928596
- 17. Segal LN, Alekseyenko AV, Clemente JC, Kulkarni R, Wu B, Chen H, et al. Enrichment of lung microbiome with supraglottic taxa is associated with increased pulmonary inflammation. Microbiome. 2013 -7-1;1:19. pmid:24450871
- 18. Segal LN, Clemente JC, Tsay JJ, Koralov SB, Keller BC, Wu BG, et al. Enrichment of the lung microbiome with oral taxa is associated with lung inflammation of a Th17 phenotype. Nat Microbiol. 2016 -4-04;1:16031. pmid:27572644
- 19. Durack J, Lynch S, Nariya S, Bhakta N, Beigelman A, Castro M, et al. Feature of the bronchial bacterial microbiome associate with atopy, asthma, and responsiveness to inhaled corticosteroid treatment. Journal of Allergy and Clinical Immunology. 2015 Jan 1,;135(1):25–30.
- 20. Sverrild A, Kiilerich P, Brejnrod A, Pedersen R, Porsbjerg C, Bergqvist A, et al. Eosinophilic airway inflammation in asthmatic patients is associated with an altered airway microbiome. Journal of Allergy and Clinical Immunology. 2016 Dec.
- 21. Einarsson GG, Comer DM, McIlreavey L, Parkhill J, Ennis M, Tunney MM, et al. Community dynamics and the lower airway microbiota in stable chronic obstructive pulmonary disease, smokers and healthy non-smokers. Thorax. 2016 Sep;71(9):795–803. pmid:27146202
- 22. Wang Z, Bafadhel M, Haldar K, Spivak A, Mayhew D, Miller BE, et al. Lung microbiome dynamics in COPD exacerbations. The European respiratory journal. 2016 Apr;47(4):1082–92. pmid:26917613
- 23. Baraket M, Oliver BG, Burgess JK, Lim S, King GG, Black JL. Is low dose inhaled corticosteroid therapy as effective for inflammation and remodeling in asthma? A randomized, parallel group study. Respir Res. 2012;13(1):11.
- 24. Miller MR, Hankinson J, Brusasco V, Burgos F, Casaburi R, Coates A, et al. Standardisation of spirometry. European Respiratory Journal. 2005 Aug 1,;26(2):319–38. pmid:16055882
- 25. Pellegrino R, Viegi G, Brusasco V, Crapo RO, Burgos F, Casaburi R, et al. Interpretative strategies for lung function tests. The European respiratory journal. 2005 Nov;26(5):948–68. pmid:16264058
- 26. Ranjan R, Rani A, Metwally A, McGee HS, Perkins DL. Analysis of the microbiome: Advantages of whole genome shotgun versus 16S amplicon sequencing. Biochemical and biophysical research communications. 2016 January 22,;469(469):967–77.
- 27. Meyer F, Paarmann D, D'Souza M, Olson R, Glass EM, Kubal M, et al. The metagenomics RAST server–a public resource for the automatic phylogenetic and functional analysis of metagenomes. BMC Bioinformatics. 2008 -9-19;9:386. pmid:18803844
- 28. Wilkerson M, Hayes DN. ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking. Bioinformatics. 2010 28 April;26(12):1572–3. pmid:20427518
- 29. Smyth GK. Linear Models and Empirical Bayes Methods for Assessing Differential Expression in Microarray Experiments. Statistical Applications in Genetics and Molecular Biology. 2004;3(1):1–25.
- 30. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biology. 2014 December 5,;15(12):550. pmid:25516281
- 31. Storey JD, Tibshirani R. Statistical significance for genomewide studies. Proc Natl Acad Sci U S A. 2003 -8-5;100(16):9440–5. pmid:12883005
- 32. Wu K, Ma J, Bai W, Cui X, Han T, Wang S, et al. Short-term intratracheal use of PEG-modified IL-2 and glucocorticoid persistently alleviates asthma in a mouse model. Sci Rep. 2016;6:31562. pmid:27527926
- 33. McGee HS, Agrawal DK. Naturally Occurring and Inducible T-Regulatory Cells Modulating Immune Response in Allergic Asthma. Am J Respir Crit Care Med. 2009 -8-1;180(3):211–25. pmid:19447898
- 34. Lai T, Wu D, Li W, Chen M, Yi Z, Huang D, et al. Interleukin-31 expression and relation to disease severity in human asthma. Scientific reports. 2016;6:22835. pmid:26956917
- 35. Green BJ, Wiriyachaiporn S, Grainge C, Rogers GB, Kehagia V, Lau L, et al. Potentially Pathogenic Airway Bacteria and Neutrophilic Inflammation in Treatment Resistant Severe Asthma. PLoS One. 2014 June 23,;9(6).
- 36. Albacker LA, Chaudhary V, Chang Y, Kim HY, Chuang Y, Pichavant M, et al. Invariant natural killer T cells recognize a fungal glycosphingolipid that can induce airway hyperreactivity. Nature Medicine. 2013 Oct 1,;19(10):1297–304. pmid:23995283
- 37. Bennett WD, Herbst M, Alexis NE, Zeman KL, Wu J, Hernandez ML, et al. Effect of inhaled dust mite allergen on regional particle deposition and mucociliary clearance in allergic asthmatics. Clinical & Experimental Allergy. 2011 Dec;41(12):1719–28.
- 38. Williams AE, José RJ, Brown JS, Chambers RC. Enhanced inflammation in aged mice following infection with Streptococcus pneumoniae is associated with decreased IL-10 and augmented chemokine production. Am J Physiol Lung Cell Mol Physiol. 2015 March 15,;308(6):539.
- 39. Sverrild A, Kiilerich P, Brejnrod A, Pedersen R, Porsbjerg C, Bergqvist A, et al. Eosinophilic airway inflammation in asthmatic patients is associated with an altered airway microbiome. Journal of Allergy and Clinical Immunology. 2016 Dec.
- 40. Taylor SL, Leong LE, Choo JM, Wesselingh S, Yang IA, Upham JW, et al. Inflammatory phenotypes in patients with severe asthma are associated with distinct airway microbiology. Journal of Allergy and Clinical Immunology. 2017 -05-04(0091–6749, 1097–6825).