In lung cancer, enrichment of the lower airway microbiota with oral commensals commonly occurs, and ex vivo models support that some of these bacteria can trigger host transcriptomic signatures associated with carcinogenesis. Here, we show that this lower airway dysbiotic signature was more prevalent in the stage IIIB–IV tumor–node–metastasis lung cancer group and is associated with poor prognosis, as shown by decreased survival among subjects with early-stage disease (I–IIIA) and worse tumor progression as measured by RECIST scores among subjects with stage IIIB–IV disease. In addition, this lower airway microbiota signature was associated with upregulation of the IL17, PI3K, MAPK, and ERK pathways in airway transcriptome, and we identified Veillonella parvula as the most abundant taxon driving this association. In a KP lung cancer model, lower airway dysbiosis with V. parvula led to decreased survival, increased tumor burden, IL17 inflammatory phenotype, and activation of checkpoint inhibitor markers.

Significance:

Multiple lines of investigation have shown that the gut microbiota affects host immune response to immunotherapy in cancer. Here, we support that the local airway microbiota modulates the host immune tone in lung cancer, affecting tumor progression and prognosis.

See related commentary by Zitvogel and Kroemer, p. 224.

This article is highlighted in the In This Issue feature, p. 211

Lung cancer has remained the leading cause of cancer death worldwide. In this past year alone, lung cancer occurred in approximately 2.1 million patients and was responsible for 1.7 million deaths (1). Targeting certain somatic mutations has improved survival, but this is applicable to only ∼30% of subjects with lung adenocarcinoma (2, 3). More recently, immunotherapy that targets inhibitory checkpoint molecules, such as programmed death 1 (PD-1), has been shown to affect the responses of T cells to neoantigens and improve survival in lung cancer (4–8). However, 40% to 60% of patients will not respond to or will develop resistance to immunotherapy (7). Recent investigations have identified gut microbiota signatures that are associated with augmenting antitumor immunity and responding to PD-1 blockade in murine models and in prospective analyses of immunotherapy-responsive cancer cohorts (9–11). For example, modulation of the microbiota in germ-free mice can enhance antitumor immunity and augment effects of checkpoint blockade (12, 13). Matson and colleagues found that in patients with melanoma, anti–PD-1 treatment responders had a higher abundance of B. longum, C. aerofaciens, and E. faecium compared with nonresponders (11). Gopalakrishnan and colleagues demonstrated that patients with higher bacterial diversity and increased relative abundance of Ruminococcaceae in the gut had enhanced systemic and antitumor immune responses (10). Routy and colleagues identified that the relative abundance of A. muciniphila was associated with a favorable clinical response to immunotherapy (9). Although most investigations have focused on the gut microbiome, no human studies have studied the lower airway microbiota and lung cancer prognosis despite growing evidence supporting the role of the lung microbiota in lower airway inflammation (14–16).

Our understanding of the role of lung microbiota in health and disease is rapidly evolving with evidence that some phenotypic characteristics of the local lung immune tone appear to be more closely correlated to the lung microbiome than to the gut microbiome (14). Culture-independent techniques show that the lower airways of normal individuals commonly harbor oral bacteria such as Prevotella and Veillonella (15, 17–19). Our group has described that lower airway dysbiosis characterized by enrichment with oral commensals is associated with increased host inflammatory tone in the lung of healthy individuals (15, 19). This same lower airway dysbiotic signature was found to differentiate between subjects with lung cancer and subjects with benign lung nodules (16). Importantly, we have shown in humans and in ex vivo experimental models that this dysbiotic signature likely triggers transcriptomic signatures (PI3K and MAPK) previously described in non–small cell lung cancer (NSCLC; refs. 16, 20), including the p53 mutation pathway (21). In order to explore the clinical implications of the lower airway microbiota in lung cancer, we utilized a prospective human cohort and a preclinical model to identify lower airway dysbiotic signatures that may affect the prognosis in this disease.

Lung Cancer Cohort

Between March 2013 and October 2018, we recruited 148 subjects with lung nodules from the NYU Lung Cancer Biomarker Center who underwent clinical bronchoscopy for diagnostic purposes and in whom lower airway brushes were obtained for research (Supplementary Fig. S1). Fifteen subjects had non–lung primary tumors (metastasis), 12 had benign lung nodules, and 38 subjects had other nonmalignant diagnosis and were excluded. The remaining 83 subjects had a final diagnosis of lung cancer and were included for this project. Among these subjects, all had microbiome 16S rRNA gene-sequencing data, 70/83 had transcriptomic data, and 75/83 had greater than six months of follow-up clinical data. Supplementary Table S1 describes the demographics and clinical characteristics of this cohort: 91% were current or former smokers with a mean history of 46 pack-years. Eighty-nine percent had a diagnosis of NSCLC, of which 65% had adenocarcinoma and 49% were found to have stage IIIB–IV. The median survival was 2.1 years; 54% received chemotherapy, 30% received radiotherapy, 24% received surgery, and 14% received immunotherapy. All biospecimens were obtained prior to treatment. Using the Cox proportional hazards model, we determined that surgical treatment and stage IIIB–IV were significantly associated with overall survival (OS; Supplementary Table S2).

Microbiomic Signatures Associated with Stage and Prognosis

In addition to lower airway brushings, we obtained buccal brushes and bronchoscope background control samples that were included in the 16S rRNA gene-sequencing analysis. As compared with background controls, the bacterial load was ∼10 times higher in lower airway brushing samples and ∼10,000 times higher in the upper airways (buccal; P < 0.001; Supplementary Fig. S2). Alpha diversity based on the Shannon index showed greater diversity among lower airway samples than upper airway and background control samples (P < 0.001; Supplementary Fig. S3A). Principal coordinate analysis (PCoA) based on Bray–Curtis dissimilarity index showed significant compositional differences across sample types (Fig. 1A; PERMANOVA P < 0.001). Across lower airway samples, there were also compositional differences between small cell lung cancer (SCLC) and NSCLC (PERMANOVA P = 0.01). Among NSCLC samples, there were no statistically significant differences in α-diversity and β-diversity between squamous cell carcinoma and adenocarcinoma.

Figure 1.

Lung microbiota in lung cancer and cancer survival. A, PCoA of airway samples shows a difference in β-diversity (P = 0.01, PERMANOVA) between small cell lung cancer and NSCLC (n = 83). B, Among patients with NSCLC (n = 74), PCoA shows a difference in β-diversity (P = 0.005, PERMANOVA) between stages IIIB–IV and I–IIIA NSCLC (left); lower airway microbiota of stages IIIB–IV were more similar to buccal microbiota than lower airway microbiota of stages I–IIIA (right, P < 0.0001, Bray–Curtis distance). C, Left, PCoA based on cancer stage and survival at six months and one year showsdifference in β-diversity (P < 0.05, PERMANOVA). C, Right, lower airway microbiota in lung cancer and worse survival at six months or one year were more similar to buccal microbiota than with better survival in both stages IIIB–IV (n = 36) and I–IIIA (n = 37) groups (P < 0.05, Bray–Curtis distance). D, Stage IIIB–IV lung cancer was associated with having a higher proportion of subjects whose lower airway microbiota were classified as enriched with oral taxa (supraglottic predominant taxa, SPT) vs. background taxa (background predominant taxa, BPT), P = 0.006. E, Enrichment of the lower airway with pneumotypeSPT was associated with better survival in stage I–IIIA cancer than enrichment with pneumotypeBPT, P < 0.05; there was no difference in stage IIIB–IV cancer. F, Bray–Curtis dissimilarity index between lower airway and buccal samples was inversely associated with delta RECIST score for stage IIIB–IV NSCLC measured at 6 to 12 months (Spearman r = −0.48, P = 0.03).

Figure 1.

Lung microbiota in lung cancer and cancer survival. A, PCoA of airway samples shows a difference in β-diversity (P = 0.01, PERMANOVA) between small cell lung cancer and NSCLC (n = 83). B, Among patients with NSCLC (n = 74), PCoA shows a difference in β-diversity (P = 0.005, PERMANOVA) between stages IIIB–IV and I–IIIA NSCLC (left); lower airway microbiota of stages IIIB–IV were more similar to buccal microbiota than lower airway microbiota of stages I–IIIA (right, P < 0.0001, Bray–Curtis distance). C, Left, PCoA based on cancer stage and survival at six months and one year showsdifference in β-diversity (P < 0.05, PERMANOVA). C, Right, lower airway microbiota in lung cancer and worse survival at six months or one year were more similar to buccal microbiota than with better survival in both stages IIIB–IV (n = 36) and I–IIIA (n = 37) groups (P < 0.05, Bray–Curtis distance). D, Stage IIIB–IV lung cancer was associated with having a higher proportion of subjects whose lower airway microbiota were classified as enriched with oral taxa (supraglottic predominant taxa, SPT) vs. background taxa (background predominant taxa, BPT), P = 0.006. E, Enrichment of the lower airway with pneumotypeSPT was associated with better survival in stage I–IIIA cancer than enrichment with pneumotypeBPT, P < 0.05; there was no difference in stage IIIB–IV cancer. F, Bray–Curtis dissimilarity index between lower airway and buccal samples was inversely associated with delta RECIST score for stage IIIB–IV NSCLC measured at 6 to 12 months (Spearman r = −0.48, P = 0.03).

Close modal

We then evaluated microbial differences in lower airway samples based on the clinical NSCLC stage, grouped as I–IIIA and IIIB–IV of tumor–node–metastasis (TNM) classification. The selection of this cutoff point for TNM classification allowed for dichotomized analyses, and we support it based on prior prognosis/survival data and cancer management guidelines related to surgical management of NSCLC (22–25). Alpha diversity was similar across staging groups of NSCLC (Supplementary Fig. S3B for comparison across individual stages and Supplementary Fig. S3C for two-group comparisons of stages IIIB–IV vs. I–IIIA). Compositional differences between the I–IIIA and IIIB–IV lung cancer groups were noted based on β-diversity analysis (Fig. 1B, left; P = 0.005), where stage IIIB–IV lung cancer was compositionally more similar to buccal samples than I–IIIA stage lung cancer samples (Fig. 1B, right). Compositional differences comparing all individual stages (I–IV) were also noted based on β-diversity analysis (Supplementary Fig. S4A; P = 0.047), where lower airway samples from more advanced stages had a greater similarity to buccal samples than lower airway samples from earlier-stage subjects (Supplementary Fig. S4B). Microbiome regression-based kernel association test (MiRKAT) analysis showed that differences noted in the microbial community profiles between stage I–IIIA and IIIB–IV NSCLC were not due to differences in location of the samples. Interestingly, subanalysis on patient samples where tumor PD-L1 expression was available (n = 39) shows that subjects with high PD-L1 expression (≥80%, n = 12) had a lower airway microbiota with greater similarity to upper airway microbiota versus the disease of similarity found among patients with lower tumor PD-L1 expression (0%, n = 16; 1%–79%, n = 11; P < 0.05; Supplementary Fig. S5).

Compositional differences based on six-month and one-year survival were also identified in β-diversity analysis (Fig. 1C, left), where samples from subjects with decreased survival were associated with greater compositional similarity to buccal samples than samples from subjects with better outcomes (Fig. 1C, right). Shannon index showed decreased α-diversity among samples from subjects with <6-month survival in both stages I–IIIA and IIIB–IV, but this difference was not statistically significant at one year (Supplementary Fig. S6). Multivariate PERMANOVA analysis demonstrated that the association between microbial community composition and six-month/one-year mortality was independent of TNM staging (Supplementary Fig. S7). No statistically significant differences were noted in α- or β-diversity analyses of buccal microbiota between subjects with different stages or mortality.

DESeq analyses were then performed to evaluate for taxonomic differential enrichment between SCLC and NSCLC and between the I–IIIA and IIIB–IV groups of NSCLC (Supplementary Fig. S8A). Importantly, lower airway samples from patients in the IIIB–IV stage group were enriched with many operational taxonomic units (OTU), which annotated to the genera Moraxella, Fusobacterium, Pseudomonas, and Haemophilus, and were decreased in abundance of Actinomycetales (Supplementary Fig. S8B; Supplementary Table S1). Using a mixed-effect model that adjusts for sample location, we report the top 20 OTUs ranked by their absolute coefficients estimates as having a differential abundance between the stage I–IIIA and IIIB–IV groups (Supplementary Table S3). Once again, stage IIIB–IV lung cancer was enriched with OTUs recognized as oral commensals, such as Haemophilus, Fusobacterium, Gemella, Prevotella, and Granulicatella.

Among stage I–IIIA and IIIB–IV subgroups, multiple OTUs were differentially enriched when worse versus better survival groups were compared (both at 6 and 12 months). Several of the OTUs annotated to the genera Veillonella, Prevotella, and Streptococcus were found to be enriched in samples from subjects with worse prognosis (Supplementary Fig. S9A–S9D; Supplementary Tables S2–S5). In order to further explore for taxonomic associations with mortality while considering TNM staging, we constructed β-diversity biplots that allowed for colocation of lower airway samples and taxa driving the spatial distribution. Using a multivariate analysis adjusted by TNM stage, Supplementary Fig. S10 shows that poor prognosis was associated with enrichment of the lower airway microbiota with oral commensals (such as Streptococcus, Prevotella, and Veillonella). When analysis was repeated only considering the lower airway samples with closest proximity to the cancer, similar results were found (Supplementary Fig. S11). Using a mixed-effect model adjusted by smoking status, stage (I–IIIA/IIIB–IV), and treatment type, we identified top OTUs associated with OS. Supplementary Table S4 reports the top 20 OTUs ranked by absolute coefficient estimates associated with OS. Poor prognosis was associated with enrichment with OTUs recognized as oral commensals that belong to the genera Prevotella, Streptococcus, Lactobacillus, and Gemella.

Utilizing a Dirichlet multinomial mixture (DMM) model, we established that samples can be divided into two clusters: cluster one consists of all the upper airway samples and ∼60% of lower airway samples and cluster two consists of all the bronchoscope background control samples and ∼40% of the lower airway samples (Supplementary Fig. S12A and S12B). Thus, similar to previously published data (15), our cohort consists of one cluster of lower airway samples enriched with background predominant taxa (BPT), such as Flavobacterium and Pseudomonas, whereas the second cluster was enriched with supraglottic predominant taxa (SPT), such as Veillonella, Streptococcus, Prevotella, and Haemophilus (Supplementary Fig. S12C; Supplementary Table S6). Supplementary Table S5 shows that we did not identify statistically significant differences in demographic or clinical characteristics, other than stage IV TNM staging (P < 0.05), between subjects with a lower airway microbiota that clustered as BPT versus SPT. Applying decontam (26) approach to these data, an analytic pipeline that accounts for taxa most likely to be contaminants, we identified Flavobacterium as a background contaminant (also most prevalent and abundant OTU in background controls) and oral commensals, such as Veillonella and Streptococcus, as most representatives of lower airway microbiota (Supplementary Fig. S13).

We then used the DMM grouping to evaluate whether the prevalence of SPT/BPT was different among stage I–IIIA and IIIB–IV NSCLC and/or associated with prognosis. The percentage of SPT was higher in lower airway samples from subjects from the stage IIIB–IV NSCLC group compared with lower airway samples from the stage I–IIIA NSCLC group (Fig. 1D; P = 0.006). Importantly, the Kaplan–Meier survival analysis shows that among subjects with stage I–IIIA NSCLC, the SPT pneumotype was associated with worse survival than the BPT pneumotype (Fig. 1E, P = 0.047). In stage IIIB–IV NSCLC, there were no statistically significant differences in survival between the SPT and BPT pneumotypes, although the overall mortality was much worse, with a median survival of less than one year as found in the above analysis. To further evaluate microbial signatures associated with treatment response, we analyzed a subset of patients with stage IIIB–IV NSCLC (thus nonsurgical) with available longitudinal imaging which allowed us to calculate the Response Evaluation Criteria in Solid Tumors (RECIST; ref. 27). Correlation analysis between delta RECIST score and β-diversity dissimilarity between upper and lower airways showed a significant inverse correlation (Fig. 1F; Spearman r = −0.48, P = 0.03). Thus, although overall mortality was not associated with pneumotype categorization in the stage IIIB–IV group, having a positive delta RECIST score, indicating tumor progression, was associated with having a lower airway microbiota more similar to that of upper airways. Taxonomic differences between a dichotomized RECIST score showed lower airway samples from patients with tumor progression (RECIST = progressive disease or stable disease) were enriched with Veillonella, Streptococcus, Prevotella, and Rothia when compared with lower airway samples from patients with tumor regression (RECIST = complete response or partial response; Supplementary Fig. S14 and Supplementary Table S7).

Transcriptomic Signatures Associated with Stage, Prognosis, and Microbiota

After quality control, RNA-sequencing (RNA-seq) data were obtained on 70 lower airway samples from 70 subjects with NSCLC. We then compared global transcriptomic differences between stage I–IIIA and IIIB–IV NSCLC with PCoA based on the Bray–Curtis dissimilarity index. In contrast to microbiota data, there were no statistically significant differences in β-diversity between these two groups. DESeq analysis showed that there were only 20 genes differentially regulated in stage IIIB–IV compared with stage I–IIIA NSCLC (Supplementary Fig. S15; Supplementary Table S8). Similarly, very few transcripts were found differentially expressed when comparing better versus worse outcomes at six-month and one-year survival (Supplementary Table S8).

We then used DESeq to compare transcriptomic signatures associated with a distinct lower airway microbiota based on DMM and found that there were 209 genes upregulated and 88 genes downregulated in airway brushes of subjects with SPT lower airway microbiota versus BPT lower airway microbiota (Fig. 2A; Supplementary Table S9, FDR < 0.25). Subanalysis of the transcriptomic data among stage I–IIIA and IIIB–IV NSCLC showed the most significant differences for SPT versus BPT within stage I–IIIA NSCLC. Functional enrichment analysis [Ingenuity Pathway Analysis (IPA)] of differentially expressed genes between SPT and BPT (all samples or stage I–IIIA NSCLC samples) showed that SPT was associated with upregulation of the following top canonical pathways: p53 mutation, PI3K/PTEN, ERK, and IL6/IL8 (Fig. 2B).

Figure 2.

Airway transcriptome in NSCLC based on lung microbiota. Comparisons between microbiome and host transcriptomic signatures were conducted using samples where paired matched data were available (n = 70). A, Volcano plot of differentially expressed genes (FDR < 0.25) between pneumotypeSPT vs. pneumotypeBPT in all, stage I–IIIA only, or stage IIIB–IV only lower airway samples. B, Unsupervised hierarchical heat map of canonical pathway analysis based on IPA (RRID:SCR_008653) using the airway transcriptome of all subjects and those with stage I–IIIA disease comparing pneumotypeSPT vs. pneumotypeBPT groups. Subanalysis using samples from patients with stage IIIB–IV disease is not presented given the paucity of differentially expressed genes between the groups. Orange shows upregulation of pathway, and blue shows downregulation of pathway. C, Network analysis based on conditional co-occurrence probability of microbiome and transcriptome data; microbiome nodes (circles) are colored red for stage IIIB–IV lung cancer and green for stage I–IIIA lung cancer (based on a gradient) and sized by relative abundance. Edges connect microbiome nodes to pathway nodes, and edge width is based on their conditional probability.

Figure 2.

Airway transcriptome in NSCLC based on lung microbiota. Comparisons between microbiome and host transcriptomic signatures were conducted using samples where paired matched data were available (n = 70). A, Volcano plot of differentially expressed genes (FDR < 0.25) between pneumotypeSPT vs. pneumotypeBPT in all, stage I–IIIA only, or stage IIIB–IV only lower airway samples. B, Unsupervised hierarchical heat map of canonical pathway analysis based on IPA (RRID:SCR_008653) using the airway transcriptome of all subjects and those with stage I–IIIA disease comparing pneumotypeSPT vs. pneumotypeBPT groups. Subanalysis using samples from patients with stage IIIB–IV disease is not presented given the paucity of differentially expressed genes between the groups. Orange shows upregulation of pathway, and blue shows downregulation of pathway. C, Network analysis based on conditional co-occurrence probability of microbiome and transcriptome data; microbiome nodes (circles) are colored red for stage IIIB–IV lung cancer and green for stage I–IIIA lung cancer (based on a gradient) and sized by relative abundance. Edges connect microbiome nodes to pathway nodes, and edge width is based on their conditional probability.

Close modal

Multiomic Analysis

To better characterize host/microbe interaction in lung cancer, we used a multiomic analytic framework that evaluates for associations between co-occurring taxa and host RNA transcriptome signatures. We estimated co-occurrence probabilities between taxa and the host transcripts adding the probability ranks for the taxa being associated with stage I–IIIA or IIIB–IV lung cancer using MMvec (27, 28). Based on the predicted microbe–transcript co-occurrences, there were two distinct clusters of taxa (Fig. 2C; interactive figure available at https://segalmicrobiomelab.github.io). The first cluster consisted of SPT-associated taxa (belonging to the genera Veillonella, Prevotella, and Streptococcus) that had a high probability of being observed in subjects with stages IIIB–IV. The second cluster consisted of BPT-associated taxa (such as Flavobacterium genus) that had a high probability of being observed in subjects with stage I–IIIA NSCLC; however, it is important to note that many of the high abundant genera in this cluster (stages I–IIIA) likely represent background taxa as identified by decontam (Supplementary Fig. S13) and not true lower airway taxa. Among SPT-associated taxa, a Veillonella taxon (OTU#585419) had the highest relative abundance and a high probability of being found in subjects with stage IIIB–IV lung cancer. This taxon was also highly associated with cell adhesion molecules, IL17, cytokines and growth factors, chemokine signaling pathway, TNF, JAK–STAT, and PI3K–AKT signaling pathway (Supplementary Table S10). Using BLAST (28), the sequence of this OTU most closely aligned with Veillonella parvula.

Lung Dysbiosis Preclinical Model

To evaluate the causal effects of lower airway dysbiosis on lung cancer progression, we tested the effects of lower airway dysbiosis induced by Veillonella parvula in a preclinical lung cancer model (KP mice; Fig. 3A). We selected this bacterium because we have found it to be a good marker for SPT, it was consistently associated with NSCLC (16), and it was the taxa with the highest relative abundance identified in our multi-omic analysis as associated with stage IIIB–IV and transcriptomic signatures. Of note, lower airway dysbiosis induced by other oral commensals, such as Streptococcus mitis and Prevotella melaninogenica, also led to increased lower airway inflammation but at a lesser degree than V. parvula (Supplementary Figs. S16 and S17A–S17B). Thus, as a proof of concept, we chose Veillonella parvula as our lower airway dysbiosis model for the KP lung cancer mice.

Figure 3.

Preclinical model of lung dysbiosis in lung cancer and cancer survival. A, Experimental condition and Kaplan–Meier survival showing decreased survival in mice with lung cancer and dysbiosis (LC + Dys, n = 22) compared with LC (n = 20) alone (P < 0.001). Dys did not affect mouse survival in WT control (n = 10 for each group). B, Quantitative data of tumor burden (measured as lumens prior to death or sacrifice normalized to baseline lumens) showing that LC + Dys mice had increased tumor burden (P < 0.05, n = 5 for each experimental condition). C, IPA was used to identify dysregulated transcriptomic pathways. D, Immune profiling of lung tissue by FACS and cytokine measurement demonstrates that lower airway dysbiosis induces Th17 and PD-1 T-cell phenotype in the lung. E, IHC analysis comparing LC and LC + dys shows increase in CD4+ T cells, CD8+ T cells, and neutrophils in the non-tumor region after dysbiosis. Minimal difference in immune response was seen within the tumor itself [n = 4 (LC) vs. n = 8 (LC + dys) mice/group; each dot represents different regions analyzed color coded by mice].

Figure 3.

Preclinical model of lung dysbiosis in lung cancer and cancer survival. A, Experimental condition and Kaplan–Meier survival showing decreased survival in mice with lung cancer and dysbiosis (LC + Dys, n = 22) compared with LC (n = 20) alone (P < 0.001). Dys did not affect mouse survival in WT control (n = 10 for each group). B, Quantitative data of tumor burden (measured as lumens prior to death or sacrifice normalized to baseline lumens) showing that LC + Dys mice had increased tumor burden (P < 0.05, n = 5 for each experimental condition). C, IPA was used to identify dysregulated transcriptomic pathways. D, Immune profiling of lung tissue by FACS and cytokine measurement demonstrates that lower airway dysbiosis induces Th17 and PD-1 T-cell phenotype in the lung. E, IHC analysis comparing LC and LC + dys shows increase in CD4+ T cells, CD8+ T cells, and neutrophils in the non-tumor region after dysbiosis. Minimal difference in immune response was seen within the tumor itself [n = 4 (LC) vs. n = 8 (LC + dys) mice/group; each dot represents different regions analyzed color coded by mice].

Close modal

Dysbiosis was induced once KP seeding was determined. Induction of lower airway dysbiosis with V. parvula in wild-type (WT) mice did not affect the mice's survival or weight gain. In contrast, within KP lung cancer mice, exposure to dysbiosis (KP + Dys) led to decreased survival, weight loss, and increased tumor burden (Fig. 3A and B; Supplementary Fig. S18A and S18B). The experiment was repeated at an early sac time point (three weeks post induction of dysbiosis) to evaluate the immune response to dysbiosis with host transcriptomics, T-cell profiling, and cytokine measurements. PCoA of host transcriptomics showed clear differences between the four experimental conditions, where dysbiosis led to greater compositional changes than lung cancer alone (Supplementary Fig. S19A). Characterization of immune cell subsets inferred from bulk transcriptomics (CIBERSORT) identified clear clustering by condition where lower airway dysbiosis led to an increase in Th1 cells and activation of dendritic cells (Supplementary Fig. S19B). IPA showed that dysbiosis led to upregulation of the PI3K/AKT, ERK/MAPK, IL17A, IL6/IL8, and inflammasome pathways (Fig. 3C). Comparisons between transcriptomic signatures induced by lower airway dysbiosis in the NSCLC mouse model and those identified in SPT among subjects with NSCLC showed concordant signatures related to IL17 signaling, chemokine, Toll-like receptor, PD-L1 signaling, and PI3K–AKT signaling, among others (Supplementary Fig. S20A and S20B). Although there are notable differences between transcriptomic signatures in human and mouse data, these findings provide a promising direction for follow-up. Lastly, lung dysbiosis induced by V. parvula led to the recruitment of Th17 cells, with increased levels of IL17 production, increased expression of PD-1+ T cells, and recruitment of neutrophils (Fig. 3D; Supplementary Fig. S21). Spatial analysis with IHC targeting CD4+ T cells, CD8+ T cells, and neutrophils shows that the increase of these inflammatory cells in response to dysbiosis occurred predominantly in tumor-spared lung tissue (Fig. 3E; Supplementary Fig. S22A). Interestingly, in the tumor there was a decrease in CD4+ T cells associated with lower airway dysbiosis.

To further assess the functional importance of dysbiotic-induced IL17 activation in lung tumorigenesis, dysbiotic-KP mice were treated with monoclonal antibodies against IL17 or isotype antibody control for two weeks after tumor initiation (Fig. 4A). Tumor luminescence data showed that IL17 blockade led to a decrease in tumor burden over the second week compared with isotype control (P = 0.0059; Fig. 4B). Immune profiling evaluated at day 14 after IL17 blockade showed that treatment with anti-IL17 antibodies was associated with decreased RORγt+ CD4+ T cells and neutrophils, and a nonstatistically significant trend toward lower IL17+ CD4+ and IL17+TCRγδ+ T cells (Fig. 4C). Histologic assessment with IHC shows that IL17 blockade led to a decrease in CD4+ T cells, CD8+ T cells, and neutrophils in the spared nontumor lung tissue but not in the tumor itself (Fig. 4D; Supplementary Fig. S22B). Overall, these data suggest that lower airway dysbiosis contributes to a tumor-inflammatory microenvironment characterized by an increase in the Th1 and Th17 phenotypes, activation of dendritic cells with potential antigen presentation capacity, and an increase in checkpoint inhibitor markers within the surrounding lung tissue.

Figure 4.

IL17 blockade during lung dysbiosis in the lung cancer preclinical model. A, Experimental conditions (anti-IL17 or isotype Ab control) were administered to KP mice with lower airway dysbiosis induced by Veillonella parvula. B, Anti-IL17 therapy was associated with decreased tumor burden change during the two weeks of antibody injections as compared with isotype control. C, Immune profiling of lung tissue by FACS demonstrates that IL17 blockade of KP mice with lower airway dysbiosis decreases RORyt+ and neutrophils (n = 6–7 for each experimental condition). D, IHC analysis shows that IL17 blockade of KP mice with lower airway dysbiosis decreases CD4+/CD8+ T cells and neutrophils in the non-tumor region (P < 0.0001 and P = 0.0002, respectively). However, a minimal difference in immune response was seen within the tumor itself [n = 8 (LC + dys) vs. n = 6 (LC + dys + anti-IL17) mice/group; each dot represents different regions analyzed color coded by mice].

Figure 4.

IL17 blockade during lung dysbiosis in the lung cancer preclinical model. A, Experimental conditions (anti-IL17 or isotype Ab control) were administered to KP mice with lower airway dysbiosis induced by Veillonella parvula. B, Anti-IL17 therapy was associated with decreased tumor burden change during the two weeks of antibody injections as compared with isotype control. C, Immune profiling of lung tissue by FACS demonstrates that IL17 blockade of KP mice with lower airway dysbiosis decreases RORyt+ and neutrophils (n = 6–7 for each experimental condition). D, IHC analysis shows that IL17 blockade of KP mice with lower airway dysbiosis decreases CD4+/CD8+ T cells and neutrophils in the non-tumor region (P < 0.0001 and P = 0.0002, respectively). However, a minimal difference in immune response was seen within the tumor itself [n = 8 (LC + dys) vs. n = 6 (LC + dys + anti-IL17) mice/group; each dot represents different regions analyzed color coded by mice].

Close modal

The lower airway microbiota, whether in health or disease state, are mostly affected by aspiration of oral secretions, and the lower airway microbial products are in constant interaction with the host immune system (15, 19, 29–31). In this study, we demonstrate that a lower airway dysbiotic signature present in patients with lung cancer affects tumor progression and clinical prognosis, likely due to alteration in stage I–IIIA immune tone promoting inflammation and checkpoint inhibition. First, patients with stage IIIB–IV NSCLC are more likely to have enrichment of the lower airway microbiota with oral commensals compared with patients with stage I–IIIA disease. In addition, this dysbiotic signature was associated with: (i) worse outcome at six months and one year (for both stage I–IIIA and IIIB–IV groups); (ii) OS in the stage I–IIIA group; and (iii) tumor progression in stage IIIB–IV disease. Our preclinical data using an NSCLC mouse support a model in which aspiration of oral commensals (identified in our human cohort) affects the lower airway inflammatory tone and promotes tumor cell proliferation. Dysbiosis in these mice led to upregulation of the ERK/MAPK, IL1, IL6, and inflammasome signaling pathways. Immune profiling showed that lung dysbiosis led to a substantial increase in Th17 and PD-1+ cells. Previous preclinical models of cancer have shown the association between lung dysbiosis and lung inflammation but have limited human microbiome data to support clinical relevance (especially considering that the human and murine microbiota differ; refs. 32–35). Our data identified that enrichment of the lower airway microbiota with human oral commensals, such as Veillonella, contributes to a local protumor immune tone leading to progression of NSCLC, suggesting that microaspiration and/or impaired airway clearance likely affect the pathogenesis of this disease (36).

Several lines of investigation have shown that increased inflammation and decreased immune surveillance, characterized by IL17 tone and checkpoint inhibition, are associated with poor prognosis in NSCLC. Increased local and systemic IL17 (37, 38), systemic IL6 (39), and higher neutrophil-to-T cell ratio (40) are associated with a poor prognosis in lung cancer. PD-L1, the ligand for PD-1, is induced in nonlymphoid cells and tumor cells under inflammatory conditions triggered by several cytokines, such as IFNγ and pathogen-associated molecular patterns (PAMP; refs. 41–43). In addition, many signaling molecules (e.g., NFκB, MAPK, PI3K, mTOR, and JAK/STAT) that affect proliferation, apoptosis, and cell survival induce PD-L1 expression (44, 45). In a bitransgenic mouse model expressing a conditional IL17A allele and a conditional KRASG12D, increased IL17 production was associated with accelerated lung tumor growth, decreased responsiveness to checkpoint inhibition, and decreased survival (46). In many cancer models (breast cancer, gastric carcinoma, and lung cancer), inflammasome activation, through IL1β signaling, leads to an inflammatory response characterized by decreased antitumor immune surveillance (47–49). In the current investigation, we show that the increase in IL17 inflammatory tone triggered by lower airway dysbiosis can be blunted by anti-IL17 blocking antibodies, which seemed to lead to a decrease in the tumor burden. More experiments are obviously needed to further characterize the phenotypic inflammatory profile in the tumor and surrounding tissue, to understand the molecular mechanisms by which lower airway inflammatory cells respond to lower airway dysbiosis, and to better characterize how these factors affect tumor burden and survival. However, the above-discussed investigations and the data presented in the current paper support that the balance between Th17 inflammation and immune surveillance affects NSCLC pathogenesis, and, thus, future investigations are warranted to explore the role of IL17 blockade in this disease.

Immune-checkpoint molecules, such as PD-1, mediate the response of T cells to neoantigens and are now first-line therapy for advanced NSCLC (4–8). However, 40% to 60% of patients will not benefit from these therapies, and existing biomarkers (e.g., expression of PD-1 ligand) have limited capacity to predict efficacy (7, 50). Different gut microbiota signatures have been identified as associated with augmenting antitumor immunity and a PD-1 blockade response (9–11). In the gut, higher α-diversity and enrichment of Ruminococcaceae were associated with a favorable response to anti–PD-1 treatment in patients with melanoma (10, 51), and modulation of the microbiota in germ-free mice can enhance antitumor immunity and augment effects of checkpoint blockade (12, 13). In germ-free or antibiotic-treated mice, lung adenocarcinoma (KRAS mutation/p53 deletion) development is decreased compared with specific pathogen-free mice (32). In this model, lung microbiota activate IL1β and IL23 cytokines from myeloid cells and induce IL17-producing γδ T cells. Thus, although most studies have focused on the effect of the gut microbiome on cancer development and progression, there is increasing evidence to suggest that the local lung microbiota play a pivotal role in lung cancer pathogenesis and treatment. Multiple lines of investigation have shown that the lower airway microbiota are major determinants of the airway immune tone in health and many disease states. For example, recent preclinical models have shown that lower airway mucosal inflammation is primarily associated with the composition of the lower airway microbiota rather than the composition of the gut or upper airway microbiota (14). In humans, we have shown that pneumotypeSPT is associated with increased local inflammatory cells and the Th17 phenotype (15, 52), and the lower airway microbial metabolism can be modulated by, for example, chronic macrolide therapy, leading to release of microbial metabolites with anti-inflammatory effects (53, 54). Anaerobes are commonly found in the lower airways and can survive oxygen stress by forming multicellular complexes within the hypoxic environment present in biofilms (55, 56). Short-chain fatty acids produced by fermentation, such as butyrate, are one energy source for anaerobes (57), and we have shown that their presence in the lower airways is higher in pneumotypeSPT and regulates IFNγ and IL17A production in CD4+/CD8+ lymphocytes (58). In NSCLC, we recently demonstrated that pneumotypeSPT is associated with several inflammatory cancer-related pathways, such as ERK/MAPK and PI3K/AKT (16), that can lead to chronic inflammation, altered Treg/Th17 balance (59–61), augmented Th17 differentiation (62, 63), and induction of PD-L1 expression (44, 45). Our current findings expand the above observations by demonstrating that a dysbiotic signature characterized by enrichment of the lower airway microbiota with oral commensals can contribute to the progression of disease.

Among the limitations pertinent to this study, we should point out that there is a significant degree of disease heterogeneity and the appropriate subanalyses could be explored only with a much larger cohort. For example, we decided to focus on NSCLC because there were few cases of small cell lung cancer. Further, within NSCLC there were several pathologic subtypes, driver mutation status, PD-L1 status, etc. The small subsample size prevents us from conducting the appropriate subanalysis. However, our analysis and models were stratified and adjusted by staging (dichotomized as I–IIIA and IIIB–IV stage groups and adjusted by individual TNM stages), which is a very significant covariate associated with prognosis and treatment modality. Interestingly, we found a few host transcriptomic signatures associated with a disease stage whereas there were many more transcriptomic signatures associated with lower airway microbiota subtype (SPT/BPT). It is possible that the histologic heterogeneity within NSCLC will affect these results, and a larger cohort may allow to control for this. Other potential confounders related to patients' clinical condition, such as swallowing and deglutition problems, cannot be fully accounted for in the current cohort but may have significant impact on our results. Given our finding of the enrichment of the lower airway microbiota with oral commensal as associated with prognosis, future investigations that include evaluation of swallowing functions are warranted. Low biomass samples are subjected to contamination with background DNA (coming from the reagents, bronchoscopy, or sequencing noise; refs. 64, 65). To address concerns regarding DNA contamination during sample collection and preparation, we applied decontam (26) analysis and showed that Flavobacterium, a taxon identified in the multiomic analysis and dominant of BPT, is likely a background contaminant. This is consistent with prior data showing no growth from lower airway samples characterized as BPT (16). We therefore induced airway dysbiosis in our mouse model with Veillonella and compared it with PBS (which 16S rRNA gene sequencing composition most resembles BPT) rather than a separate bacterium as control. Our investigation supports the hypothesis that the lower airway microbiota contribute to a local protumor immunity; however, we did not investigate the systemic inflammatory response in this model. Further support for the relevance of this mechanism will need to focus on blocking the immune response to the microbial exposure in the setting of lung cancer and evaluating the effects of induced lower airway dysbiosis during immunotherapy. In the current investigation, we did not explore the association between lung microbiota and response to immunotherapy because this treatment was applied in a relatively small fraction of patients (16%) and the vast majority of the samples were collected before this therapy became standard of care. Also, although we identified a taxonomic signature associated with inflammatory tone and prognosis in lung cancer, we cannot determine the molecular signatures present in the microbial community that may be responsible for this association. Future investigations that exploit novel functional microbiomic approaches (e.g., metagenome, metatranscriptome, and metabolome) should focus on molecular markers with significant immunomodulatory activity. In our preclinical model, we tested whether Veillonella parvula was sufficient to induce lower airway inflammation and worsening of tumor progression. Other oral commensal, when present in the lower airways, may also be contributing to this process and may need to be further evaluated as key components of lower airway dysbiosis in isolation or in complex microbial communities. Although the lower airway microbiota were associated with staging and survival, other dysbiotic signatures in other mucosae could also have significant associations. Even though we did not identify significant microbiota signatures in the buccal samples, future investigations should include gut samples as well to establish the relative role of the microbiota of different mucosae niches to the pathogenesis of lung cancer. Finally, further validation of the results presented here will require a second cohort where sampling approach and design are customized to overcome some of the limitations here described.

This study has broad clinical implications regarding lung cancer pathogenesis and treatment response. Identification of lower airway dysbiotic signatures associated with lung cancer prognosis may be important to personalize approaches for lung cancer treatment and prognosis. Fecal microbiota transplant, a strategy with proven efficacy in difficult-to-treat Clostridium difficile infection and inflammatory bowel disease (66, 67), can influence the susceptibility to anti–PD-1 cancer immunotherapy (9, 10), and its clinical impact is now being tested in humans within ongoing clinical trials. Despite the evidence that the local microbiota affect the local inflammatory tone of the lung, there are no human trials aiming to modify the lung microbiome in the setting of malignancy. The data presented here suggest that lower airway dysbiosis induced by microaspiration of oral commensals affects lung tumorigenesis by promoting an IL17-driven inflammatory phenotype, a pathway amenable for targeted therapy that may have potential implications in this disease. A better understanding of the microbial host interaction in the lower airways will be needed to uncover how the lung cancer–associated microbiota could be modulated to affect prognosis and response to immunotherapies.

Subjects

All subjects signed written informed consent to participate in this study that was approved by the Institutional Review Board of New York University. Participants included patients who had suspicious nodules on chest imaging and who underwent clinical bronchoscopy. Lung cancer subtype, somatic mutation, and stage were recorded after histopathologic confirmation. We excluded subjects with a prior history of cancer or recent (less than one month) antibiotic use. RECIST (27) score was analyzed at the 6- to 12-month time point after diagnosis of lung cancer, where these data were most consistently available.

Bronchoscopic Procedure

Both background and supraglottic (buccal) samples were obtained prior to the procedure as previously described (16). The background samples were obtained by passing sterile saline through the suctioning channel of the bronchoscope prior to the procedure. For this project, we obtained multiple lower airway samples from different locations, including 82 from the right mainstem, 59 from the airways leading to the lung cancer lesion (involved segments), and 69 from the airways spared of disease on the contralateral lung. A detailed description of the number of samples and the analyses performed in them is provided in Supplementary Table S6.

Bacterial 16S rRNA-Encoding Gene Sequencing

High-throughput sequencing of bacterial 16S rRNA-encoding gene amplicons (V4 region; ref. 68) was performed. Reagent control samples and mock mixed microbial DNA were sequenced and analyzed in parallel (Supplementary Fig. S23). The obtained 16S rRNA gene sequences were analyzed with the Quantitative Insights Into Microbial Ecology (QIIME; RRID:SCR_008249) 1.9.1 package (69). OTUs were not removed from upstream analysis. PERMANOVA testing was used to compare the compositional differences of groups. A prevalence-based method using the R package decontam (v1.6.0; ref. 26) was used to identify potential contaminants from the sequencing data sets. In this process, all reads from background bronchoscope control samples were identified as negative controls and thus possible sources of contaminants. No OTU was removed from the analyses performed and data from the 16S microbiome for this manuscript are available (Sequence Read Archive, RRID:SCR_001370: #PRJNA592147).

Sample clustering of meta-communities was based on Dirichlet multinomial mixtures (DMM) modeling (70).

Transcriptome of Bronchial Epithelial Cells

RNA-seq was performed on bronchial epithelial cells obtained by airway brushing, as described (71–73), using the Hi-Seq/Illumina platform at the NYU Langone Genomic Technology Center (data available at Sequence Read Archive: # PRJNA600487). KEGG (74, 75) annotation was summarized at levels 1 to 3. Genes with a false discovery rate (FDR)–corrected adjusted P value <0.25 were considered significantly differentiated, unless otherwise specified. Pathway analysis using differentially regulated genes (FDR < 0.25) was done using IPA (RRID:SCR_008653; QIAGEN Inc.; ref. 76). Gene set enrichment analysis was performed with differential genes (FDR < 0.25) for data set comparison, R package fgsea v1.4.1 (77).

Experimental Mouse Model

The mice utilized in this experiment were five-week-old females at the time of use. The strain was B6(Cg)-Tyrc-2J/J mice purchased via vendor (Jackson Laboratory; cat. #000058). The mice were kept in Skirball Animal Facility and were kept under controlled conditions with cycles of 12-hour daylight and 12-hour darkness. Mice were euthanized by carbon dioxide asphyxiation followed by cardiac puncture. Blood, skin swabs, oral swabs, lung lavage, lung tissue, humerus bone marrow, cecum, terminal ileum, and fecal pellets were collected for study. The Institutional Animal Care and Use Committee of the New York University School of Medicine approved all procedures, and experiments were carried out following their guidelines (IACUC# s16-00032).

KP Model Lung Cancer

The KP model of lung cancer histopathologically resembles that of human cancers and has been used to study translational models of lung cancer in mice (78). The KP model of lung cancer is based on KRASLSL-G12D/+;p53fl/fl. NSCLC models require induction by use of replication-deficient adenoviruses expression Cre (Ad-Cre) to induce transient Cre expression in the lungs of mice. Once tumor burden is increased in the mice, the lungs are harvested and the KP lung cancer cells grown in cell culture (79). Cell culture lines of KP lung cancer cells are grown in DMEM 10% FBS plus gentamicin under aerobic conditions with 5% carbon dioxide at 37°C. Cells were harvested from the cell culture when 90% congruent. The goal was to grow cells to 3,000,000 KP cells/mL (or 150,000 cells/50 μL). To detect in vivo luminescence, images were acquired using the IVIS spectrum (PerkinElmer) after intraperitoneal injection of Luciferin (Promega). We then proceeded to intratracheal inoculation of KP cells. The mice were anesthetized utilizing isoflurane until sedated. The mice were then placed on an intubation platform, and with blunt forceps, their tongue was gently pulled ventrally until the pharynx was exposed (78). Then, an Exel Safelet catheter (Exel International Inc.; cat. # 26746) was introduced through the vocal cords of the mice, and a 50-μL inoculum of lung cancer (1.5 × 105 KP cells) was placed into the catheter. The mice were then removed from the intubation platform to recover from anesthesia on a heat pad.

Creation of Veillonella Parvula Inoculum

The following human oral commensals were obtained: Veillonella parvula, Prevotella melaninogenica, and Streptococcus mitis (ATCC). These bacteria were grown in anaerobic conditions (Bactron 300, Shel Labs) and then stored in 20% glycerol tryptic soy broth at −80°C. To prepare oral commensal challenges, bacteria strains were thawed and streaked on anaerobic PRAS-Brucella Blood agar plates (Anaerobe Systems). The plates were incubated at 37°C in an oxygen-free environment (tri-mix: 5% carbon dioxide, 5% hydrogen, and 90% nitrogen) in an anaerobic chamber for 24 to 48 hours. The colonies were collected from the plate and resuspended in 1 mL of sterile PBS. The OD620 was measured to calculate the approximate concentration prior to use.

Intratracheal Microbial and Control Challenge

Mice were assigned to receive the microbial challenge with Veillonella parvula twice a week via intratracheal inoculation starting two weeks after the inoculation with lung cancer. First, mice were sedated with the use of isoflurane anesthesia. The mice were then suspended by their dorsal incisors upon an elastic cord; a blunt pair of forceps was used to ventrally pull the tongue forward to expose the larynx. Then, a pneumatic otoscope (Welch-Allyn; cat. #71000C) with a 2-mm ear specula was advanced until the vocal cords were visualized. Using a gel-loading tip, a 50-μL volume of the Veillonella parvula was deployed into the trachea of the mouse. These exposures occurred twice a week, spaced three to four days/week apart. Mice were monitored during this process; no mice died due to the inoculation procedure. A control procedure to inoculate mice with PBS was performed in the same manner.

Immune Inhibition Experiment

Two weeks after KP cell inoculation, mice were challenged intratracheally with Veillonella parvula similar to above. At this time, mice were randomized 1:1 to receive anti-IL17 (1 mg/mL; Bio X Cell), anti-IL17 isotype control (2 mg/mL; Bio X Cell). Antibody dose was diluted in 100 μL and given via intraperitoneal injection twice a week for a total of two weeks.

Organization and Measurements on Mice

Once lung tumor development was detected by IVIS (two weeks post inoculation), mice were randomized according to tumor burden to receive either PBS or dysbiosis with V. parvula while maintaining co-house conditions. For the KP mice, those with median lumens of 8 × 105 to 7 × 106 photon-flux (photons/s/cm2/steradian) at two weeks were utilized for the experiments. Wild-type mice from the same strain and no KP exposure were used as control mice and were exposed to sterile PBS or V. parvula. Thus, in all experiments, mice were organized to the following groups: (i) wild-type with PBS control (WT), (ii) wild-type with dysbiosis with V. parvula (Dys), (iii) KP lung cancer with PBS control (LC), and (iv) KP lung cancer with V. parvula (LC + Dys). Imaging the mice utilizing luciferins expression (lumens) occurred two weeks after inoculation with KP lung cancer cells. The platform we used to image the mice was PerkinElmer IVIS Spectrum (PerkinElmer, cat. # 124262). Luciferin (1.5 mg; PerkinElmer, Xeno-Light D-Luciferin Potassium Salt; cat. # 122799) was given intraperitoneally. Mice received 50 μL of their respective inoculum with the Veillonella condition receiving 1.5 × 106 cfu/mL. The mice were organized into groups based upon their median lumens to establish experimental groups of mice with the same luminosity for a baseline. The imaging of the mice occurred twice every week on the day prior to inoculation. For the survival experiment, we utilized 60 mice that were followed for six weeks after initiation of microbial challenge or PBS control. Forty additional mice were divided in the same four conditional groups for immune phenotyping on lung homogeneate, including lung transcriptomics, flow cytometry, and cytokine measurement. For this experiment, mice were sacrificed after two weeks post initiation of microbial or PBS exposure. For host RNA transcriptome, flash-frozen lung samples were defrosted and then homogenized utilizing a hand TissueRuptor II on the second lowest setting (Qiagen). Then samples were spun down on a tabletop centrifuge 14,000 rpm for two minutes, and the pellet was collected and sent for RNA processing. RNA was extracted from collected supernatant using the Qiagen miRNeasy Mini Kit (Qiagen; cat. #74135). Quality control was established with RNA integrity number cutoff >6. RNA-seq was performed using Hi-Seq (Illumina) at the NYU Genomic Technology Center. RNA-seq library preps were made using the Illumina TruSeq Stranded mRNA LT Kit (Illumina; cat. #RS-1222-2101) on a Beckman Biomek FX instrument, using 250 ng of total RNA as input, amplified by 12 cycles of PCR, and run on an Illumina 2500 (v4 chemistry), as single-read 50 bp. Sequences from the murine lung homogenate were aligned against the murine ensemble reference genome utilizing STAR, RRID:SCR_015899 (v2.5) aligner (80). Gene counting of each sample was performed using featureCounts, RRID:SCR_012919 (v1.5.3; refs. 81, 82). FACS was performed on single-cell suspension derived from lung homogenate. First, lung samples were minced and dissociated utilizing Liberase (Hoffmann-La Roche) for 35 minutes in a 37°C water bath and followed by mechanical disruption through a 70-μm filter. Liberase was used at a concentration of 0.5 mg/mL in DMEM supplemented with 10% FBS. For intracellular cytokine staining, the cells were treated with a cell stimulation and protein transport inhibition cocktail containing PMA, Ionomycin, Brefeldin A, and Monensin (500x eBioscience Affymetrix) for four hours. The cells were surface stained, fixed in 2% paraformaldehyde, and permeabilized with 0.5% saponin. Cell staining with fluorochrome-conjugated antibodies was performed targeting CD3+, CD4+, CD8+, CD69+, PD1+, and IL17+ (Thermo-Fisher), and measurements were performed on a BD LSR II flow cytometer (BD Biosciences). Acquired data were analyzed using FlowJo, RRID:SCR_008520 version 10.3 (Tree Star Inc.). Cytokines and chemokines were measured using Luminex (Murine Cytokine Panel II, EMD Millipore). Lung homogenates were thawed and processed according to the recommended protocol using the Murine Cytokine/Chemokine Magnetic Bead Panel (# MCYTMAG-70K-PXkl32). All cytokines/chemokines concentrations were normalized by the gram of lung homogenate and included those with dynamic range: G-CSF, Eotaxin, IFNg, IL1a, Il1b, IL3, IL4, IL5, IL6, IL7, IL9 IL10, IL12p40, IL12p70, LIF, IL17, IP10, KC, MCP1, MIP1a, MIP1b, M-CSF, MIP2, MIG, RANTES, VEGF, and TNFa.

Multiplex Immunostaining

Five-micrometer sections of paraffin-embedded preserved lung tissue were stained with Akoya Biosciences Opal multiplex automation kit reagents unless stated otherwise. Automated staining was performed on Leica BondRX autostainer. The protocol was performed according to the manufacturers' instructions with the antibodies specified in Supplementary Table S7. Briefly, all slides underwent sequential epitope retrieval with Leica Biosystems epitope retrieval 1 (ER1; citrate based, pH 6.0; cat. AR9961) and 2 solution (ER2; EDTA based, pH9; cat. AR9640), primary and secondary antibody incubation and tyramide signal amplification with Opal fluorophores (Supplementary Table S7). Primary and secondary antibodies were removed during epitope retrieval steps while fluorophores remain covalently attached to the epitope.

Image Acquisition and Analysis

Semiautomated image acquisition was performed on a Vectra Polaris multispectral imaging system. After whole-slide scanning at 20×, the tissue was manually outlined to select fields for spectral unmixing and analysis using InForm version 2.4.10 software from Akoya Biosciences. Fields of view for analysis were separated as containing tumor only or areas of pulmonary parenchyma where tumor was not apparent. For each field of view, cells were segmented based on nuclear signal (DAPI). Cells were phenotyped after segmentation using InForm's trainable algorithm based on glmnet (83) package in R. Four algorithms were created to classify cells as Ly6g+ (neutrophils) or “other,” CD4+ or “other,” CD8+ or “other” and F4/80+ or “other.” Phenotypes were reviewed for different samples during training iterations. Data were exported as text containing sample names, field of acquisition coordinates, individual cell information including coordinates and identified phenotype. Each image was analyzed with all four algorithms, so that every cell was classified four times. Concatenation of all phenotyping information was performed in R using the Phenoptr Reports package (Kent S. Johnson 2020). phenoptr: InForm Helper Functions. R package version 0.2.7. https://akoyabio.github.io/phenoptr/) in RStudio software [RStudio Team (2015). RStudio: Integrated Development for R. RStudio, Inc.; URL http://www.rstudio.com/.] Statistical analysis (Mann–Whitney U test) was run for the following groups: lung cancer vs. lung cancer + dysbiosis (n = 4 and 8 mice, respectively, Fig. 3E), and lung cancer + dysbiosis vs. lung cancer + dysbiosis + anti-IL17 (n = 8 and 6 mice, respectively; Fig. 4D), taking each field as an independent value.

Statistical and Multiomic Analyses

In Supplementary Table S2, categorical variables were presented as frequencies and percentages and their distribution difference between groups with dead or alive OS status were assessed by the Fisher exact test. The Cox proportion hazards models (84) were used to evaluate each variable's marginal association with the time to death. Hazard ratio and P value were reported.

The MiRKAT (85) was used to investigate whether the community-level microbial profile among lower airway samples was different between any paired samples from right main, involved, or noninvolved locations, and between stage I–IIIA and IIIB–IV while adjusting for smoking status within each location samples. The survival version of MiRKAT test: MiRKAT-S (86) was used to investigate whether the community-level microbial profile is associated with the OS while adjusting for smoking status, stage, and surgery within each location sample. The paired Bray–Curtis dissimilarity was used in all tests.

For the taxonomic level analysis, we used the linear mixed-effect model on the arcsine square root–transformed relative abundance at genus level for their associations with stage (I–IIIA/IIIB–IV; Supplementary Table S3). In the model, the subject was set as the random effect to take care of the correlation among three location samples from the same subjects. The stage was set as fixed effect while adjusting for smoking status. We used the two-stage linear mixed-effect model (87) on the arcsine square root–transformed relative abundance at genus level for their associations of the OS (Supplementary Table S4) while adjusting for smoking status, stage, and surgery. In the first stage, the linear mixed-effect model was used to take care of the correlation among three location samples from the same subjects. The random intercept estimates from the first stage were used in the Cox proportional hazards model in the second stage to investigate their association with the OS.

Because the distributions of microbiome data are nonnormal, and no distribution-specific tests are available, we used nonparametric tests of association. For association with discrete factors, we used either the Mann–Whitney test (in the case of two categories) or the Kruskal–Wallis ANOVA (in the case of >2 categories). For tests of association with continuous variables, we used the nonparametric Spearman correlation tests. FDR was used to control for multiple testing (88). To evaluate for taxonomic or transcriptomic differences between groups, we utilized DESeq2 (89).

Differential abundance of microbes related to lung cancer stage (IIIB–IV vs. I–IIIA) was calculated using Songbird as previously described (90). We then computed the microbe-transcript co-occurrence probability (probability of observing a transcriptomic pathway when a microbe is observed) using mmvec (91). A probability matrix of the top 10 transcriptome-related pathways for each microbe was generated and used to create a network based on the Fruchterman–Reingold force-directed algorithm using R package ggnet v 0.1.0. (https://cran.r-project.org/web/packages/GGally/index.html). Microbe nodes were colored based on differential analysis of stage IIIB–IV versus I–IIIA NSCLC.

Data Storage

Sequencing data are available at Sequence Read Archive (92, 93) under accession number 16S Microbiome PRJNA592147, Human RNA-seq PRJNA600487, and Murine RNA-seq PRJNA600489. Codes utilized for the analyses presented in the current article are available at https://github.com/segalmicrobiomelab/lung_cancer_prognosis_microbiome.

J. Carpenito reports grants from NIH during the conduct of the study. J. Pass reports grants from NCI EDRN U01 and AACR/J&J Grant during the conduct of the study. K.-k. Wong reports being a founder and equity holder of G1 Therapeutics. K.-k. Wong has sponsored research agreements with Pfizer, Janssen, Takeda, Targimmune, AstraZeneca, Novartis, Merck, Ono, and Array. L. Segal reports grants from NCI/NIH, FNIH, and AACR Grant during the conduct of the study. No disclosures were reported by the other authors.

J.-C.J. Tsay: Conceptualization, resources, data curation, formal analysis, funding acquisition, validation, investigation, methodology, writing–original draft, writing–review and editing. B.G. Wu: Investigation, methodology, writing–original draft. I. Sulaiman: Conceptualization, data curation, software, formal analysis, investigation, visualization, and methodology. K. Gershner: Data curation. R. Schluger: Data curation. Y. Li: Resources, data curation, and supervision. T.-A. Yie: Data curation. P. Meyn: Data curation and software. E. Olsen: Data curation, software, and formal analysis. L. Perez: Data curation. B. Franca: Data curation, software, and formal analysis. J. Carpenito: Data curation and investigation. T. Iizumi: Investigation. M. El-Ashmawy: Data curation, software, and formal analysis. M. Badri: Conceptualization, data curation, software, formal analysis, visualization, methodology, and writing–original draft. J.T. Morton: Conceptualization, data curation, software, investigation, visualization, methodology, and writing–original draft. N. Shen: Data curation, software, and formal analysis. L. He: Data curation, software, and formal analysis. G. Michaud: Data curation. S. Rafeq: Data curation. J.L. Bessich: Data curation. R.L. Smith: Data curation. H. Sauthoff: Data curation. K. Felner: Data curation. R. Pillai: Investigation. A.-M. Zavitsanou: Investigation. S.B. Koralov: Methodology. V. Mezzano: Software, investigation, and visualization. C.A. Loomis: Conceptualization and resources. A.L. Moreira: Software, investigation, and visualization. W. Moore: Formal analysis. A. Tsirigos: Conceptualization, supervision, and methodology. A. Heguy: Conceptualization, supervision, and methodology. W.N. Rom: Conceptualization and supervision. D.H. Sterman: Conceptualization and supervision. H.I. Pass: Conceptualization, data curation, supervision, funding acquisition, and investigation. J.C. Clemente: Conceptualization, software, formal analysis, and supervision. H. Li: Conceptualization, formal analysis, supervision, and validation. R. Bonneau: Conceptualization and supervision. K.-K. Wong: Conceptualization, resources, supervision, funding acquisition, methodology, writing–review and editing. T. Papagiannakopoulos: Conceptualization, supervision, funding acquisition, validation, investigation, methodology, writing–review and editing. L.N. Segal: Conceptualization, resources, data curation, software, formal analysis, supervision, funding acquisition, investigation, visualization, methodology, writing–original draft, project administration, writing–review and editing.

This study was supported by R37 CA244775 (L.N. Segal, NIH/NCI); PACT grant (L.N. Segal, FNIH); K23 AI102970 (L.N. Segal, NIH/NIAD); EDRN 5U01CA086137-13 (W.N. Rom); DoD W81XWH-16-1-0324 (J.-C.J. Tsay); the 2018 AACR–Johnson and Johnson Lung Cancer Innovation Science Grant Number 18-90-52-ZHAN (H.I. Pass/L.N. Segal); A Breath of Hope Foundation (J.-C.J. Tsay); Simons Foundation (R. Bonneau); and CTSI Grant #UL1 TR000038 (L.N. Segal). The Genome Technology Center and ExpPath Core are partially supported by the Cancer Center Support Grant P30CA016087 at the Laura and Isaac Perlmutter Cancer Center (A. Heguy and A. Tsirigos); T32 CA193111 (B.G. Wu); UL1TR001445 (B.G. Wu); FAMRI Young Clinical Scientist Award (B.G. Wu); Stony Wold-Herbert Fund Grant-in-Aid/Fellowship (B.G. Wu, I. Sulaiman, and K. Gershner); R01 HL125816 (L.N. Segal and S.B. Koralov, NIH/NHLBI); and R01 DK110014 (H. Li and L. He). We would like to thank the Genome Technology Center (GTC) for expert library preparation and sequencing, the Applied Bioinformatics Laboratories (ABL) for providing bioinformatics support and helping with the analysis and interpretation of the data, and Experimental Pathology Research Laboratory for histopathology services and imaging. GTC, ExpPath Core, and ABL are shared resources partially supported by the Cancer Center Support Grant P30CA016087 at the Laura and Isaac Perlmutter Cancer Center. This work has used computing resources at the NYU School of Medicine High Performance Computing Facility. Financial support for the PACT project is possible through funding support provided to the FNIH by AbbVie Inc., Amgen Inc., Boehringer-Ingelheim Pharma GmbH and Co. KG, Bristol-Myers Squibb, Celgene Corporation, Genentech Inc., Gilead, GlaxoSmithKline plc, Janssen Pharmaceutical Companies of Johnson and Johnson, Novartis Institutes for Biomedical Research, Pfizer Inc., Sanofi, and Shared Instrumentation Grant S10 OD021747 award for the Vectra Imaging system in Experimental Pathology.

1.
Siegel
RL
,
Miller
KD
,
Jemal
A
. 
Cancer statistics, 2019
.
CA Cancer J Clin
2019
;
69
:
7
34
.
2.
Cancer Genome Atlas Research N
. 
Comprehensive molecular profiling of lung adenocarcinoma
.
Nature
2014
;
511
:
543
50
.
3.
Rosell
R
,
Karachaliou
N
. 
Large-scale screening for somatic mutations in lung cancer
.
Lancet
2016
;
387
:
1354
6
.
4.
Patel
SP
,
Kurzrock
R
. 
PD-L1 expression as a predictive biomarker in cancer immunotherapy
.
Mol Cancer Ther
2015
;
14
:
847
56
.
5.
Rizvi
NA
,
Hellmann
MD
,
Snyder
A
,
Kvistborg
P
,
Makarov
V
,
Havel
JJ
, et al
Cancer immunology. Mutational landscape determines sensitivity to PD-1 blockade in non-small cell lung cancer
.
Science
2015
;
348
:
124
8
.
6.
Herbst
RS
,
Soria
JC
,
Kowanetz
M
,
Fine
GD
,
Hamid
O
,
Gordon
MS
, et al
Predictive correlates of response to the anti-PD-L1 antibody MPDL3280A in cancer patients
.
Nature
2014
;
515
:
563
7
.
7.
Gandhi
L
,
Rodriguez-Abreu
D
,
Gadgeel
S
,
Esteban
E
,
Felip
E
,
De Angelis
F
, et al
Pembrolizumab plus chemotherapy in metastatic non-small-cell lung cancer
.
N Engl J Med
2018
;
378
:
2078
92
.
8.
Reck
M
,
Rodriguez-Abreu
D
,
Robinson
AG
,
Hui
R
,
Csoszi
T
,
Fulop
A
, et al
Pembrolizumab versus chemotherapy for PD-L1-positive non-small-cell lung cancer
.
N Engl J Med
2016
;
375
:
1823
33
.
9.
Routy
B
,
Le Chatelier
E
,
Derosa
L
,
Duong
CPM
,
Alou
MT
,
Daillere
R
, et al
Gut microbiome influences efficacy of PD-1-based immunotherapy against epithelial tumors
.
Science
2018
;
359
:
91
7
.
10.
Gopalakrishnan
V
,
Spencer
CN
,
Nezi
L
,
Reuben
A
,
Andrews
MC
,
Karpinets
TV
, et al
Gut microbiome modulates response to anti-PD-1 immunotherapy in melanoma patients
.
Science
2018
;
359
:
97
103
.
11.
Matson
V
,
Fessler
J
,
Bao
R
,
Chongsuwat
T
,
Zha
Y
,
Alegre
ML
, et al
The commensal microbiome is associated with anti-PD-1 efficacy in metastatic melanoma patients
.
Science
2018
;
359
:
104
8
.
12.
Sivan
A
,
Corrales
L
,
Hubert
N
,
Williams
JB
,
Aquino-Michaels
K
,
Earley
ZM
, et al
Commensal bifidobacterium promotes antitumor immunity and facilitates anti-PD-L1 efficacy
.
Science
2015
;
350
:
1084
9
.
13.
Vetizou
M
,
Pitt
JM
,
Daillere
R
,
Lepage
P
,
Waldschmitt
N
,
Flament
C
, et al
Anticancer immunotherapy by CTLA-4 blockade relies on the gut microbiota
.
Science
2015
;
350
:
1079
84
.
14.
Dickson
RP
,
Erb-Downward
JR
,
Falkowski
NR
,
Hunter
EM
,
Ashley
SL
,
Huffnagle
GB
. 
The lung microbiota of healthy mice are highly variable, cluster by environment, and reflect variation in baseline lung innate immunity
.
Am J Respir Crit Care Med
2018
;
198
:
497
508
.
15.
Segal
LN
,
Clemente
JC
,
Tsay
JC
,
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
;
1
:
16031
.
16.
Tsay
JJ
,
Wu
BG
,
Badri
MH
,
Clemente
JC
,
Shen
N
,
Meyn
P
, et al
Airway microbiota is associated with upregulation of the PI3K pathway in lung cancer
.
Am J Respir Crit Care Med
2018
;
198
:
1188
98
.
17.
Charlson
ES
,
Bittinger
K
,
Haas
AR
,
Fitzgerald
AS
,
Frank
I
,
Yadav
A
, et al
Topographical continuity of bacterial populations in the healthy human respiratory tract
.
Am J Respir Crit Care Med
2011
;
184
:
957
63
.
18.
Dickson
RP
,
Erb-Downward
JR
,
Freeman
CM
,
McCloskey
L
,
Falkowski
NR
,
Huffnagle
GB
, et al
Bacterial topography of the healthy human lower respiratory tract
.
MBio
2017
;
8
:
e02287
16
.
19.
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
;
1
:
19
.
20.
Gustafson
AM
,
Soldi
R
,
Anderlind
C
,
Scholand
MB
,
Qian
J
,
Zhang
X
, et al
Airway PI3K pathway activation is an early and reversible event in lung cancer development
.
Sci Transl Med
2010
;
2
:
26ra5
.
21.
Greathouse
KL
,
White
JR
,
Vargas
AJ
,
Bliskovsky
VV
,
Beck
JA
,
von Muhlinen
N
, et al
Interaction between the microbiome and TP53 in human lung cancer
.
Genome Biol
2018
;
19
:
123
.
22.
Yoon
SM
,
Shaikh
T
,
Hallman
M
. 
Therapeutic management options for stage III non-small cell lung cancer
.
World J Clin Oncol
2017
;
8
:
1
20
.
23.
Fan
H
,
Shao
ZY
,
Xiao
YY
,
Xie
ZH
,
Chen
W
,
Xie
H
, et al
Incidence and survival of non-small cell lung cancer in Shanghai: a population-based cohort study
.
BMJ Open
2015
;
5
:
e009419
.
24.
Goldstraw
P
,
Crowley
J
,
Chansky
K
,
Giroux
DJ
,
Groome
PA
,
Rami-Porta
R
, et al
The IASLC lung cancer staging project: proposals for the revision of the TNM stage groupings in the forthcoming (seventh) edition of the TNM classification of malignant tumours
.
J Thorac Oncol
2007
;
2
:
706
14
.
25.
Ettinger
DS
,
Wood
DE
,
Aggarwal
C
,
Aisner
DL
,
Akerley
W
,
Bauman
JR
, et al
NCCN guidelines insights: non-small cell lung cancer, version 1.2020
.
J Natl Compr Canc Netw
2019
;
17
:
1464
72
.
26.
Davis
NM
,
Proctor
DM
,
Holmes
SP
,
Relman
DA
,
Callahan
BJ
. 
Simple statistical identification and removal of contaminant sequences in marker-gene and metagenomics data
.
Microbiome
2018
;
6
:
226
.
27.
Schwartz
LH
,
Litiere
S
,
de Vries
E
,
Ford
R
,
Gwyther
S
,
Mandrekar
S
, et al
RECIST 1.1–update and clarification: from the RECIST committee
.
Eur J Cancer
2016
;
62
:
132
7
.
28.
Johnson
M
,
Zaretskaya
I
,
Raytselis
Y
,
Merezhuk
Y
,
McGinnis
S
,
Madden
TL
. 
NCBI BLAST: a better web interface
.
Nucleic Acids Res
2008
;
36
:
W5
9
.
29.
Dickson
RP
,
Erb-Downward
JR
,
Freeman
CM
,
McCloskey
L
,
Beck
JM
,
Huffnagle
GB
, et al
Spatial variation in the healthy human lung microbiome and the adapted island model of lung biogeography
.
Ann Am Thorac Soc
2015
;
12
:
821
30
.
30.
Bassis
CM
,
Erb-Downward
JR
,
Dickson
RP
,
Freeman
CM
,
Schmidt
TM
,
Young
VB
, et al
Analysis of the upper respiratory tract microbiotas as the source of the lung and gastric microbiotas in healthy individuals
.
mBio
2015
;
6
:
e00037
.
31.
Morris
A
,
Beck
JM
,
Schloss
PD
,
Campbell
TB
,
Crothers
K
,
Curtis
JL
, et al
Comparison of the respiratory microbiome in healthy non-smokers and smokers
.
Am J Respir Crit Care Med
2013
;
187
:
1067
75
.
32.
Jin
C
,
Lagoudas
GK
,
Zhao
C
,
Bullman
S
,
Bhutkar
A
,
Hu
B
, et al
Commensal microbiota promote lung cancer development via gammadelta T cells
.
Cell
2019
;
176
:
998
1013
.
33.
Gui
QF
,
Lu
HF
,
Zhang
CX
,
Xu
ZR
,
Yang
YH
. 
Well-balanced commensal microbiota contributes to anti-cancer response in a lung cancer mouse model
.
Genet Mol Res
2015
;
14
:
5642
51
.
34.
Cheng
M
,
Chen
Y
,
Wang
L
,
Chen
W
,
Yang
L
,
Shen
G
, et al
Commensal microbiota maintains alveolar macrophages with a low level of CCL24 production to generate anti-metastatic tumor activity
.
Sci Rep
2017
;
7
:
7471
.
35.
Cheng
M
,
Qian
L
,
Shen
G
,
Bian
G
,
Xu
T
,
Xu
W
, et al
Microbiota modulate tumoral immune surveillance in lung through a gammadeltaT17 immune cell-dependent mechanism
.
Cancer Res
2014
;
74
:
4030
41
.
36.
Lee
SH
,
Sung
JY
,
Yong
D
,
Chun
J
,
Kim
SY
,
Song
JH
, et al
Characterization of microbiome in bronchoalveolar lavage fluid of patients with lung cancer comparing with benign mass like lesions
.
Lung Cancer
2016
;
102
:
89
95
.
37.
Chen
X
,
Wan
J
,
Liu
J
,
Xie
W
,
Diao
X
,
Xu
J
, et al
Increased IL-17-producing cells correlate with poor survival and lymphangiogenesis in NSCLC patients
.
Lung Cancer
2010
;
69
:
348
54
.
38.
Xu
C
,
Hao
K
,
Yu
L
,
Zhang
X
. 
Serum interleukin-17 as a diagnostic and prognostic marker for non-small cell lung cancer
.
Biomarkers
2014
;
19
:
287
90
.
39.
Liao
C
,
Yu
Z
,
Guo
W
,
Liu
Q
,
Wu
Y
,
Li
Y
, et al
Prognostic value of circulating inflammatory factors in non-small cell lung cancer: a systematic review and meta-analysis
.
Cancer Biomark
2014
;
14
:
469
81
.
40.
Tomita
M
,
Shimizu
T
,
Ayabe
T
,
Nakamura
K
,
Onitsuka
T
. 
Elevated preoperative inflammatory markers based on neutrophil-to-lymphocyte ratio and C-reactive protein predict poor survival in resected non-small cell lung cancer
.
Anticancer Res
2012
;
32
:
3535
8
.
41.
Loke
P
,
Allison
JP
. 
PD-L1 and PD-L2 are differentially regulated by Th1 and Th2 cells
.
Proc Natl Acad Sci U S A
2003
;
100
:
5336
41
.
42.
Liu
J
,
Hamrouni
A
,
Wolowiec
D
,
Coiteux
V
,
Kuliczkowski
K
,
Hetuin
D
, et al
Plasma cells from multiple myeloma patients express B7-H1 (PD-L1) and increase expression after stimulation with IFN-gamma and TLR ligands via a MyD88-, TRAF6-, and MEK-dependent pathway
.
Blood
2007
;
110
:
296
304
.
43.
Qian
Y
,
Deng
J
,
Geng
L
,
Xie
H
,
Jiang
G
,
Zhou
L
, et al
TLR4 signaling induces B7-H1 expression through MAPK pathways in bladder cancer cells
.
Cancer Invest
2008
;
26
:
816
21
.
44.
Lee
SK
,
Seo
SH
,
Kim
BS
,
Kim
CD
,
Lee
JH
,
Kang
JS
, et al
IFN-gamma regulates the expression of B7-H1 in dermal fibroblast cells
.
J Dermatol Sci
2005
;
40
:
95
103
.
45.
Chen
J
,
Feng
Y
,
Lu
L
,
Wang
H
,
Dai
L
,
Li
Y
, et al
Interferon-gamma-induced PD-L1 surface expression on human oral squamous carcinoma via PKD2 signal pathway
.
Immunobiology
2012
;
217
:
385
93
.
46.
Akbay
EA
,
Koyama
S
,
Liu
Y
,
Dries
R
,
Bufe
LE
,
Silkes
M
, et al
Interleukin-17A promotes lung tumor progression through neutrophil attraction to tumor sites and mediating resistance to PD-1 blockade
.
J Thorac Oncol
2017
;
12
:
1268
79
.
47.
Guo
B
,
Fu
S
,
Zhang
J
,
Liu
B
,
Li
Z
. 
Targeting inflammasome/IL-1 pathways for cancer immunotherapy
.
Sci Rep
2016
;
6
:
36107
.
48.
Zhong
FL
,
Mamai
O
,
Sborgi
L
,
Boussofara
L
,
Hopkins
R
,
Robinson
K
, et al
Germline NLRP1 mutations cause skin inflammatory and cancer susceptibility syndromes via inflammasome activation
.
Cell
2016
;
167
:
187
202
.
49.
Kolb
R
,
Phan
L
,
Borcherding
N
,
Liu
Y
,
Yuan
F
,
Janowski
AM
, et al
Obesity-associated NLRC4 inflammasome activation drives breast cancer progression
.
Nat Commun
2016
;
7
:
13007
.
50.
Maleki Vareki
S
,
Garrigos
C
,
Duran
I
. 
Biomarkers of response to PD-1/PD-L1 inhibition
.
Crit Rev Oncol Hematol
2017
;
116
:
116
24
.
51.
Blacher
E
,
Levy
M
,
Tatirovsky
E
,
Elinav
E
. 
Microbiome-modulated metabolites at the interface of host immunity
.
J Immunol
2017
;
198
:
572
80
.
52.
Pradhan
D
,
Segal
LN
,
Kulkarni
R
,
Chung
S
,
Rom
WN
,
Weiden
MD
, et al
Bronchial reactivity in early emphysema may be associated with local neutrophilic inflammation
.
Am J Respir Crit Care Med
2013
:
A1110
.
53.
Segal
LN
,
Clemente
JC
,
Wu
BG
,
Wikoff
WR
,
Gao
Z
,
Li
Y
, et al
Randomised, double-blind, placebo-controlled trial with azithromycin selects for anti-inflammatory microbial metabolites in the emphysematous lung
.
Thorax
2017
;
72
:
13
22
.
54.
Dickson
RP
,
Morris
A
. 
Macrolides, inflammation and the lung microbiome: untangling the web of causality
.
Thorax
2017
;
72
:
10
2
.
55.
Lone
AG
,
Atci
E
,
Renslow
R
,
Beyenal
H
,
Noh
S
,
Fransson
B
, et al
Staphylococcus aureus induces hypoxia and cellular damage in porcine dermal explants
.
Infect Immun
2015
;
83
:
2531
41
.
56.
Williamson
KS
,
Richards
LA
,
Perez-Osorio
AC
,
Pitts
B
,
McInnerney
K
,
Stewart
PS
, et al
Heterogeneity in Pseudomonas aeruginosa biofilms includes expression of ribosome hibernation factors in the antibiotic-tolerant subpopulation and hypoxia-induced stress response in the metabolically active population
.
J Bacteriol
2012
;
194
:
2062
73
.
57.
Bourriaud
C
,
Robins
RJ
,
Martin
L
,
Kozlowski
F
,
Tenailleau
E
,
Cherbut
C
, et al
Lactate is mainly fermented to butyrate by human intestinal microfloras but inter-individual variation is evident
.
J Appl Microbiol
2005
;
99
:
201
12
.
58.
Segal
LN
,
Clemente
JC
,
Li
Y
,
Ruan
C
,
Cao
J
,
Danckers
M
, et al
Anaerobic bacterial fermentation products increase tuberculosis risk in antiretroviral-drug-treated HIV patients
.
Cell Host Microbe
2017
;
21
:
530
7
.
59.
Barbi
J
,
Pardoll
D
,
Pan
F
. 
Metabolic control of the Treg/Th17 axis
.
Immunol Rev
2013
;
252
:
52
77
.
60.
Okkenhaug
K
,
Patton
DT
,
Bilancio
A
,
Garcon
F
,
Rowan
WC
,
Vanhaesebroeck
B
. 
The p110delta isoform of phosphoinositide 3-kinase controls clonal expansion and differentiation of Th cells
.
J Immunol
2006
;
177
:
5122
8
.
61.
Sauer
S
,
Bruno
L
,
Hertweck
A
,
Finlay
D
,
Leleu
M
,
Spivakov
M
, et al
T cell receptor signaling controls Foxp3 expression via PI3K, Akt, and mTOR
.
Proc Natl Acad Sci U S A
2008
;
105
:
7797
802
.
62.
Kurebayashi
Y
,
Nagai
S
,
Ikejiri
A
,
Ohtani
M
,
Ichiyama
K
,
Baba
Y
, et al
PI3K-Akt-mTORC1-S6K1/2 axis controls Th17 differentiation by regulating Gfi1 expression and nuclear translocation of RORgamma
.
Cell Rep
2012
;
1
:
360
73
.
63.
Liu
H
,
Yao
S
,
Dann
SM
,
Qin
H
,
Elson
CO
,
Cong
Y
. 
ERK differentially regulates Th17- and Treg-cell development and contributes to the pathogenesis of colitis
.
Eur J Immunol
2013
;
43
:
1716
26
.
64.
Erb-Downward
JR
,
Falkowski
NR
,
D'Souza
JC
,
McCloskey
LM
,
McDonald
RA
,
Brown
CA
, et al
Critical relevance of stochastic effects on low-bacterial-biomass 16S rRNA gene analysis
.
mBio
2020
;
11
:
e00258
20
.
65.
Salter
SJ
,
Cox
MJ
,
Turek
EM
,
Calus
ST
,
Cookson
WO
,
Moffatt
MF
, et al
Reagent and laboratory contamination can critically impact sequence-based microbiome analyses
.
BMC Biol
2014
;
12
:
87
.
66.
Moayyedi
P
,
Surette
MG
,
Kim
PT
,
Libertucci
J
,
Wolfe
M
,
Onischi
C
, et al
Fecal microbiota transplantation induces remission in patients with active ulcerative colitis in a randomized controlled trial
.
Gastroenterology
2015
;
149
:
102
9
.
67.
Paramsothy
S
,
Kamm
MA
,
Kaakoush
NO
,
Walsh
AJ
,
van den Bogaerde
J
,
Samuel
D
, et al
Multidonor intensive faecal microbiota transplantation for active ulcerative colitis: a randomised placebo-controlled trial
.
Lancet
2017
;
389
:
1218
28
.
68.
Caporaso
JG
,
Lauber
CL
,
Walters
WA
,
Berg-Lyons
D
,
Huntley
J
,
Fierer
N
, et al
Ultra-high-throughput microbial community analysis on the Illumina HiSeq and MiSeq platforms
.
ISME J
2012
;
6
:
1621
4
.
69.
Caporaso
JG
,
Kuczynski
J
,
Stombaugh
J
,
Bittinger
K
,
Bushman
FD
,
Costello
EK
, et al
QIIME allows analysis of high-throughput community sequencing data
.
Nat Methods
2010
;
7
:
335
6
.
70.
Holmes
I
,
Harris
K
,
Quince
C
. 
Dirichlet multinomial mixtures: generative models for microbial metagenomics
.
PLoS One
2012
;
7
:
e30126
.
71.
Mortazavi
A
,
Williams
BA
,
McCue
K
,
Schaeffer
L
,
Wold
B
. 
Mapping and quantifying mammalian transcriptomes by RNA-Seq
.
Nat Methods
2008
;
5
:
621
8
.
72.
Wilhelm
BT
,
Marguerat
S
,
Watt
S
,
Schubert
F
,
Wood
V
,
Goodhead
I
, et al
Dynamic repertoire of a eukaryotic transcriptome surveyed at single-nucleotide resolution
.
Nature
2008
;
453
:
1239
43
.
73.
Sultan
M
,
Schulz
MH
,
Richard
H
,
Magen
A
,
Klingenhoff
A
,
Scherf
M
, et al
A global view of gene activity and alternative splicing by deep sequencing of the human transcriptome
.
Science
2008
;
321
:
956
60
.
74.
Tanabe
M
,
Kanehisa
M
. 
Using the KEGG database resource
.
Curr Protoc Bioinformatics
2012
;
Chapter 1:Unit1.12
.
75.
Kanehisa
M
,
Goto
S
,
Sato
Y
,
Furumichi
M
,
Tanabe
M
. 
KEGG for integration and interpretation of large-scale molecular data sets
.
Nucleic Acids Res
2012
;
40
:
D109
14
.
76.
Kramer
A
,
Green
J
,
Pollard
J
,
Tugendreich
S
. 
Causal analysis approaches in Ingenuity Pathway Analysis
.
Bioinformatics
2014
;
30
:
523
30
.
77.
Subramanian
A
,
Tamayo
P
,
Mootha
VK
,
Mukherjee
S
,
Ebert
BL
,
Gillette
MA
, et al
Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles
.
Proc Natl Acad Sci U S A
2005
;
102
:
15545
50
.
78.
DuPage
M
,
Dooley
AL
,
Jacks
T
. 
Conditional mouse lung cancer models using adenoviral or lentiviral delivery of Cre recombinase
.
Nat Protoc
2009
;
4
:
1064
72
.
79.
Romero
R
,
Sayin
VI
,
Davidson
SM
,
Bauer
MR
,
Singh
SX
,
LeBoeuf
SE
, et al
Keap1 loss promotes Kras-driven lung cancer and results in dependence on glutaminolysis
.
Nat Med
2017
;
23
:
1362
8
.
80.
Dobin
A
,
Davis
CA
,
Schlesinger
F
,
Drenkow
J
,
Zaleski
C
,
Jha
S
, et al
STAR: ultrafast universal RNA-seq aligner
.
Bioinformatics
2013
;
29
:
15
21
.
81.
Liao
Y
,
Smyth
GK
,
Shi
W
. 
The Subread aligner: fast, accurate and scalable read mapping by seed-and-vote
.
Nucleic Acids Res
2013
;
41
:
e108
.
82.
Liao
Y
,
Smyth
GK
,
Shi
W
. 
featureCounts: an efficient general purpose program for assigning sequence reads to genomic features
.
Bioinformatics
2014
;
30
:
923
30
.
83.
Friedman
JH
,
Hastie
T
,
Tibshirani
R
. 
Regularization paths for generalized linear models via coordinate descent
.
J Stat Softw
2010
;
33
:
22
.
84.
Cox
DR
. 
Regression models and life-tables
.
J R Stat Soc Series B Methodol
1972
;
34
:
187
220
.
85.
Zhao
N
,
Chen
J
,
Carroll
IM
,
Ringel-Kulka
T
,
Epstein
MP
,
Zhou
H
, et al
Testing in microbiome-profiling studies with MiRKAT, the microbiome regression-based kernel association test
.
Am J Hum Genet
2015
;
96
:
797
807
.
86.
Plantinga
A
,
Zhan
X
,
Zhao
N
,
Chen
J
,
Jenq
RR
,
Wu
MC
. 
MiRKAT-S: a community-level test of association between the microbiota and survival times
.
Microbiome
2017
;
5
:
17
.
87.
Sayers
A
,
Heron
J
,
Smith
A
,
Macdonald-Wallis
C
,
Gilthorpe
MS
,
Steele
F
, et al
Joint modelling compared with two stage methods for analysing longitudinal data and prospective outcomes: a simulation study of childhood growth and BP
.
Stat Methods Med Res
2017
;
26
:
437
52
.
88.
Reiner
A
,
Yekutieli
D
,
Benjamini
Y
. 
Identifying differentially expressed genes using false discovery rate controlling procedures
.
Bioinformatics
2003
;
19
:
368
75
.
89.
Love
MI
,
Huber
W
,
Anders
S
. 
Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2
.
Genome Biol
2014
;
15
:
550
.
90.
Morton
JT
,
Marotz
C
,
Washburne
A
,
Silverman
J
,
Zaramela
LS
,
Edlund
A
, et al
Establishing microbial composition measurement standards with reference frames
.
Nat Commun
2019
;
10
:
2719
.
91.
Morton
JT
,
Aksenov
AA
,
Nothias
LF
,
Foulds
JR
,
Quinn
RA
,
Badri
MH
, et al
Learning representations of microbe-metabolite interactions
.
Nat Methods
2019
;
16
:
1306
14
.
92.
Kodama
Y
,
Shumway
M
,
Leinonen
R
,
International Nucleotide Sequence Database C
. 
The sequence read archive: explosive growth of sequencing data
.
Nucleic Acids Res
2012
;
40
:
D54
6
.
93.
Leinonen
R
,
Sugawara
H
,
Shumway
M
,
International Nucleotide Sequence Database C
. 
The sequence read archive
.
Nucleic Acids Res
2011
;
39
:
D19
21
.