Anti-aging effects of a functional food via the action of gut microbiota and metabolites in aging mice

WS improves the weight, food intake, and general physical conditions of aging mice

The hair color and body weight of the mice were recorded to note the base conditions during the dietary intervention process. As shown in Figure 1A, the hair color of the mice in the control group was snow-white, whereas the hair color of the aging mice in the SD, WS, and FMT groups was deep yellow. The mice in the SD group were thinner than the mice in the WS group. During the experimental intervention, the body weight and food intake of mice in the control group were found to be significantly higher than that of mice in the other three groups, and the administration of the WS and FMT resulted in a slight recovery in body weight (Figure 1B, 1C).

Figure 1. The effects of WS on the general nutritional status of aging mice. (A) The hair color and bodyweight of mice in each group. (B) The dynamic weight of mice in each group during the 0-15-week intervention process. (C) The food intake of mice in each group was recorded at days 0, 20, 40, 60, 80, 100, and 120. (D) The level of ALB (g/L) in the serum of mice. (E) The muscle circumference (cm) of the hind limbs of mice. (F) The lean body mass of mice. (G) The fat content in mice. (H) Analysis of body composition of mice. *P < 0.05, **P < 0.01, ***P < 0.001 in a comparison with the labeled group.

Further, the influence of WS on the physical fitness of the mice was examined. We found that the level of serum ALB was markedly decreased in the experimental group compared to the control group (P < 0.001), and the ALB levels were restored in mice from the WS group (Figure 1D). The measurements for the hind limb (Figure 1E) showed a distinct decrease in the muscle circumference in the SD group relative to the control group, whereas the administration of WS resulted in a significant increase in muscle circumference than in the SD group (P < 0.01). Compared to the control group, the lean body mass and fat content significantly decreased in the other three groups, while treatment with WS and FMT caused an increase in the lean body mass and fat percentage (Figure 1F, 1G). Magnetic resonance imaging was used to identify the fat distribution and content in mice from the SD and WS groups, and the subcutaneous fat content was lower in the SD group compared to the WS group (Figure 1H). These results suggested that WS might effectively alleviate the aging-induced reduction in food intake, body weight, serum ALB level, hind limb muscle circumference, lean body mass, and fat content.

WS alleviates cognitive impairments and improves neurological function in SAMP8 mice

The learning and memorization abilities of mice were analyzed using the shuttle-box test. The results showed that the ratio of the successful avoidance response of mice in the control and WS groups was significantly higher than in the other experimental group from the 4^th day of training (Figure 2A, P < 0.05). The retention of the passive avoidance response of mice in the SD group was significantly higher than that in the WS and FMT groups from the 3^rd day of training (Figure 2B, P < 0.05). In addition, mice that were administered FMT and WS showed no significant difference compared with the control group. These results indicated that WS might ameliorate the learning and memorization abilities of SAMP8 mice.

Figure 2. Effects of WS on the learning memory and the neuronal molecules in mice. (A) The ratio of active avoidance response (%). (B) Passive avoidance reaction time (s). *P < 0.05, **P < 0.01, ***P < 0.001 relative to the control group, and ^#P < 0.05, ^##P < 0.01, ^###P < 0.001 relative to the SD group. (C) Images of the HE stained hippocampus in the mouse cerebrum at 100× and 400× magnification. (D) Immunohistochemical images of Aβ expression in tissues from the mouse cerebrum at 20× and 400× magnification. (E) Immunohistochemical images of GFAP expression in tissues from the mouse cerebrum at 50× and 400× magnification.

The hematoxylin-eosin (HE) staining of the hippocampus (Figure 2C) showed that the features of cells in the CA1 region of the hippocampus in mice from the control group were abundant in neurons, and showed a normal morphological structure, clear cytoplasmic boundaries, and an evident nucleolus. However, the cells of mice in SD group were characterized by multiple shriveled neurons, reduced numbers of neurons and cell layers, nuclear hyperchromatism, unclear nucleolus, and rare cytoplasm. The cells from the CA1 region in the WS group showed deeper nuclear staining and no detectable shriveled neurons. The number of neuronal layers in WS group significantly increased relative to that in the SD group. Moreover, several shriveled neurons, a deeper nuclear staining, unclear cytoplasmic and nuclear boundaries, and relatively tightly packed neurons were observed in the cells of the CA2 region from the WS group. Notably, the characteristics of cells in the FMT group were similar to that in the WS group, except for the presence of a greater number of shriveled neurons. These results indicated that WS might improve the aging-induced morphological changes in the hippocampal cells in the cerebrum of mice.

Immunohistochemical (IHC) assays for the hippocampus showed that the expression of glial fibrillary acidic protein (GFAP, Figure 2E) and amyloid β-protein (Aβ, Figure 2D) considerably increased in mice from the SD group than those from the control group, while this effect was ameliorated in mice from the WS and FMT groups. These results demonstrated the effects of WS in the improvement of declining learning memory and neurological impairment in aging mice.

Effects of WS on inflammatory cytokines, oxidation, and antioxidation indicators

The levels of inflammatory cytokines IL-2, IL-4, IL-6, IL-10, TNF-α, and IFN-γ in the serum of mice are shown in Table 1. The results indicated that the level of IL-2 in the control group was significantly higher than that in the SD, WS, and FMT groups (P < 0.001). In contrast, the levels of IL-6, IL-10, and TNF-α in the control group were significantly decreased relative to that in the SD, WS, and FMT groups (P < 0.001). Compared with the SD group, the IL-2 level significantly increased in the WS group (P < 0.001), whereas the IL-6 and TNF-α levels were obviously reduced in the WS group (P < 0.05). However, the levels of IL-2, IL-6, and TNF-α did not show significant differences between the WS and FMT groups. Further, no significant changes in the levels of IL-4 and IFN-γ were found among the four groups. These findings showed that WS might exert a partial anti-inflammatory effect by regulating several inflammation-related cytokines.

Additionally, the indicators for oxidation and antioxidation were also measured in the cerebrum of mice. The results demonstrated that the malondialdehyde (MDA) content in the cerebrum significantly decreased, and superoxide dismutase (SOD) and glutathione peroxidase (GSH-Px) was significantly higher in the control group than in the SD, WS, and FMT groups (Table 2, P < 0.001). Interestingly, the content of MDA in the cerebrum was significantly reduced in the WS group relative to the SD group (P < 0.001). Moreover, a significant activation of SOD activity in the cerebrum was detected in the WS group compared to the SD group. Notably, no significant differences of SOD and MDA activities were observed between the WS and FMT groups. Additionally, no significant change for GSH-Px activity was found in the SD, WS, and FMT groups. Together, these results demonstrated that WS had an antioxidative effect on aging mice.

WS improves immune functions in SAMP8 mice

We investigated the immune cell subsets in the spleen of the mice in the four groups. The results showed that the percentage of lymphocytes, natural killer (NK) cells, CD3⁺CD8⁺T cells, and CD4⁺IFNγ⁺T cells in the aging mice (SD, WS, and FMT groups) were significantly decreased relative to that of the control group. As expected, WS increased the number of these cells in the aging mice, and the percentages of these four cell types were not significantly different between the WS and FMT groups (Table 3). In addition, no significant change in the numbers of B cells, T cells, CD3⁺CD4⁺T cells, CD4⁺IL-4⁺T cells, and CD4⁺IL-17⁺T cells were detected among the four groups.

Influence of WS on bacterial diversity and composition in the gut

We profiled the gut microbiota in fecal samples from each mouse via 16S rRNA gene sequencing. The α-diversity (Chao1 index) analysis revealed that the types of bacterial species and their abundance in the gut microbiota from SD group were significantly decreased relative to mice from the other three groups (Figure 3A, 3B). No significant differences were observed in the gut microbiota among control, WS, and FMT groups, indicating the potential capability of WS to restore the altered imbalance of gut microbial distribution in the aging mice. Mice in the FMT group mimicked the status of intestinal microbiota from the WS group to some extent. In addition, the species accumulation curves demonstrated that the number of microbial species detected showed a leveling off, indicating that sufficient sampling had been achieved and our results were reliable (Figure 3C).

Figure 3. The α- and β-diversity of gut microbiota from feces samples. (A) Chao1 index of gut microbes in mice from the four groups. (B) The species index of gut microbes in mice from four groups. (C) Species accumulation curves. The X-axis represents the number of microbial samples sequenced, whereas the Y-axis represents the number of OTU detected. (D) PCoA plot showing the distribution of feces samples from SAMP8 and SAMR1 mice at the experiment baseline. (E) PCoA plot showing the distributions of plasma samples detected in four groups at the experiment’s metaphase. (F) PCoA plot showing the distributions of plasma samples from SAMP8 and SAMR1 mice at the experimental endpoint. The different colors represent samples from different groups. The closer the samples are, the more similar the microbial composition and structure between the samples. (G) The relative abundance of the top 15 differential OTUs between the WS and SD groups. (H) The heatmap for the top 20 differential OTUs. (I) The boxplot for the top 10 differential OTUs.

Furthermore, the β-diversity of the bacterial community in each group of mice during the baseline, metaphase, and end-point periods were compared. The principal coordinates analysis (PCoA) plots showed a distinct gut microbial community structure in a comparison between the control and SD groups, accounting for 43.62% of the total variability in the data (PC1=22.14%, PC2=10.96%, and PC3=10.52%; Figure 3D). After the administration of WS and FMT, the bacterial community structure in the group of aging mice that were administered SD were found to be increasingly similar to the control group as the intervention time increased (Figure 3E, 3F). In addition, we found that the bacterial community structure of FMT group was highly similar to that of WS group; however, these results were not consistent, indicating that the status of intestinal microbes in mice from the FMT group was simulated from the WS group to a certain extent.

Subsequently, the top 15 relatively abundant microbes in SD and WS groups were identified (Figure 3G), and the top 20 as well as top 10 differential genera are respectively presented in the heatmap (Figure 3H) and boxplot (Figure 3I). Results indicated that the abundance of bacteria from the family Prevotellaceae, and the genera Ruminococcus and Butyrivibrio markedly increased in the WS group relative to the SD group, whereas the abundance of Lachnoclostridium, Ruminiclostridium, Eubacterium coprostanoligenes, Intestinimonas, Clostridium sp. ASF356 and microbes from the family Lachnospiraceae significantly decreased in the WS group compared to the SD group.

Comparison of the number of differential metabolites at different experimental times

The differences in the metabolites level between the SAMP8 and SAMR1 mice at baseline, 8 weeks, and 15 weeks were analyzed. PCA plot revealed that the distribution of metabolites in SAMP8 and SAMR1 mice was clearly separated at the baseline, indicating that the SAMP8 samples had distinct metabolic characteristics relative to the SAMR1 controls. However, the distance between the two sample types gradually approached during the experiment. In brief, the amount of the different metabolites between the SAMP8 and the SAMR1 mice significantly, decreased from 1131 at baseline, to 917 at 8 weeks, and to 111 at 15 weeks (end-point) (Supplementary Figure 1A). Further, no separate distributions for metabolites were observed between the SD and WS groups at 8 weeks; however, distinct distributions were presented at 15 weeks. Specifically, the number of differentially expressed metabolites between the SD and WS groups increased at 15 weeks (887) compared with 8 weeks (56) (Supplementary Figure 2B). Furthermore, for the samples from the WS and FMT groups, the PLS-DA plots showed an approaching distance at 15 weeks compared to 8 weeks, and the number of differentially expressed metabolites between the WS and FMT groups was reduced at 15 weeks (632) than at 8 weeks (375) (Supplementary Figure 2C). These results indicated that the FMT group simulated the alterations of metabolites from the WS group at 15 weeks. However, there was a certain difference in the number of metabolites between the WS and FMT groups, suggesting that the simple intestinal fecal transplantation in this study could not completely replicate the efficacy of oral WS.

Pathways analysis of differentially produced metabolites between the WS and SD groups

For the fecal samples, compared with the SD group, 887 differentially expressed metabolites (eg., vitexin, melatonin, and histamine) were identified in the WS group. The volcano plot of the metabolites is shown in Figure 4A, 4B shows the heatmap for the top 50 differentially expressed metabolites. Based on the VIP values, the top metabolites including 2E, 13Z-octadecadienal, alpha-artemisic acid, and D-glucuronic acid were significantly upregulated in the WS group, whereas 4-Pyridoxic acid, 17-hydroxy-linolenic acid, 5-Hydroxyindoleacetic acid, and (S)-(-)-Perillyl alcohol were significantly downregulated. The log₂ fold-change (FC) values indicated that Tetrahydrofolyl-[Glu](n), diltiazem, candletoxin A, formyl-5-hydroxykynurenamine, ponasterone A, quillaic acid, and vitexin were the top significantly upregulated metabolites in the WS group, whereas propinol adenylate, synaptolepis factor K1, and dehydrosoyasaponin I were the top significantly downregulated metabolites in the WS group (Supplementary Table 3). The results for the pathway enrichment analysis indicated that the differentially expressed metabolites were significantly involved in choline metabolism in cancer, linoleic acid metabolism, and neuroactive ligand-receptor interaction pathways (Figure 4C, 4D).

Figure 4. The enriched functional pathways for the differential metabolites between the WS and SD groups from feces samples. (A) The volcano plot for different metabolites. (B) The heatmap for the top 50 differentially enriched microbes. (C) The bubble map for the top 20 significantly enriched KEGG pathways. (D) The bar graph for the top 10 significantly enriched KEGG pathways.

In addition, 121 differentially produced metabolites were screened from the serum samples of the WS and SD groups. The volcano plot and top 50 metabolites are shown in Supplementary Figure 2A, 2B. Based on the VIP values, the top three metabolites were PA(17:1(9Z)/22:2(13Z,16Z), LysoPC(18:1(11Z), and PC(14:0/22:4(7Z,10Z,13Z,16Z) (Supplementary Table 4). These differential metabolites were significantly associated with protein digestion and absorption, aminoacyl-tRNA biosynthesis, and glycerophospholipid metabolism pathways (Supplementary Figure 2C, 2D).

Correlation analysis between microbial and metabolic profiles

Next, we explored the relationships between the differential microbes and altered metabolites based on Person’s correlation coefficient. A total of 535 microbe-metabolite relationship pairs involving 21 microbes and 94 metabolites were obtained. The heatmap and network of the top 20 correlation pairs are shown in Figure 5A–5C. The results showed that Lactobacillus positively correlated with calenduloside B (r=0.659, P=0.006), soyasaponin III (r=0.591, P=0.016), as well as soyasaponin IV (r=0.625, P =0.010), and negatively correlated with ganoderic acid A (r=-0.506, P=0.046); Butyrivibrio positively correlated with D-glucuronic acid (r=0.499, P=0.049), N-acetyl-a-neuraminic acid (r=0.498, P=0.049), and 5-acetamidovalerate (r=0.867, P<0.001); Ruminococcus positively correlated with D-glucuronic acid (r=0.745, P=0.001) and 5-acetamidovalerate (r=0.513, P=0.042).

Figure 5. The relationship between changed microbiotas and altered metabolites. (A) Correlation square heatmap, (B) circular heatmap, and (C) the network of associations between the differential OTUs of gut microbes and the differentially expressed metabolites. Red color indicates a positive correlation and the blue color represents a negative correlation. A deeper color indicates a greater degree of correlation (*P < 0.05 and **P < 0.01).

Preparation of WS

WS was produced by Yantai Wushen Food Science and Technology Co. Ltd. (Yantai, Shandong, China) using a series of special manufacturing processes, including spray drying, microwave treatment, and low-temperature baking. The proportions of protein, fat, carbohydrate, and fiber in WS were 22.1%, 12.7%, 36.36%, and 21.7%, respectively (Supplementary Table 1). In addition, WS is composed of 55 ingredients, and contains several antioxidant and immune-related nutrients such as crude polysaccharides (12,050 mg/100 g), squalene (164 mg/100 g), procyanidins (244 mg/100 g), anthocyanin (15.87 mg/100 g), vitamin C (24.4 mg/100 g), vitamin E (6.72 mg/100 g), zeaxanthin (236 mg/100 g), selenium (15.87 mg/100 g), and taurine (15.87 mg/100 g) (Supplementary Table 2).

Animals and experimental design

All animal procedures were approved by the Institutional Animal Care and Committee of Second Military Medical University. A total of 60 male SAMP8 mice (aged 3-4 months) and 20 male senescence-accelerated mouse resistant 1 (SAMR1, aged 3-4 months) were purchased from Zhishan (Beijing) Institute of Health Medicine (Beijing, China). SAMP8 mice were allocated into three groups (n=20 per group): standard diet (SD) group, mice had ad libitum access to SD; WS group, mice had ad libitum access to processed WS laboratory feed; and FMT group, mice were fed with SD and were injected in the anus with 200 μL of fecal suspension derived from the WS group. SAMR1 mice fed with SD were used as control group. The SD consisted of 9% water, 19% protein, 4% fat, 5% fiber, 8% ash, and 55% carbohydrate. All mice were allowed to have free water, and maintained at standard temperatures (25 ± 2° C) and relative humidity (40 ± 5 %) conditions under 12 h light/dark cycle. The interventions were continued for 15 consecutive weeks.

Shuttle-box test

Shuttle-box test was performed in accordance with a previous study [48], and it was used to evaluate the learning and memorization abilities of mice. Before the experiment, the mice were placed into chambers (25 ×18.5 ×30 cm) and acclimatized to the new situation for 5 min. Each training consisted of a conditioned stimulus (60 dB, 5 s) followed by an electrical shock (100 V, 50 Hz, AC, 10 s). If the mice escaped to the other side of the shuttle box before electrical stimulation, it is recorded as active avoidance. If the mice completed the shuttle after the electrical shock, it was considered as passive avoidance. Moreover, if no escape occurred after both stimuli, no avoidance behavior was recorded. The test was conducted for 5 consecutive days and 10 times per day with an inter-trial interval of 20 s. After the experiment, the response rates for active avoidance and passive avoidance were calculated to evaluate the memorization ability of mice. Additionally, the time needed for active and passive avoidance evaluated positively with the learning ability of the mice.

Detection of basic indices

The body weight of each mouse was measured every week, and the daily food intake of mice was calculated. In addition, the activity and condition of the hair of the mice in each group were observed and recorded. At the end of the trial, the lean body mass and fat percentage of the mice in each group were examined by using the Awake animal body composition analyzer (MesoQMR23-060H; Shanghai Electronic Technology Co., Ltd, China). Meanwhile, the fat content and its distribution were detected by using the magnetic resonance imaging analyzer MesoMR23-060H-I.

Sample collection

After the experiments, the mice were anesthetized using isoflurane and blood samples were collected. Serum was extracted after centrifugation of the blood samples at 3000 rpm for 15 min and stored at −80° C. Subsequently, the mice were sacrificed and immediately anatomized on ice to obtain the cerebrum, spleen, and hind limbs as described previously [4]. After washing with saline and drying with filter paper, the weights of the harvested organs were recorded. The muscle circumference of the hind limbs of the mice were measured and recorded. The whole hippocampus was isolated from the cerebrum for HE staining and IHC analysis. The feces (50 mg) of the mice at days 0, 8, and 15 weeks were collected and stored at −80° C for the analysis of short-chain fatty acids (SCFA) and the intestinal microflora.

Hematoxylin-eosin staining and immunohistochemical analysis

The isolated hippocampus was immediately fixed upon isolation with 4% paraformaldehyde (Wuhan Servicebio Technology Co., Ltd., Wuhan, China), and subsequently dehydrated in ethanol and embedded in paraffin. Tissue sections of 5 μm thickness were prepared. For HE staining, the sections were stained with hematoxylin (Sigma-Aldrich) for 10 min, followed by incubation with eosin (Sigma-Aldrich) for 30 s at 25° C. The histological morphology of neurons in the CA1 area of the hippocampus was observed using an optical microscope (NIKON Eclipse Ci, Japan) at a magnification of 100× and 400×, respectively. For IHC staining, the sections were subjected to 3 min of microwave heating for antigen retrieval and were subsequently blocked with 5% goat serum (Sigma-Aldrich) at 37° C for 1 h. Next, the sections were incubated with rabbit anti-glial fibrillary acidic protein (GFAP; 1:500; Abcam) primary antibody at 4° C overnight, followed by treatment with horseradish peroxidase (HRP)-conjugated sheep anti-rabbit IgG (1:1500; Abcam) at 37° C for 1 h. The visualization of the immunoreaction was performed using 3, 3′-diaminobenzidine (Sigma-Aldrich). The stained images were photographed using a microscope at 400× magnification.

Quantification of immune cells via flow cytometric analysis

Flow cytometry was used to identify the immune cell subsets from the fresh spleen tissues of mice. Briefly, the spleen tissues were dissociated into single cells by trituration and cell lysis with red blood cell lysis buffer (Cat. No: 420301; Biolegend, San Diego, CA, USA). After centrifugation, the density of the single-cell suspension was adjusted to 1 × 10⁶ cells/mL. The CD3⁺CD8⁺T cells, CD3⁺CD4⁺T cells, CD4⁺IFNγ⁺T cells, and Treg cells were stained with antibodies for labeling. Next, the percentage of the different immune cell subsets was identified by flow cytometry (FACSFortessa X20). The antibodies used were as follows: CD3 eF450 (Cat. No: 11-0042-82, Thermo), CD4 PE-Cy7 (Cat. No: 100528; BD Biosciences, San Jose, CA, USA), CD8a PerCP-Cy5.5 (Cat. No: 45-0081-82; Thermo), CD4 FITC (Cat. No: 11-0042-82; Thermo), and IFN eF450 (Cat. No: 48-7311-82; Thermo).

Enzyme-linked immunosorbent assay (ELISA) analysis

The levels of inflammatory cytokines IL-2 (Cat. No: SU-B20010), IL-4 (Cat. No: SU-B20011), IL-6 (Cat. No: SU-B20012), IL-10 (Cat. No: SU-B20005), TNF-α (Cat. No: SU-B20220), and IFN-γ (Cat. No: SU-B20707) in the mouse serum were detected using ELISA kits from Kenuodi Biotechnology Co., Ltd (Fujian, China) as per manufacturers’ instructions. The oxidation marker MDA, and the antioxidant activities of SOD and GSH-Px in cerebrum samples were detected using ELISA kits from Sigma. Serum albumin (ALB) was detected using a kit from the Nanjing Jiancheng Bioengineering Institute (Nanjing, China).

16S rRNA sequencing

Total microbial genomic DNA (gDNA) from the fecal samples (5 g) of each mouse was extracted using the QIAamp DNA Stool Mini Kit (50) (51504, Qiagen) as per the manufacturer’s protocol. The V3-V4 regions of the bacterial 16S rRNA gene were amplified using the primers pairs 356F (5'-CCTACGGGNGGCWGCAG-3') and 803R (5'-GACTACHVGGGTATCTA ATCC-3') [49]. Next, the PCR amplification, the library preparation and sequencing were performed as described previously [49]. The Illumina Miseq (Illumina, San Diego, CA, USA) sequencing platform was used to generate paired-end sequencing reads (2×300 bp). Subsequently, the FLASH software (Trimmomatic) was used to merge the matched pairs and perform quality control for the reads obtained.

Bioinformatics analysis

USEARCH (version 7.0) was used for Operational Taxonomic Units (OTU) clustering at a similarity level of 97%. The Ribosomal Database Project classifier v2.2 was applied to perform the taxonomy for each representative OTU against the Silva database [50]. Then, the obtained OTU data were used for taxonomical assignments using RDP Classifier algorithm (http://rdp.cme.msu.edu/classifier/classifier.jsp). The analyses for the richness index, alpha (Chao index [51]), and beta diversity (PCoA) were conducted through QIIME (Version 1.9, http://qiime.sourceforge.net/). The t-test and the Wilcoxon sum-rank test were used to assess the difference between the groups. The OTUs with P<0.05 were considered to be significantly different between groups. A principal component analysis (PCA) and UniFrac-based PCoA were performed using the R statistical software to perform the clustering of different samples.

Metabolomic analysis

A total of 100 μL serum and 60 mg stool samples from each mouse in different groups (n=8) were collected. Subsequently, 10 μL and 20 μL of the internal standard (0.3 mg/mL L-2-chlorophenyl alanine and 0.01 mg/mL 17:0 Lyso PC) was added to the serum and stool samples, respectively, and 400 μL mixture of methanol and water (V:V=1:4) was added to remove impurities. The supernatant of each sample was used for liquid chromatography-mass spectrometry (LC-MS) analysis using Dionex UltiMate 3000 UHPLC, (Thermo Scientific). Chromatographic separation was conducted using a 100 × 2.1 mm, 1.8 μm ACQUITY UPLC HSS T3 column at 50° C with a flow rate of 0.35 mL/min. Water with 0.1% formic acid was used as the solvent A of the mobile phase, and acetonitrile with 0.1% formic acid was used as solvent B.

Data pre-processing, multivariate analysis, and metabolic annotation

The Progenesis QI v2.3 software was used for data processing, including baseline filtering, peak identification, integration, correction of retention time, peak alignment, and normalization. Next, multivariate statistical analyses including PCA and orthogonal partial least-squares discriminant analysis (OPLS-DA) were used to observe the distribution of samples and the differential metabolites among groups. The degree of difference between the two groups was evaluated using the t-test, and the metabolites showing significant differences were selected using thresholds of VIP > 1 (VIP = variable importance in the projection) and P-values < 0.05. Volcano plots were created to identify the P and the fold-change values for metabolites between the groups. Next, hierarchical clustering was performed to identify the expression differences of the top 50 metabolites in different samples. The functions of differential metabolites were annotated using the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis.

Correlation analysis for microbial taxonomy and metabolites

To evaluate the potential associations between the different microbes and metabolites identified in the WS and SD groups, a Spearman’s correlation analysis was conducted. The Spearman’s correlation coefficient was calculated for the relative abundance of OTU and the response intensity of the metabolites, and the relationship between microbial taxonomy and metabolites. A network showing the top 20 results was plotted.

Statistical analysis

Data were expressed as mean ± standard deviation (mean ± SD) and were analyzed using the SPSS 19.0 software (IBM, New York, NY, USA). Data from more than three groups were analyzed using a one-way analysis of variance (ANOVA) and the least significant difference (LSD) test. A value of P < 0.05 was considered as statistical significant.