A clinical trial with 30 ccRCC patients (RCC group) and 30 healthy controls (Control group) was conducted. All ccRCC patients in the RCC group were diagnosed exactly by pathological examinations after surgery, and all healthy controls were recruited from Beijing Friendship Hospital, Capital Medical University (Beijing, China). The healthy controls at this stage were free of renal cell carcinoma by medical examinations and had no history before enrollment. All subjects were Chinese Han ethnic population to ensure their living and eating habits were basically similar. All experimental subjects conformed to the following standards. Inclusion Criteria: (1) no infectious disease or abnormal renal function, no neurological symptoms, no upper gastrointestinal disease, (2) performance status (PS) score ranging from 0 to 1.5, proper functioning of major organs, no other active malignancies (except non-melanoma skin cancer), (3) no antibiotics, probiotics, vitamins, minerals, NSAIDs, prokinetics, bismuth, antacids, H2-receptor antagonists, omeprazole, Sucralfate, or misoprostol intake within the past 4 weeks, no recent hormone therapy taken, and (4) no history of gastroduodenal ulcer or major gastrointestinal surgery. Exclusion criteria: (1) patients with gastrointestinal discomfort such as acid reflux and nausea, diabetes or other serious systemic diseases, (2) HIV-infected or hepatitis B surface antigen positive patients, (3) patients with incomplete heart function, unstable angina pectoris or recent (within 6 months) myocardial infarction, (4) patients with history of major organ transplant, radiotherapy or chemotherapy, and (5) females who are pregnant or breast-feeding.

The study was approved by the Ethics Committee of Beijing Friendship Hospital Affiliated to Capital Medical University (Beijing, China) (Document number: 2021-P2-072-01). All subjects have signed the informed consent.

Patients in the RCC group were all initially diagnosed as renal carcinoma by CT/MRI scan. After being informed about the purpose and procedures of the study, sufficient, middle, and internal fecal samples were freshly collected from both of the RCC group and Control group by the laboratory personnel with special fecal sampling tubes, and samples were transported to the laboratory on ice within 1 h (Xing Kang Medical, Jiangsu). Fecal samples of patients in the RCC Group were all collected before surgery, and urine and urinal pollution to the fecal samples were avoided during collection processes. Samples were stored at −80°C until further processing.

A fasting blood sample of all subjects were taken before receiving any medical treatments or surgery. Serum specimens were separated by centrifugation for a thorough analysis, including blood routine test, C-reactive protein test, liver and kidney function test, as well as electrolyte detection.

Note that the fecal and blood sample used for experiments and analysis are from patients in RCC group that have postoperative pathological diagnosis confirmation of ccRCC.

Total microbial DNA was extracted from the fecal samples using the OMEGA Soil DNA Kit (Omega BioTek, Norcross, GA, USA) according to the manufacturer’s instructions. The DNA concentration were evaluated in a NanoDrop ND-1000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA), respectively, and all DNA samples were stored at −80°C prior to further analysis. The A260:A280 ratios were ranging from 1.8 to 2.0.

PCR amplification of the microbial 16S rRNA V3–V4 hypervariable regions was performed. PCR primers of the V3–V4 hypervariable region were designed as previously described (Wang et al., 2016). The forward primer used was 5′-TCGTCGGCAGC GTCAGATGTGTATAAGAGACAGCCTACGGGNGGCWGCAG-3′, and the reverse primer was 5′-GTCTCGTGGGCTCGGAG ATGTGTATAAGAGACAGGACTACHVGG GTATCTAATCC-3′. Amplification was performed on an ABI 9600 instrument according to the instruction successively with a denaturation at 94°C for 4 min, 25 cycles of denaturation at 94°C for 45 s, annealing at 60°C for 45 s, elongation at 72°C for 45 s and a final extension cycle at 72°C for 8 min.

Gut microbial 16S rRNA gene V3–V4 region was sequenced using Illumina MiSeq platform as previously described (Wang et al., 2016). The high-quality data was obtained after quality control processing. Bioinformatics analysis of sequencing data was performed using QIIME2 (Bolyen et al., 2019).

After removal of primers with Cutadapt, the sequences were denoised, merged and chimera were removed with DADA2. The remaining high-quality sequences were clustered into amplicon sequence variants (ASVs). A representative sequence from each ASV was assigned taxonomically against the SILVA database. Next, an original ASV table was created which contained a readable matrix of the ASV abundance for each sample. Alpha diversity based on Shannon, Simpson, Chao 1, Observed species, Faith’s PD, Pielou’s evenness, and Good’s coverage were estimated. Beta-diversity was calculated based on Euclidean distance, and the significance of difference between the two groups was assessed by permutational multivariate analysis of variance (PERMANOVA) test in MATLAB R2018b (The MathWorks, Inc., Natick, MA, USA). Random Forest and linear discriminant analysis effect size (LEfSe) were contrasted based on the relative abundance of ASVs. The top 20 most important ASVs were shown in Random Forest. The ASVs with significant changes were found out in LEfSe when alpha value of the factorial Kruskal–Wallis test was <0.05 and the logarithmic LDA score was >2.0. The associations between the differential ASVs and the blood biochemical and metabolic indexes were performed by Spearman correlation.

Blood samples were collected in Vacutainer tubes containing the chelating agent ethylene diamine tetraacetic acid (EDTA). The samples were then centrifuged at 1,500 g, 4°C for 15 min. Plasma samples were stored at −80°C prior to further analysis.

A total of 200 μl of plasma samples were mixed with 400 μl of cold methanol/acetonitrile (1:1, v/v) to remove the protein. The mixture was centrifuged 14,000 g, 4°C for 15 min. and the supernatants were dried in a vacuum centrifuge. The samples were re-dissolved in 100 μl acetonitrile/water (1:1, v/v) solvent.

Analysis was performed using an UHPLC (1290 Infinity LC, Agilent Technologies) coupled to a quadrupole time-of-flight (AB Sciex TripleTOF 6600). Samples were analyzed by HILIC separation. A 2.1 mm × 100 mm ACQUIY UPLC HSS T3 1.8 μm column (waters, Ireland) was used. In ESI positive mode, the mobile phase contained A = water with 0.1% formic acid and B = acetonitrile with 0.1% formic acid; and in ESI negative mode, the mobile phase contained A = 0.5 mM ammonium fluoride in water and B = acetonitrile. The gradient was 1% B for 1.5 min and was linearly increased to 99% in 11.5 min and kept for 3.5 min. then it was reduced to 1% in 0.1 min and a 3.4 min of re-equilibration period was employed. The gradients were at a flow rate of 0.3 ml/min, and the column temperatures were kept constant at 25°C. 2 μl aliquot of each sample was injected.

The ESI source conditions were set as follows: Ion Source Gas1 (Gas1): 60, Ion Source Gas2 (Gas2): 60, curtain gas (CUR): 30, source temperature: 600°C, IonSpray Voltage Floating (ISVF): ±5,500 V; TOF MS scan m/z range: 60–1,000 Da, product ion scan accumulation time: 0.20 s/spectra; TOF MS/MS scan m/z range: 25–1,000 Da, product ion scan accumulation time: 0.05 s/spectra. The product ion scan was acquired using information dependent acquisition (IDA) with high sensitivity mode selected. The parameters were set as follows: collision energy (CE): 35 V ± 15 eV; declustering potential (DP): ±60 V; exclude isotopes within 4 Da; candidate ions to monitor per cycle: 10.

The raw MS data were converted from wiff.scan files to MzXML files using ProteoWizard MSConvert, then data were imported into freely available SCMS software for peak picking and grouping. Collection of Algorithms of Metabolite pRofile Annotation (CAMERA) were used for annotation of isotopes and adducts. Compound identification of metabolites was performed by comparing of accuracy m/z value (<25 ppm), and MS/MS spectra with an in-house database established with available authentic standard.

After normalized to total peak intensity, the data were then analyzed suing R package (v3.2.0) (Chen et al., 2018). Univariate analysis and orthogonal partial least-squares discriminant analysis (OPLS-DA) were performed to find out the differential metabolites between the two groups. The sevenfold cross-validation and response permutation testing was used to evaluate the robustness of the model. The variable importance in the projection (VIP) value of each variable in the OPLS-DA model was calculated to indicate its contribution to the classification. Metabolites with the VIP value >1 was further applied to Student’s t-test at univariate level to measure the significance of each metabolite, the P-values less than 0.05 were considered as statistically significant. After combining the differential metabolites screened by positive and negative ion modes, the KEGG pathway was annotated and analyzed.

In brief, public gene expression data of ccRCC patients was downloaded from the Cancer Genome Atlas (TCGA) database. A well-established metabolic signature gene set was adopted in a previous study (Peng et al., 2018), which included genes for amino acid metabolism (348 genes), carbohydrate metabolism (286 genes), energy metabolism (110 genes), lipid metabolism (766 genes), nucleotide metabolism (90 genes), TCA cycle (148 genes), and vitamin cofactor metabolism (168 genes). A total number of 200 genes that related to inflammatory response were also collected. Differentially expression genes between tumor and normal tissues were detected by “limma” R package with fold change >2 and false discovery rate <0.05. In addition, we performed the prognostic analysis for selected genes by using univariate Cox regression model. Finally, all the data was visualized by different R packages.

Differences between two groups were compared by Student’s t-test or Mann–Whitney U test. All statistical tests were performed using GraphPad Prism Software. Data were considered significant at P < 0.05.

Sixty fecal samples were sequenced from the RCC group (n = 30) and the Control group (n = 30). Baseline clinical parameters of subjects were listed in Table 1, with no significant difference between the RCC group and the Control group (P > 0.05).

1,800,840 valid sequences were obtained from the 60 qualified fecal samples through sequencing and data processing. The Shannon and Simpson indices were shown in Figure 1A. Results indicated that there was no difference in the richness and evenness of gut microbiota between the RCC group and the Control group. More alpha diversity indices such as Chao 1, Observed species, Faith’s PD, Pielou’s evenness, and Good’s coverage confirmed the same result (Supplementary Figure 1). Beta-diversity analysis revealed there existed a significant difference in the overall structure of the gut microbiota between the RCC group and Control group (Figure 1B).

To further clarify the features of gut microbiota in patients with ccRCC, we analyzed the relative abundance of specific microbiota, and classified them at different taxonomic levels. The abundance distribution of 10 main gut microbiota, including Firmicutes, Bacteroidetes, Proteobacteria, Actinobacteria, Verrucomicrobia, Tenericutes, Fusobacteria, Cyanobacteria, TM7, and Synergistetes were described at phyla level (Figure 1C). Detailed taxonomic abundance difference of gut microbiota at family and phyla levels was shown according to the sequencing dataset (Supplementary Figure 2). Moreover, significant difference of 18 genera between the RCC group and the Control group is identified by Random Forest analysis (Figure 1D). Results indicated that 11 species were enriched in ccRCC patients, including Bradyrhizobiaceae, Desulfovibrionaceae, two kinds of Enterobacteriaceae, Streptophyta, Peptostreptococcaceae, Allobaculum, Candidatus Aquiluna rubra, Ruminococcus bromii, Akkermansia muciniphila, and Collinsella aerofaciens, On the other hand, nine species including Wolbachia, four kinds of Lactobacillus, Bacteroides, Clostridium perfringens, Prevotella copri, and Shigella, were down-regulated notably. Moreover, it can be concluded from the results that Wolbachia, Bradyrhizobiaceae, and Lactobacillus iners differed most typically based on the Importance score. LEfSe was then adopted to analyze the significant microbial variations between the two groups, too (Supplementary Figure 3). Our data preliminarily suggested that patients with ccRCC had alterations of specific bacterial functional species compared with those of the Control group.

To pursue our primary hypothesis that changes of gut microbiota are associated with ccRCC tumorigenesis, we examined the potential related differential genes involved in serum metabolites of the two groups. Except for the undefined superclass, the “organic acids and derivatives” accounted for 19.832% of all relevant metabolites with the largest proportion, according to the serum metabolomic analysis, followed by “lipids and lipid-like molecules” (Supplementary Figure 4A). Volcano plots indicated that 256 altered metabolites in positive ion mode and 102 altered metabolites in negative ion mode are found in the comparison between the RCC group and the Control group, and most of these metabolites fell into “organic acids and derivatives” and “lipids and lipid-like molecules” classes (Supplementary Figures 4B, C).

A superior OPLS-DA model (Q² > 0.5 both in positive and negative ion mode) was generated to analyze the difference of serum metabolomic profiles (Figures 2A, B). As shown in Figures 2C, D, majority of the 18 major metabolites altered significantly in positive ion mode were classified as “organic acids and derivatives” based on superclass, and majority of the 7 metabolites in negative ion mode were classified as “lipids and lipid-like molecules.” Heatmap was also used to intuitively illustrate the abundance distribution of the 25 metabolites in all fecal samples screened by OPLS-DA, respectively, in positive ion mode (Figure 2E) and negative ion mode (Figure 2F). In general, 15 metabolites (such as taurine) were downregulated and 10 metabolites (for example, arachidonic acid) were upregulated in the RCC group. In addition, we also found from Figures 1D, 2E that the abundance of Desulfovibrionaceae in high-grade ccRCC (T3) patients was significantly higher than that in low-grade (T1) patients, while the level of serum taurine was lower.

In order to explore the potential pathways by which the above altered metabolites function, a gene-KEGG pathways dot plot was drawn (Figure 2G). The results showed that the identified metabolites involved 20 metabolic pathways, among which protein digestion and absorption, aminoacyl-tRNA biosynthesis, and ATP binding cassette (ABC) transporter pathways showed the strongest correlation. Meanwhile, a spearman analysis was conducted to investigate the patterns of interactions between metabolites that might play a role in positive and negative ion mode (Supplementary Figure 5). In the heatmap, the abundance of particular metabolites, such as taurine and arachidonic acid was further analyzed and demonstrated (Supplementary Figure 6).

Furthermore, clinical indices including blood routine, C-reactive protein, liver and kidney function, and electrolyte analysis were performed to evaluate the clinical status difference between the RCC group and the Control group (Figure 3A and Supplementary Figure 7). As shown in Figure 3, inflammation related indices such as WBC, CRP, neutrophil, and MO, were significantly upregulated in the RCC group. We investigated whether the 20 microbial species differing between the RCC group and the Control group were correlated with well-established clinical indices and altered serum metabolites (Figure 3B). Remarkably in the two groups, Lactobacillus, known as protective bacteria, was negatively correlated with inflammatory indices such as WBC, MO, and CRP. Until now, we could draw the conclusion that there exists a correlation between the altered gut microbiota and inflammation level in patients with ccRCC.

Considering that the significant variation in a subset of gut microbiota might be closely associated with inflammation in the RCC group, we search for the possible link between inflammation encoding genes and ccRCC development. Result showed that 22 genes, including Gabbr1 and etc., were identified as putative targets, which were overlapped across the univariate Cox regression and different expression analysis (Figure 4A and Supplementary Table 1). We then estimated the influence of the 22 inflammatory factors on the prognosis of ccRCC (Figure 4B). Strikingly, Gabbr1 was also a risk factor associated with poor prognosis of ccRCC.

We also explored serum metabolites associated with ccRCC through TCGA database, in order to conduct further combined analysis with previous metabolomics results. The Venn diagram illustrated the overlap of 202 metabolism-related transcripts both identified by the univariate Cox regression and different expression analysis (Figure 4C and Supplementary Table 2). A total of 202 metabolism-related genes were filtrated from different expression analysis with significant prognosis (Figure 4D). Finally, the 202 selected metabolism-related genes with significant different expression and prognosis were clustered into 7 metabolism processes, including lipid (42%), amino acid (19%), carbohydrate (17%), vitamin cofactor (8%), nucleotide (5%), energy (5%), and TCA cycle (4%) (Figure 4E). Results confirmed the crucial role of lipid metabolism in ccRCC tumorigenesis and development, which showed strong consistency with our metabolomics results.