The PM hybrid (Broussonetia kazinoki × B. papyrifera) used in this study was cultivated at the mountain experimental field of Guizhou Qianchang Shenghe modern agriculture Co., Ltd. (106^o37′E longitude and 25^o73′N latitude, Daihua, Guizhou, China). The station had a mid-subtropical monsoon humid climate with a mean annual temperature of 15.5°C (extreme maximum temperature is 34°C–36°C, minimum temperature is −6°C to −9°C), annual relative humidity of 81%, and annual rainfall of 1,250–1,400 mm. The growth height of PM is about 1.8 m, and it is harvested 4–5 times a year, leaving 20 cm of stubble each time. Fresh young shoots and leaves were collected from the second cutting, and three field replicates were performed on 12 June 2020.

The CB, RB, and WB materials were purchased from Guilian Animal Feed Co., Ltd. (Guiyang, China). First use a chopper (92-2S, Sida Agri-Machine Co., Ltd., Luoyang, China) to cut PM into a length of 1–2 cm, and then silages were prepared from PM unmixed. The mixture silages were also prepared from PM and brans (CB, RB, or WB) at a percentage of 70 + 30, based on the fresh matter (FM). Approximately 20 kg of PM and mixture materials for each treatment were packed into the polyethylene drum silos (30 L, Huafang Plastic Co., Ltd., Jiangsu, China) with three replicates, air in the material was tightly compressed and discharged, and then sealed with a cover to maintain anaerobic conditions. The silos were stored in a ventilated warehouse at ambient temperature (21°C–35°C). After 128 days of fermentation, the silages were opened, and the samples of each treatment were divided into three parts for microbial and chemical analysis. The fresh parts were immediately collected in a sterilization bag and kept in an icebox, and then immediately transported to the laboratory for analysis of microbial population and chemical composition. At the same time, samples for analyzing microbial community and metabolites are quickly stored in a freezer at −80°C for future extraction of DNA and liquid.

The dry matter (DM), organic matter (OM), crude protein (CP), ether extract (EE), neutral detergent fiber (NDF), acid detergent fiber (ADF), and acid detergent lignin (ADL) were analyzed in accordance with the methods of the Association of Official Analytical Chemists (AOAC International, 2000). WSC, pH, ammonia nitrogen (NH₃-N) and organic acids were determined in accordance with methods as described by Cai (2004). The pH was determined by a pH meter (D-71, Horiba, Kyoto, Japan). The NH₃-N content was measured using the Kjeltech auto distillation system (2200, Foss Tecator, Hoganas, Sweden). The organic acid was analyzed by a high-performance liquid chromatography system (HPLC, LC-2000, Jasco, Tokyo, Japan) by following conditions: detector, UV-2070 (JASCO Corporation, Tokyo, Japan), 450 nm; column, Shodex RS pak KC-811 (8.0 mm 30 cm; Shoko, Tokyo, Japan); oven temperature, 60°C; eluent, 3 mmol/L HClO₄; flow rate, 1.0 ml/min. Lactic acid buffer capacity (LBC) was analyzed by a titration method (Muck et al., 1991).

The microbial population was determined by a plate count method as described by Cai (2004). The LAB, aerobic bacteria, coliform bacteria and yeasts (including molds) were analyzed by the mediums of de Man, Rogosa, and Sharpe agar (Difco Laboratories, Detroit, MI, United States), nutrient agar (Nissui-Seiyaku Co., Ltd., Tokyo, Japan), blue light broth agar (Nissui-Seiyaku Co., Ltd.) and potato dextrose agar (Nissui-Seiyaku Co., Ltd.), respectively. LAB were incubated in an anaerobic box (ANX-1; Hirosawa Ltd., Tokyo, Japan) at 30°C for 2–3 days, and other microorganisms were in an anaerobic incubator at 30°C for 2–5 days. The growth status of microbial colonies were observed every day during the incubation period, and the colonies were counted as viable numbers in colony-forming unit per gram (cfu/g) of FM.

For SMRT sequencing analysis, the triplicate samples (10 g) were mixed with 90 ml of sterilized 0.85% NaCl solution and shaken at the speed of 250 rpm in a 4°C refrigerator for 45 min. The liquid mixture was filtered through a four-layer cheesecloth pre-autoclaved, and then the filtrate was centrifuged at 10,000 rpm for 10 min at 4°C to obtain microbial precipitate (Du et al., 2021b). Total genomic DNA from each sample was extracted using the TIANamp Bacteria DNA isolation kit (DP302-02, Tiangen, Beijing, China). The quality of the extracted DNA was assessed in 1% agarose gels using electrophoresis and spectrophotometry (optical density at 260/280 nm ratio). All DNA samples were purified through a DNA kit column and stored at −20°C until further analyses (Hou et al., 2015; Du et al., 2021b).

The PCR product of the full length of 16S rRNA regions were amplificated as described by Du et al. (2021b). The 16S rRNA gene was amplified by PCR for SMRT sequencing using the 27F (AGRGTTTGATYNTGGCTCAG) and 1492R (TASGGHTACCTTGTTASGACTT) primers. The PCR program was: 95°C for 5 min; 30 cycles of 95°C for 30 s, 50°C for 30 s, and 72°C for 60 s, with a final extension of 72°C for 7 min. The Agencourt AMPure XP Beads (Beckman Coulter, Indianapolis, IN) were used for purification of the total of PCR amplicons. The Qubit dsDNA HS Assay Kit and Qubit 3.0 Fluorometer (Invitrogen, Thermo Fisher Scientific, Oregon, United States) was used for quantification analysis. Library construction used PacBio SMRTbell™ Template Prep Kit 2.0 (Menlo Park, CA, United States) on normalized pooled PCR products (Earl et al., 2018). The amplicons were sequenced on a PacBio RS II platform (Pacific Biosciences) using P6-C4 chemistry (Hou et al., 2015). The quality control of PCR amplifications and sequence pre-processing were performed according to the method described by Mosher et al. (2013).

Raw data were processed using the protocol RS_Readsofinsert.1 available in SMRT Portal version 2.7 software (PacBio). Firstly, the high-quality sequences were performed using the Quantitative Insights In to Microbial Ecology (QIIME) package (version 1.7.0; Caporaso et al., 2010). Under 100% clustering of sequence identity, using the analyses of Python nearest alignment space termination (PyNAST) and clustering and classification inference with U-statistics (UCLUST) to align the extracted high-quality sequences to obtain representative sequences (Edgar, 2010). After the selection of the representative sequences, the unique sequences were classified into operational taxonomic unit (OTU) based on 99% threshold identity using the UCLUST algorithm (Lozupone and Knight, 2005). Potential chimeric sequences in the representative set of OTU were removed using the Chimera Slayer tool (Haas et al., 2011). According to the analysis of the OTU clustering, and the research requirements to analyze both the core and unique information for material and silage samples. The SILVA database version 132 was implemented to classify different OTU and annotate the taxonomic information for each OTU representative sequence based on Bergey’s taxonomy at the genus, family, order, class, and phylum levels, according to classification at an 80% minimum bootstrap threshold (Quast et al., 2013). OTU that occurred only once or twice were discarded. Unweighted Pair-group Method with Arithmetic Means (UPGMA) clustering was performed to use average linkage to explain the distance matrix in QIIME software (version 1.7.0; Chen et al., 2020).

The silage sample (5 g) was placed in 20 ml of cold extraction liquid (methanol:water = 4:1) for extraction in a microcentrifuge tube. Each treatment performs six biological replicates for metabolomic profiling detection analysis. The untargeted metabolomics analyses were carried out by ultra HPLC tandem time-of-flight mass spectrometry (UHPLC/TOF-MS). The liquid chromatography with mass spectrometry (LC–MS) untargeted metabolomics analyses were performed on an Agilent 1100 series HPLC system (Agilent Technologies, Palo Alto, CA, United States) with a UPLC BEH Amide column (1.7 μm, 2.1 mm × 100 mm, Waters, Milford, United States) coupled to Triple TOF 5600 (Q-TOF, AB Sciex, Framingham, Massachusetts, United States). The detailed information of the UHPLC/TOF-MS program was following the previous study by Wang et al. (2019). In order to evaluate the stability of the analytical system during the determination process, every five samples were designed as a quality control (QC) sample.

Before statistical analysis, the files of LC–MS raw data were converted for pre-processing by the software of Progenesis QI (version 2.0, Waters Corporation, Milford, MA, United States). The result of pre-processing generates a data matrix with the following contents: baseline filtration, peak identification, integration, retention time correction, peak alignment, and finally a data matrix of retention time, mass-to-charge ratio, and peak intensity. In all the samples, variables with non-zero values more than 50% were retained; the missing values with one-half of the minimum value in the original matrix were filled; the total peaks were normalized, and the QC samples with relative standard deviation (RSD) ≥ 30% were deleted. The LC-MS analytical system had appropriate stability and repeatability and that the acquired data were sufficient quality for subsequent determination.

The metabolite annotation was performed by software Progenesis QI (version 2.0, Waters Corporation, Milford, MA, United States), integrated with commercial databases of human metabolome database (HMDB). Based on the detected peaks, the metabolites were identified by the interquartile range denoising method. Missing raw data values were set to half of the minimum value of the detection limit. In addition, an internal standard normalization method was used in these data analyses.

According to metabolite clustering, the Venn diagram package (version 1.2) of R packages statistical tools was draw (Shade and Handselman, 2012). The value of variable importance in projection (VIP) score was analyzed by ropls (R packages, version 1.6.2) and scipy (Python, version 1.0.0; Chong et al., 2018). The principal component analysis (PCA), partial least-square discriminant analysis (PLS-DA), and orthogonal partial least squares discriminant analysis (OPLS-DA) were analyzed by the ropls (R packages, version 1.6.2) software package (Wehrens and Mevik, 2007). Input data were the total mass of the signal integration area of each sample, and the signal integration area was normalized with method of internal standard normalization for each sample. Correlation analyses was also performed using a Python-based tool (Langfelder and Horvath, 2008).

The ANOVA was performed using the general linear model (GLM) procedure of Statistical Package for the Social Sciences (SPSS version 19.0, SPSS Inc., Chicago, IL, United States) to examine the differences between samples, and significance was declared at p < 0.05. The LBC, chemical composition and fermentation quality data for the materials and silages were subjected to one-way ANOVA. Tukey’s honest significant difference test (HSD) was employed for different sample means (Steel and Torrie, 1980).

Welch’s two-sample T-test was used to identify biochemicals that differed significantly between experimental groups. We defined significantly different compounds between treatment groups by the criteria of VIP > 1 and p < 0.05. The correlation analyses of the bacterial community with metabolites with biofunctional activity at species level, respectively. The metabolites are displayed horizontally, respectively, and the bacterial community information is displayed vertically. The corresponding value of the middle heat map is the Spearman correlation coefficient r, which ranges between −1 and 1, r < 0 indicates a negative correlation (blue), r > 0 indicates a positive correlation (red), and “*,” “**,” and “***” represent p < 0.05, p < 0.01, and p < 0.001, respectively.

In the tropics, with the increase in population and decrease in arable land, the competition between humans and animals for food and feed has become increasingly prominent. Conventional feeds such as forage crop, grass, and grain cannot meet the needs of ruminants. In response to the shortage of feed in the dry season, natural feed resources such as WP are under active development for use in animal production. PM is a popular WP in tropical regions and has good livestock palatability. Its nutritional value is similar to the nutritional value of legume grasses, and it has the reputation of woody alfalfa (Dong et al., 2020). In this study, the CP content of fresh stems and leaves of PM was 24% of DM, indicating that it can be used as a source of high-protein livestock feed. Because of the high moisture of WP, like other herbaceous forages, it is generally difficult to make high-quality silage (Cai et al., 1999; Du et al., 2021b). As shown in Table 1, the moisture content of PM was >80%, indicating the need for adjustment during the silage process. In addition, WSC, LAB, and LBC affect silage fermentation (Cai et al., 1999). The addition of 30% CB, RB, and WB to PM reduced the moisture content to <65%, within the ideal range (60%–70%) for preparing high-quality silage (McDonald et al., 1991). A suitable moisture content inhibits fermentation by harmful clostridia (Wan et al., 2021), while the addition of exogenous additives increases the fermentation substrate for LAB and alleviates LBC, thereby improving silage fermentation quality (Cai et al., 2020).

In this study, the addition of CB and WB improved silage fermentation quality (Table 2). However, contrary to our expectations, the addition of RB resulted in butyric acid fermentation and proteolysis in PM silage, and the underlying mechanism of this process is unclear. RB is rich in EE and can promote the proliferation of clostridia during ensiling (Kiyoshi and Choseki, 1975), and lead to protein degradation. Clostridia are commonly found in anaerobic environments such as paddy soils, and their spores are incorporated into bran during rice harvesting and rice production. Generally, the anaerobic environment formed by irrigation in paddy fields during rice cultivation is much higher than that in upland field, such as corn and wheat cultivation (Du et al., 2021b). Therefore, we speculate that absolute anaerobic microorganisms such as Clostridium species in paddy fields are higher than those in ordinary upland field. In this study, the results of SMRT sequencing confirmed the presence of Clostridium spp. in RB, whose silage contained a large amount of NH₃-N. Therefore, RB will affect the fermentation quality of PM silage, resulting in nutrient loss.

Various microorganisms are present in the stems and leaves of forage crops, grasses, and other plants. Their community structures, species diversities, and metabolic functions are important determinants of silage fermentation (McEniry et al., 2010; Guan et al., 2018; Keeling et al., 2019). In this study (Figure 1A), the genera Pantoea, Enterobacter, and Weissella were the major epiphytic bacteria in PM material. Pantoea is a Gram-negative pathogen, present in various environments; however, it is uncommon in silage. Pantoea decompose and metabolize organic compounds, leading to the deterioration of fermented food and feed (Avila and Carvalho, 2020). Enterobacter and Weissella are common microorganisms in silage. Enterobacter is a Gram-negative genus that dominates the bacterial community of fresh forage. In the early stage of silage fermentation, Enterobacter and Weissella compete with LAB, causing nutrient loss. However, as fermentation progresses, their proliferation is inhibited; they are subsequently replaced by LAB (Heron et al., 1993; Pahlow et al., 2003). In the silage process, D-lactic acid and carbon dioxide are produced, influencing the proportion of lactic acid isomers and increasing the loss of DM in silage (Cai et al., 1998). Oryza sativa Japonica Group Japanese rice was significantly enriched in the CB and WB materials; Methylobacterium, Sphingomonas, and Clostridium were predominant in the RB materials. Except for Clostridium, these bacteria are rare in silage materials, and their influence on silage fermentation is unclear. They are generally considered to be unfavorable bacteria for silage production (Du et al., 2021b). Clostridium affects the fermentation of silage by degrading proteins and causing nutrient loss.

After ensiling, Enterobacter was predominant in the PM and PR silages; it competed with LAB for nutrients in the anaerobic environment (Santos et al., 2016). The genus Clostridium was also a major taxon in PR silage. The addition of RB may increase the EE and clostridial spore contents, thus promoting the proliferation of Clostridium spp. (Queiroz et al., 2018). The effects of adding CB and RB to PM silage on the moisture and fermentation substrate contents may promote the growth of LAB and inhibit the growth of harmful bacteria. Indeed, the dominant genus in both PC and PW silages was Lactobacillus. Notably, the epiphytic microbial community in the materials comprised <1% of genera, indicating that these bacteria were rapidly suppressed after ensiling.

The dominant species in the materials were P. agglomerans, E. asburiae, and W. paramesenteroides in PM; they were Methylobacterium aquaticum and Sphingomonas endophytica in RB (Figure 1B). Except for W. paramesenteroides, these bacteria are considered harmful and presumed to reduce silage fermentation quality. W. paramesenteroides is a heterofermentative LAB typically found in silage. In the initial stage of fermentation, it grows vigorously and produces carbon dioxide, creating an anaerobic environment that promotes fermentation by other acid-resistant LAB. L. plantarum was predominant in high-quality PC and PW silage, whereas E. asburiae and C. tyrobutyricum were predominant in low-quality PM and PR silage. Such changes in microbial population structure indicate that these microbes are determinants of silage fermentation. The preparation of silage affects fermentation quality, and the microbial community in silage is important for the physiology and health of livestock.

The combination of microbiome and metabolomics profiling can provide microbial community and metabolite information; it can also expand the functions of silage (Xu et al., 2019). The materials before ensiling were aerobic but not acidic, and the high moisture content of fresh PM promoted the growth of microorganisms with corresponding metabolite production. However, the low moisture content of CB, RB, and WB suppressed this growth and metabolite production (Figure 2A).

Microorganisms are typically sensitive to anaerobic and acidic environments (Du et al., 2021b). Some Gram-negative aerobic bacteria have thin cell walls and weak acid resistance; they are rapidly inactivated by the anaerobic and acidic silage fermentation environment. However, LAB grow rapidly and metabolize both amino acids and carbohydrates in anaerobic and acidic environments. The increased LAB population in PC and PW silages increased their metabolite contents (Figure 2B). In the PM and PR silages, Enterobacter and Clostridium spp. were more abundant than LAB, reducing the metabolic capacity of LAB and their corresponding metabolite production.

Several functional organic acids (Figure 3A) and derivatives (carboxylic acids and derivatives) and flavoring agents (benzoic acids and derivatives; Figure 3B) were detected in the PM and mixture silages. These metabolites were mainly lipids and lipid-like molecules, organic acids and derivatives, and organic oxygen compounds (Table 3). Lipids and lipid-like molecules comprise a superclass of organic compounds, including fatty acids and their derivatives. LAB produce multiple metabolites of fatty acids during fermentation (Sun et al., 2012). Fatty acids improve the growth and milk fat components of dairy cattle, thereby improving their health. The greater proportion of LAB in PC and PW silages led to higher fatty acid content than in other silages. In the PM and PR silages, the harmful bacteria inhibited the growth of LAB, suppressing fatty acid synthesis. Organic acids (e.g., fumaric acid, malic acid, citric acid, and succinic acid) in silage can increase livestock performance and feed efficiency; they are used as substitutes for feed antibiotics (Ke et al., 2017). In this study, the greater proportion of LAB in PM and PR silages increased production of organic acids and derivatives, thereby improving silage fermentation quality. Organic oxygen compounds comprise a superclass of organooxygen compounds, including carbohydrates and carbohydrate conjugates. Compared with the materials, the organic oxygen compound content was significantly reduced in the silages. This is because during the ensiling process, Lactobacillus and Clostridium use glycosyl compounds to produce lactic acid and spores, thereby decreasing the carbohydrate content.

Microbial fermentation products including organic acids and 1,2-propanediol have been used as indicators of silage fermentation quality and aerobic stability (Xu et al., 2019; Cai et al., 2020). However, these indicators cannot fully reflect the complex metabolic processes of silage. Metabolomics analysis enables evaluation of the fermentation, nutritional, and functional characteristics of silage (Guo et al., 2018). The quality of feed protein depends on the balance among multiple amino acids (Tedeschi et al., 2015). In this study (Figures 4A–C), greater abundances of Leu, Try, Met, Phe, Ser, His, Val, Asp, and Ile were associated with better ensiling quality in PC and PW silages, compared with PM silage. The increased amino acid contents in PC and PW silages may be related to the unique sugar fermentation and amino acid biosynthesis pathways of LAB, such as L. plantarum. These data indicate that the high-quality PM silage generated by addition of CB and WB suppresses protein degradation. The opposite pattern was observed in PR silage (Figure 4D). Clostridium spp. in PR silage could promote the degradation of protein to NH₃-N. Feeding silage can provide fatty acids to ruminants, which can increase the content of unsaturated fatty acids in milk and reduce the ratio of saturated fatty acids; it also increases the levels of essential fatty acids. Therefore, the fatty acids produced by LAB in silage can be used to promote animal health and welfare (Salsinha et al., 2018). The CB and WB had significantly increased levels of unsaturated fatty acids (linoleic acid, isopropyl linoleate, indole carboxylic acid, and nitro-octadecenoic acid), which exert an antioxidant effect, eliminate free radicals, and reduce cancer incidence. However, RB inhibits the synthesis of unsaturated fatty acids. Some of the aromatic compounds detected are sweeteners or flavoring agents. Lactobacillus plantarum in PC and PW silages uses sugars to produce aromatic compounds (allyl butyrate, α-carboxy-delta-decalactone, sclareol, and ethylsuberic acid), improving feed flavor and quality. However, Clostridium spp. in PR silage inhibited aromatic compound biosynthesis, explaining its poor aroma. In addition, antimicrobial metabolites were detected. Benzopyran and ferulic acid exhibit antioxidant properties; they promoted the proliferation of LAB and inhibited harmful bacteria in PC and PW silages. However, acetylononin is an inhibitor of isoflavonoid o-glycosides, and it may inhibit organic compound synthesis.

In PC and PW silages, the addition of CB and WB caused Lactobacillus plantarum to predominate, thereby altering metabolite levels. Compared with PM silage, the addition of RB did not affect the bacterial community, likely because Enterobacter and Clostridium were predominant and led to similar metabolite profiles (Figures 5A,B). The results of permutation testing supported this finding (Figure 5C).

Spearman correlation analysis regarding the main metabolites of bacteria in silage showed that during fermentation, metabolites were positively correlated with LAB and negatively correlated with other bacteria (Xu et al., 2019). Citric acid and L-malic acid are key intermediates in the citric acid cycle; they were produced by most LAB species LAB can use D-(+)-cellobiose as a substrate for energy production, accelerating their growth and inhibiting the growth of yeasts, molds, and Clostridium and Enterobacter spp. (Ke et al., 2017). L-Malic acid and citric acid may suppress protein hydrolysis by reducing pH, which enhances mixed ruminal microorganism fermentation in vitro and improves animal performance (Wang et al., 2009). In addition, L-malic acid and citric acid are antioxidants and inhibit lipoxygenases. In this study (Figure 6), L. plantarum and W. paramesenteroides were positively correlated with the L-malic acid, citric acid, and D-(+)-cellobiose contents. In contrast, C. tyrobutyricum and E. asburiae were negatively correlated with these metabolites. Overall, our findings indicate that the silage microbiota and metabolome are closely related and exhibit important interactions. Succession of the microbial population affects metabolite levels; changes in metabolite levels also influence the microbial population. The metabolic activities of microorganisms affect silage fermentation quality and have important roles in the production of healthy livestock.