Parishin was provided by Qi Jianhua’s research team at Zhejiang University (Zhejiang University Patent No. CN201610061288.2). It was dissolved in normal saline and prepared as a solution with a concentration of 2 mg/ml.

SPF-grade C57BL/6 male adult mice (12 weeks, weight 20 ± 2 g) and aged mice (19 months, 20–30 g) were provided by Zhejiang Experimental Animal Center, and all the mice were fed a normal diet and housed at room temperature with 12/12 h day/night alternation. Due to the lack of in vivo studies on the antiaging effects of parishin, the dose setting of parishin in this study was based on the results of in vitro experiments (Lin et al., 2016) and with reference to the reported study of gastrodin in the long-term treatment of natural aging-related disease (Huang et al., 2015). After 2 weeks of acclimatization, the mice were randomly divided into four groups (N = 10, 5 per cage) as follows (Charan and Kantharia, 2013): the young group (YC), adult mice were gavaged with saline; the old group (AC), aged mice were gavaged with saline; the low-dose group (LPar), aged mice were gavaged with parishin (10 mg/kg/d); and the high-dose group (HPar), aged mice were gavaged with parishin (20 mg/kg/d). The treatment lasted for 8 weeks, once daily, and the feces of all the groups were collected at the end of 8 weeks. After completion of the treatment, the mice were anesthetized with pentobarbital sodium before sacrifice. Eye blood was taken for hematological testing, and heart and lung tissues were taken for histopathological evaluation. All the experimental procedures followed the National Institutes of Health (NIH) Guide for the Care and Use of Laboratory Animals. The study protocol was approved by the Animal Experimentation Ethics Committee of Zhejiang University (no. 1446).

The blood samples were allowed to stand for 30 min at room temperature, and the serum was separated by centrifugation at 3,000 rpm for 10 min. The serum levels of GDF15, p16^Ink4a, and IL-6 were measured by ELISA (Cusabio Biotech, Wuhan, China).

Fecal DNA was extracted using the QIAamp Fast DNA Stool Mini Kit (Qiagen, Valencia, United States). PCR amplification of the 16S rDNA variable region (V3 + V4) was performed using primers 338F (5′-ACTCCTACGGGAGGCAGCAG-3′) and 806R (5′-GGACTACHVGGGTWTCTAAT-3′; Lv et al., 2021b). The PCR products were purified using AMPure XPbeeds (Agencourt, Beckman Coulter, United States), and the library was quantified by real-time quantitative PCR. Paired-end sequencing (2 × 300 bp) was performed using the Illumina MiSeq platform (Illumina, San Diego, CA). FLASH (v1.2.8) software was used to clean, filter and merge the raw reads. Vsearch (v2.3.4) was then used to select operational taxonomic units (OTUs) with > 97% sequence similarity. OTU clustering and identification based on the RDP database and the NCBI-16S database as well as microbial diversity and differential enrichment were analyzed using QIIME (v1.8.0). PICRUSt2 was used to predict the gene function of the obtained OTUs.

Sample preparation for metabolic assays was performed as described previously (Lv et al., 2021a). Briefly, a mixture of 20 mg feces and 800 μL methanol was homogenized three times using a Precellys Evolution instrument (Bertin Technologies, United States). After centrifugation at 14,000 rpm for 15 min, the supernatant was filtered through a 0.22 μm membrane, and 20 μL of heptadecanoic acid (1 mg/mL, Sigma–Aldrich, St. Louis, MO, United States) was added as an internal reference and then dried under nitrogen at room temperature. After drying, the samples were methoxymated with methoxypyridine (20 mg/mL, Sigma–Aldrich, St. Louis, MO, United States) and trimethylsilylated with N,O-bistrifluoroacetamide containing 1% trimethylsilyl chloride. The pretreated samples were analyzed with an Agilent 7890A-5975C GC–MS system (Agilent, United States). The downstream data were compared with NIST 17 databases to identify the corresponding metabolites (matching score ≥ 80%).

The mouse heart and lung tissues were fixed in 4% paraformaldehyde (Beijing Solarbio Technology Co., Ltd.) at room temperature for 24 hours and then embedded in paraffin and cut into 4-μm-thick sections. Hematoxylin and eosin (H&E) and Masson’s trichrome staining were used to visualize the architecture and fibrosis of cardiac and lung tissues, respectively. The slices were examined under an inverted light microscope (Leica, Berlin, Germany). The degree of fibrosis was quantified by ImageJ software and expressed as a percentage of the fibrotic area of the whole region.

For the comparison of serum GDF15, p16^Ink4a, and IL-6 levels, quantitative cardiac and lung histopathology results and intestinal bacterial alpha diversity and fecal metabolites between groups, the Shapiro–Wilk method was first used to test whether the data of each group conformed to a normal distribution. Linear mixed models was used to compare any two data sets that were normally distributed; otherwise, the Mann–Whitney U test was used. The Wilcoxon rank sum test combined with the Benjamini–Hochberg method was used to compare the relative abundance of each taxonomic level of intestinal bacteria between groups. Correlation analysis between variables was performed using Spearman’s rank correlation test. P < 0.05 was considered a statistically significant difference.

To investigate the ameliorative effect of parishin on the aging phenotype in mice, we examined changes in some important aging-related physiopathological indicators after parishin treatment. The ELISA results showed that the serum p16^Ink4a, GDF15, and IL-6 levels were significantly lower in the parishin-treated mice (especially in the HPar group) than in the AC group mice (Figure 1A). In the lung tissue, H&E staining showed that the alveoli were collapsed, the disturbed alveolar structure, inflammatory cells infiltration and fibroblasts proliferated were observed in the AC group mice compared with the YC group mice; In the LPar and HPar groups, the morphological structure of the lung tissue was improved, the number of infiltrating inflammatory cells was decreased, and alveolar structure was restored. Masson staining showed a large number of blue-stained collagen fibers in the AC group. In the parishin treatment groups, collagen deposition was reduced (Figure 1B). In the heart tissue, H&E staining showed that cardiomyocytes presented a disordered arrangement in the AC group compared to the YC group. The cardiomyocyte arrangement became significantly ordered in the LPar and HPar mice, especially in the HPar treatment group. Masson staining showed a large area of collagen accumulation in the AC group, which indicated a high level of fibrosis, while parishin treatment alleviated collagen accumulation in the LPar and HPar groups. Additionally, the fibrosis degree analysis confirmed that parishin treatment improved the cardiopulmonary fibrosis in aged mice (Figure 1C).

Changes in structure and function of gut microbiota accompany the aging of gut tissue. To investigate the regulatory effect of parishin on the gut microbiota, we quantitatively determined and analyzed the gut bacteria of each group of mice based on 16S rDNA sequencing. A total of 1,304,626 reads from 40 fecal samples were obtained by 16S rDNA sequencing. Alpha diversity analysis showed that the community biodiversity indicated by the Shannon index and the community richness indicated by the Chao 1 index were not significantly different between the YC and AC groups, while the community biodiversity and community richness in the HPar and LPar groups were significantly changed compared to those in the YC and AC groups (Figure 2A). For beta diversity, microbial community dissimilarities among the four groups of mice were visualized by principal coordinate analysis (PCoA) and non-metric multidimensional scaling (NMDS) plots (Figure 2B), and the significance of intergroup differences between groups was determined by analysis of similarities (ANOSIM) and Adonis. The results showed significant differences in all comparisons except HPar vs. LPar. (Table 1).

We further analyzed the changes in the relative abundance of individual bacterial taxa between the YC, AC, HPar and LPar groups. The phylum Firmicutes and genera Limosilactobacillus and Saccharimonas were enriched in the AC group compared with the YC group. The abundance of bacterial taxa from Firmicutes or Limosilactobacillus showed an increasing trend after parishin treatment compared to the AC group, while the abundance of Saccharimonas in the LPar group decreased significantly compared to that in the AC group (P_adj = 0.016). In addition, the abundance of bacterial taxa from families Atopobiaceae and Peptostreptococcaceae and genera Prevotellaceae NK3B31 group and Rombousia were enriched in the LPar and/or HPar group compared to the AC group, while the families Tannerellaceae, Lachnospiraceae, and Beijerinckiaceae and the genera Erysipelatoclostridium, Defluviitaleaceae UCG-011, Turicibacter, Rikenellaceae RC9 gut group, Parabacteroides, Oscillospiraceae sp., Oscillibacter, Lachnospiraceae sp., Desulfovibrio and an uncultured Clostridum were depleted (Figure 2C).

To explore the possible health-promoting role of gut microbial alterations induced by parishin teatment, we predicted the functions of the differential gut microbes between the YC, AC, LPar and HPar groups by using Picrust 2 software and then mapped the differential functional genes to the MetaCyc pathway database and clustered them according to abundance. As shown in Figure 3, these differentially expressed functional genes were mainly enriched in five classes of cell wall biosynthesis, generation of precursor metabolites and energy and degradation/utilization/assimilation, fatty acid and lipid biosynthesis, amino acid biosynthesis, and nucleoside and nucleotide biosynthesis for a total of 35 pathways (P < 0.05).

We further analyzed the functional differences in the differentially enriched gut microbes in each comparison. In the comparison of AC and YC, differential functional genes located in nucleoside and nucleotide biosynthesis-related pathways, including the superpathway of pyrimidine deoxyribonucleotide de novo biosynthesis, superpathway of pyrimidine ribonucleosides salvage, superpathway of pyrimidine deoxyribonucleoside salvage, pyrimidine deoxyribonucleotide de novo biosynthesis II and inosine-5’-phosphate biosynthesis III pathways, and genes located in the generation of precursor metabolites and energy and degradation/utilization/assimilation-related pathways, including bifidobacterial shunts, pentose phosphate pathway and heterolactic fermentation pathway, were significantly enriched; however, genes located in the superpathway of β-D-glucuronide and D-glucuronide degradation, superpathway of hexuronide and hexuronate degradation, TCA cycle VIII and L-1,2-propanediol degradation pathway, and an amino acid biosynthesis-related pathway, L-arginine biosynthesis III, were significantly reduced (Figure 4A).

After high-dose parishin treatment, all nucleoside and nucleotide biosynthesis-related pathways in which the differentially expressed functional genes that located in the AC vs. YC were no longer significant. In addition, the differentially expressed genes that were significantly downregulated in the HPar group were mainly located in fourteen cofactor, carrier and vitamin biosynthesis-, precursor metabolite and energy production and degradation/utilization/assimilation-, and sugar nucleotide biosynthesis-related pathways, such as biotin biosynthesis I, L-tyrosine degradation I, and CMP-legionaminate biosynthesis I, while the significantly upregulated genes in the HPar group were located in the O-antigen building blocks biosynthesis pathway and fatty acid and lipid biosynthesis-related pathways, such as CDP-diacylglycerol biosynthesis I and II, and superpathways of phospholipid biosynthesis I (Figure 4B).

However, the LPar group exhibited upregulated differentially expressed genes in four nucleoside and nucleotide biosynthesis-related pathways, including the superpathways of adenosine nucleotide de novo biosynthesis I and II, guanosine deoxyribonucleotide de novo biosynthesis II, and adenosine deoxyribonucleotide de novo biosynthesis II, compared to the AC group, while the downregulated genes were located in the superpathway of pyrimidine deoxyribonucleotide de novo biosynthesis (Figure 4C). On the other hand, similar to the those associated with the high-dose parishin treatment, the upregulated genes associated with the low-dose parishin treatment were also located in the O-antigen building block biosynthesis pathway and in three fatty acid and lipid biosynthesis-related pathways, namely, phosphatidylglycerol biosynthesis I and II and the superpathway of phospholipid biosynthesis I. Additionally, the altered genes associated with low-dose parishin treatment were mainly located in seventeen cell wall biosynthesis-, cofactor, carrier, and vitamin biosynthesis-, generation of precursor metabolites and energy and degradation/utilization/assimilation-, amino acid biosynthesis-, and sugar nucleotide biosynthesis-related pathways, such as peptidoglycan biosynthesis IV, biotin biosynthesis II, TCA cycle, L-arginine biosynthesis III and UDP-N-acetyl-D-glucosamine biosynthesis I (Figure 4C).

Abnormal gut metabolism is one of the common aging phenotypes, and the co-metabolites between gut microbes and the host are the bridge for their interactions. To clarify the effect of parishin treatment on aging-induced gut metabolic alterations, we quantified fecal metabolites in mice using gas chromatography–mass spectrometry (GC–MS). A total of 66 metabolites were identified from 40 fecal samples. The results of orthogonal projections to latent structures discriminant analysis (OPLS-DA) showed a clear separation between the YC, AC and LPar groups and between the YC, AC and HPar groups (Figure 5A), indicating significant differences in their metabolomic profiles. When using the values of variable importance for projection (VIP) > 1.5 as a threshold, the metabolites with significant contributions to differentiate the fecal metabolome profiles of the YC, AC and LPar groups were batyl alcohol, 3-(3-hydroxyphenyl)-propanoic acid, uracil and β-sitosterol (Figure 5B); the metabolites that contributed significantly to the fecal metabolome profiles of the YC, AC and HPar groups were ethanolamine, glycine, 2-α-mannobiose palmitelaidic acid, L-valine, sebacic acid, β-sitosterol and β-alanine (Figure 5C).

Next, we compared the fecal metabolite levels between each group. The levels of 3-(3-hydroxyphenyl)-propanoic acid, L-proline, 4-hydroxybenzeneacetic acid, 5-aminovaleric acid and β-sitosterol were significantly lower in the AC group than in the YC group. After treatment with high-dose parishin, the levels of these five metabolites were restored and no longer significantly different from those in the YC group; however, only the levels of 3-(3-hydroxyphenyl)-propanoic acid were recovered after low-dose parishin treatment (Figure 5D).

To explore the potential relationship between the gut microbiota and metabolome, Spearman’s rank correlation test was used to analyze the correlations between significantly altered bacterial taxa and metabolites in the feces. Both the absolute value of the correlation coefficient r > 0.4 and P < 0.05 were used as the screening threshold. The results showed that Turicibacter and Saccharimonas were negatively associated with 4-hydroxybenzeneacetic acid; Atopobiaceae and Prevotellaceae NK3B31 group were positively associated with β-sitosterol; Tannerellaceae and Parabacteroides were negatively correlated with β-sitosterol; Limosilactobacillus was negatively correlated with 5-aminovaleric acid; and Oscillospiraceae sp. was negatively correlated with L-proline (Figure 6A).

Next, we analyzed the association of gut microbiota and metabolites with aging-related indicators in mice. The results showed that GDF15 was positively correlated with Tannerellaceae, Saccharimonas, Turicibacter, Parabacteroides and Oscillospiraceae sp. and negatively correlated with 3-(3-hydroxyphenyl)-propanoic acid and 4-hydroxybenzeneacetic acid; p16^Ink4a was negatively correlated with the Prevotellaceae NK3B31 group and positively correlated with Defluviitaleaceae UCG-011; IL-6 was negatively correlated with Atopobiaceae; lung fibrosis degree was positively correlated with Tannerellaceae, Lachnospiraceae, Turicibacter, Parabacteroides, Oscillospiraceae sp. and Defluviitaleaceae UCG-011 and negatively correlated with Prevotellaceae NK3B31 group and Atopobiaceae; and the cardiac fibrosis degree was positively correlated with Tannerellaceae, Turicibacter, Rikenellaceae RC9 gut group and Parabacteroides and negatively correlated with Atopobiaceae (Figure 6B).