Tumor microbiota of renal cell carcinoma affects clinical prognosis by influencing the tumor immune microenvironment

Abstract

Despite reported influences of the intratumoral microbiome on cancer progression, its role in this subtype remains unclear. This study aimed to characterize the microbial landscape and signatures of kidney renal clear cell carcinoma using RNA-Seq data from The Cancer Genome Atlas. Following microbial decontamination, differential microbial analysis was conducted between tumorous and adjacent non-tumorous samples. Compared to non-tumorous samples, tumorous microbiota exhibited reduced α and β diversity and distinct phylum-level communities. Differential microbial analysis between patients exhibiting long and short overall survival revealed ten significant differential microbial genera, with six genera correlating with a positive prognosis (Plasmodium, Babesia, Toxoplasma, Cytobacillus, Alicyclobacillus, Verrucomicrobium) and four with a negative prognosis (Colletotrichum, Leuconostoc, Gluconobacter, and Parabacteroides). Employing Cox regression analysis and support vector machines, a prognosis-related microbiome risk signature was developed, achieving an AUC of 0.809. Based on this risk signature, two microbiome-based subtypes were found to be significantly associated with distinct clinical prognoses and immune microenvironments. These findings were corroborated by significant correlations between prognostic-relevant microorganisms and 30 immune-related differentially expressed genes. Specifically, microbial genera associated with a negative prognosis were linked to a pro-tumor acute inflammatory immune response, whereas genera related to a positive prognosis were associated with an anti-tumor adaptive immune response. In conclusion, microbiome-based subtyping revealed correlations between tumor microbiome, clinical prognosis, and tumor microenvironment, indicating intratumoral microbiota as a promising prognostic biomarker for kidney renal clear cell carcinoma.

Keywords: Biomarkers; Intratumor microbiota; Kidney renal clear cell carcinoma; Tumor microenvironment.

Fig. 1

Overall profile of the comparison between cancerous and paracancerous tissue microbiotas Abbreviations: Tumor, cancerous samples; Normal, paracancerous samples; PCoA, principal coordinate analysis; KM curve, Kaplan–Meier survival curve; Cox analysis, Cox's proportional hazards regression model analysis. (A) Overlay of the microbial compositions at the phylum and genus levels. The vertical coordinates represent the relative abundance, and colors represents the top 10 microorganisms. (B) Box plots showing the α-diversity (Shannon and Observed species indices) of cancerous and paracancerous samples, Wilcoxon test. (C) PCoA plots showing the β-diversity of cancerous and paracancerous samples, Adonis test. (D) Cladistic map of microbial taxa in cancerous and paracancerous tissues, marked with significantly different genus. (E) KM curves based on risk level. P values were determined using the Cox proportional hazards risk model.

Fig. 2

Analysis of differential microbial groups between patients with different survival times Abbreviations: Long, samples from long-survival patients; Short, samples from short-survival patients; KIRC, kidney renal clear cell carcinoma; PCoA, principal coordinate analysis; ∗, p < 0.05; ∗∗, p < 0.01; ∗∗∗, p < 0.001; ∗∗∗∗, p < 0.0001. (A) Overlay of the microbial compositions at the phylum and genus levels. The vertical coordinates represent the relative abundance, and colors represents the top 10 microorganisms. (B) Box plots showing the Shannon α-diversity of KIRC samples with differences in long- and short-term survival, Wilcoxon test. (C) PCoA plots showing the Jaccard index β-diversity of the KIRC samples with survival differences, Adonis test. (D) Heat map showing the differential expression of microorganisms between different survival groups. (E) Box plots showing the differential expression of microorganisms between different survival groups, Wilcoxon test.

Fig. 3

Survival risk scores based on differentially abundant microorganisms Abbreviations: KM curve, Kaplan–Meier survival curve; Cox analysis, Cox's proportional hazards regression model analysis; ROC, receiver operating characteristic; SVM, support vector machine; AUC, area under the ROC curve; glm, generalized linear model; rpart, decision-making tree. (A) KM curves showing univariate Cox regression analysis results (Plasmodium, Babesia, Toxoplasma, Cytobacillus, Alicyclobacillus, and Verrucomicrobium). p values were determined using Cox proportional hazards risk models. (B) ROC curve demonstrating SVM training for 10 differential microorganisms. (C) ROC curves from three prediction methods: SVM, glm, and rpart. (D) KM curve showing multiple Cox analyses, in which the dataset was randomly divided into training and prediction sets at a 7:3 ratio.

Fig. 4

Differentially abundant microorganisms were related to immune function genes Abbreviations: High: patients in a high-risk group; low: patients in a low-risk group; p-adj, adjusted p value. (A) Volcano plot showing differential gene expression in patients with different risk group. Threshold: log2 fold change >1, p-adj <0.05. (B) Heat map showing 30 differential immune gene expression (24 upregulated and 6 downregulated) in patients with different risk groups. (C) Bubble plots showing the correlation between 30 differential immune gene expression and 10 differential microorganisms. Spearman correlation significant (p < 0.05) coefficients are labeled in the figure.

Fig. 5

Survival risk scores distinguished high- and low-risk patients with distinct immune-related characteristics Abbreviations: p-adj, adjusted p value; GO, Gene Ontology; GSEA, gene set enrichment analysis; NES, normalized enrichment score; FDR, false discovery rate. (A) GO functional annotation pathways of significant upregulated genes (p < 0.05). (B) GSEA graph illustrating the results of upregulated genes in the high-risk group. Screening thresholds: |NES| > 1, p value < 0.05, and FDR (p-adj) < 0.25. (C) GO functional annotation pathways of significant downregulated genes (p < 0.05). (D) GSEA graph illustrating the results of downregulated genes in the high-risk group. (E) GSEA ridge plot with NES as the horizontal coordinate. (F) Heat map showing correlations between survival-related microorganisms and immune cells. Pearson correlation significant (p < 0.05) coefficients are labeled.