Metabolomics reveals potential biomarkers in the rumen fluid of dairy cows with different levels of milk production

Objective In the present study, an liquid chromatography/mass spectrometry (LC/MS) metabolomics approach was performed to investigate potential biomarkers of milk production in high- and low-milk-yield dairy cows and to establish correlations among rumen fluid metabolites. Methods Sixteen lactating dairy cows with similar parity and days in milk were divided into high-yield (HY) and low-yield (LY) groups based on milk yield. On day 21, rumen fluid metabolites were quantified applying LC/MS. Results The principal component analysis and orthogonal correction partial least squares discriminant analysis showed significantly separated clusters of the ruminal metabolite profiles of HY and LY groups. Compared with HY group, a total of 24 ruminal metabolites were significantly greater in LY group, such as 3-hydroxyanthranilic acid, carboxylic acids, carboxylic acid derivatives (L-isoleucine, L-valine, L-tyrosine, etc.), diazines (uracil, thymine, cytosine), and palmitic acid, while the concentrations of 30 metabolites were dramatically decreased in LY group compared to HY group, included gentisic acid, caprylic acid, and myristic acid. The metabolite enrichment analysis indicated that protein digestion and absorption, ABC transporters and unsaturated fatty acid biosynthesis were significantly different between the two groups. Correlation analysis between the ruminal microbiome and metabolites revealed that certain typical metabolites were exceedingly associated with definite ruminal bacteria; Firmicutes, Actinobacteria, and Synergistetes phyla were highly correlated with most metabolites. Conclusion These findings revealed that the ruminal metabolite profiles were significantly different between HY and LY groups, and these results may provide novel insights to evaluate biomarkers for a better feed digestion and may reveal the potential mechanism underlying the difference in milk yield in dairy cows.


INTRODUCTION
Nutritional requirements plays an important role in dairy cows' productivity, it is necessary to find an effectively digested feed and to maintain healthy physiological functions for efficient milk production. Interestingly, our previous study revealed that the compositions and structures of ruminal bacterial communities between high-yield (HY) and low-yield (LY) dairy cows were significantly different, which advanced the understanding underlying regulation response to the different milk yield [1]. Although the rumen microbial population has an explicit relationship with volatile fatty acids (VFAs), amino acids as well as fatty acids, ruminal metabolites can comprehensively reflect rumen digestive and systemic health [2]. However, the nutrient metabolism in the ruminal system between dairy cows having different ability in milk production remains unclear. Therefore, further research exploring the metabolic profile is necessary to elucidate the mechanism of metabolic regulation of dairy cows.
Metabolomics is an innovative global assessment and quantitatively measures the small endogenous metabolites in biological sample using high-throughput techniques, such as nuclear magnetic resonance [3], liquid chromatography/ mass spectrometry (LC/MS) [4], gas chromatography−mass spectrometry [5] and gas chromatography time-of-flight/mass spectrometry (GC-TOF/MS) [6]; moreover, in recent years, metabolomics approaches have been widely utilized for their high resolution and detection sensitivity for ruminal metabolite biomarkers. Most previous studies employed interesting metabolomics applications to reveal the impact of diet on animal health and performance [7], to investigate potential metabolite biomarkers of disease or physiology [4,8] and to examine diet-induced metabolic alterations in the ruminal microbial metabolic profile and their effects on milk production-related traits in ruminants [5,9,10]. Thus, the identification and integrative analysis of ruminal metabolites may facilitate the comprehensive characterization of metabolic mechanisms at the molecular and cellular levels in response to internal or external stimuli.
Therefore, it is interesting whether differences in the milk yields of cows were directly associated with their ruminal bacterial metabolomics. The objective of this study was to illustrate the metabolic mechanisms underlying the differences between high and low milk production, especially the rumen bacterial community, by examining the metabolite profiles of ruminal fluid in dairy cows.

Animals and experimental design
The experimental procedures were approved by the Animal Care Committee, Beijing University of Agriculture (Beijing, China). The experiment was designed as previously described [1]. Briefly, sixteen Holstein lactating dairy cows were divided into two groups according to the characteristics of parity (2.6 ±0.4), average body weight (670±24 kg), average dry mater intake (24.2±2.7 kg/d), number of lactation days (114.6±7.5), and these cows were assigned to either HY or LY group (31.90 ±1.76 kg/d and 19.30±1.76 kg/d, respectively; n = 8). The experiment was conducted over 21 days, with the first 14 days serving as an adaptation period. Cows were individually housed in tethered stalls in a barn with good ventilation and were fed and milked using a pipeline milking system three times per day at 0630, 1400, and 2000 h. The methods of samples collection and analysis were described in our previous study [1].

Liquid chromatography/mass spectrometry analysis
The rumen fluid samples were thawed at room temperature. Samples (100 μL) were used for metabolomics analysis. L-2-chlorophenylalanine (10 μL, 0.3 mg/mL), added, mixed and shaken for 10 s. Then, 300 μL of methanol-acetonitrile (2:1, v/v) added, mixed and shaken for 60 s. The solution was ultrasonically extracted on ice for 5 min, incubated at -20°C for 30 min and then centrifuged for 15 min at 3,850 g at 4°C. Two hundred microliters of supernatant was used for LC/MS analysis.
The LC/MS was carried out using an Ultimate 3000-Velos Pro system equipped with a binary solvent delivery manager and a sample manager, coupled with an LTQ Orbitrap mass spectrometer equipped with an electrospray interface (Thermo Fisher Scientific, Sunnyvale, CA, USA). The LC was an Acquity BEH C18 column (100 mm×2.1 mm i.d., 1.7 μm; Waters, Milford, MA, USA). Separation was achieved with solvent B (acetonitrile) and solvent A (aqueous 0.1% (v/v) formic acid) with the following gradient at a flow rate of 0.40 mL/min: 5% B-25% B over 0-1.5 min; 25% B-100% B over 1.5-10.0 min; 100% B-100% B over 10.0-13.0 min; 100% B-5% B over 13.0-13.5 min; and 13.5-14.5 min holding at 5% B. The injection volume was 3.0 μL, and the column temperature was set at 45.0°C. The mass spectrometric data were collected using an LTQ Orbitrap mass spectrometer equipped with an electrospray ionization source operating in either positive or negative ion mode. The capillary and source temperatures were set at 350°C, with a desolvation gas flow of 45 L/h. Centroid data were collected from 50 to 1,000 m/z with a 30,000 resolution.
The quality control (QC) sample was prepared by mixing aliquots of all samples into a pooled sample. Then, the QC sample was analysed using the same method as the experimental samples. The QC samples were injected at regular intervals (every 10 samples) throughout the analytical run to provide a set of data to assess repeatability.

Statistical analysis
The positive and negative data were combined to obtain a combined dataset that was imported into the SIMCA-P+ 14.0 software package (Umetrics, Umeå, Sweden). Principle component analysis (PCA) [11] and (orthogonal) partial least squares discriminant analysis ((O)PLS-DA) were carried out to visualize the metabolic alterations among experimental groups after mean centring and unit variance scaling. The variable importance in the projection (VIP) ranks the overall contribution of each variable to the (O)PLS-DA model, and those variables with VIP>1.0 are considered relevant for group discrimination. In this study, the default 7-round cross-validation was applied, and one-seventh of the samples were excluded from the mathematical model in each round to guard against overfitting. Significant differences in the metabolites of the two groups were analysed using the Wilcoxon rank-sum test. A heatmap of the key metabolites analyze in this process with methods similar to those previously published [12]. The correlation matrix between metabolites and ruminal bacterial www.ajas.info

81
Zhang et al (2020) Asian-Australas J Anim Sci 33:79-90 species was generated using Spearman correlation coefficient and visualized by using R language.

RESULTS
Ruminal pH, volatile fatty acids, dry matter intake, milk yield and milk composition Data describing the ruminal parameters and performance of the experimental cows have been reported in Tong's prior study [1]. Briefly, the ruminal pH was similar in both HY and LY groups. Compared with LY group, the proportion of propionate (p = 0.08) and total VFA concentration (p<0.05) in HY group were increased, whereas the acetate to propionate ratio (p = 0.06) and the proportion of acetate (p = 0.06) in HY group showed decreasing trends. In addition, the dry matter intake was higher in HY group than LY group (p = 0.03). Compared with LY group, milk production, 4% fatcorrected milk yield and energy-corrected milk yield were significantly increased in HY group. However, there was an increasing trend in milk fat content (p = 0.08) and milk protein content (p<0.01) in LY group compared to HY group. No differences were observed in milk lactose content and somatic cell count between LY and HY groups.

Differences in the ruminal metabolites of the two groups
In total, 10,544 practicable peaks that were unique and nonoverlapping were obtained from the rumen fluid samples. After strict QC and identification, 366 metabolites were dif-ferentially detected in LY and HY groups. These metabolites were mainly glycerophosphoethanolamines, fatty acids, amino acids, peptides, carbohydrates and pyrimidines.
To characterize the variations in the metabolic profiles of HY and LY dairy cows, PCA and OPLS-DA were conducted. As presented on the PCA score plots (Figure 1), HY and LY groups were clearly separated; PCA 1 and PCA 2 accounted for 25.8% and 15.9% of the total variation, respectively. Furthermore, in our study, the PCA of the LC/MS data showed that all the QC samples obtained in the positive and negative ion modes mostly overlapped in the same area, which indicated that this model was stable, reproducible and consistent for all the samples. As shown in Figure 2, not surprisingly, OPLS-DA model presented that the rumen fluids of HY and LY dairy cows had significantly different and distinct metabolite compositions. Both PCA and OPLS-DA indicated that metabolites found in the rumen fluid samples of HY and LY dairy cows were markedly distinct.
After the results of statistical analysis and the VIP value obtained from OPLS-DA were considered, 54 significantly different metabolites (p<0.05 and VIP>1) between the rumen fluids of HY and LY dairy cows were identified, as shown in Table 1. These identified metabolites were diverse and included the following: benzene and substituted derivatives; carboxylic acids and derivatives; coumarins and derivatives; diazines; fatty acyls; glycerophospholipids; indoles and derivatives; organooxygen compounds; pyrimidine nucleosides; steroids and steroid derivatives. As shown in Figure 3, each class of compound had multiple differentially abundant me-     VIP, variable importance in the projection. 1) FC represents the fold change, the ratio of the mean value of the peak area obtained from low-yield cows and the mean value of the peak area obtained from high-yield cows. FC > 1 indicates that this metabolite is more abundant in the low-yield cows than in the high-yield cows. www.ajas.info

83
Zhang et al (2020) Asian-Australas J Anim Sci 33:79-90 tabolites in the rumen fluids of HY and LY dairy cows. Compared with HY group, 24 ruminal metabolites were significantly increased in LY group, such as 3-hydroxyanthranilic acid, carboxylic acids, carboxylic acid derivatives (Lisoleucine, L-valine, L-tyrosine, L-lysine, dimethylglycine), diazines (uracil, thymine, cytosine), palmitic acid, and imidazopyrimidines. The remaining 30 metabolites were decreased in LY group, such as gentisic acid, caprylic acid, myristic acid, linoleic acid, and glycerophospholipids (p<0.05). Interestingly, remarkable alterations in the content of carboxylic acids, carboxylic acid derivatives, coumarins, coumarin derivatives, diazines, imidazopyridines, indoles, indole derivatives, organonitrogen compounds, phenols, and pyrroles in LY and HY groups were observed in the present study. Specifically, the VIP, variable importance in the projection. 1) FC represents the fold change, the ratio of the mean value of the peak area obtained from low-yield cows and the mean value of the peak area obtained from high-yield cows. FC > 1 indicates that this metabolite is more abundant in the low-yield cows than in the high-yield cows.
pyrrole metabolites and fluvastatin levels were increased by 6.06-fold (p<0.01) in LY group compared with HY group. As shown in Figure 4, hierarchical clustering analysis (HCA) shows significantly changed metabolite responses to different milk yields in dairy cows, which provides a further understanding of the metabolites that differ between HY and LY dairy cows. In the present study, HCA revealed that the significantly increased metabolites in LY group compared with HY group were mainly gathered into two subclusters. As shown in the lower part of Figure 4, compared with HY group, the metabolites either increased or decreased in LY group were separated clearly and significantly different subclusters were located. One subcluster consisted of the significantly upregulated metabolites in LY group compared with HY group, such as 3-hydroxyanthranilic acid, L-lysine, thymine and uracil. The other subcluster consisted of the remarkably upregulated metabolites from HY group, including glycerophospholipids, prenol lipids, and organooxygen compounds.

Metabolic pathway analyses
To further understand how multiple pathways alternated in response to HY and LY milk production, analysis of the functions of the pathways association with differential metabolites was subjected to the Kyoto encyclopedia of genes and genomes (KEGG). The dramatically impacted pathways identified with this enrichment analysis are shown in Table 2. ABC transporters, 2-oxocarboxylic acid metabolism, biosynthesis of amino acids, biosynthesis of unsaturated fatty acids and protein digestion and absorption were determined to be the main pathways that were different between HY and LY in the present study.
As determined by pathway topology analysis and as shown in Figure 5, 13 main metabolic pathways were enriched in our study, including the followings: ABC transporters; biosynthesis of unsaturated fatty acids; protein digestion and absorption; 2-oxocarboxylic acid metabolism; biosynthesis of amino acids; fatty acid biosynthesis; aminoacyl-tRNA biosynthesis; pyrimidine metabolism; cyanoacetic acid metabolism; glucose-inolate biosynthesis; sphingolipid metabolism; valine, leucine and isoleucine biosynthesis; and mineral absorption. Among these metabolic pathways, protein digestion and absorption, ABC transporters and biosynthesis of unsaturated fatty acids were the three most impacted pathways between the two groups, as determined by the richness factor values.

Correlation analysis between the ruminal microbiome and metabolites
For further understanding of the functional correlation between the perturbations in the ruminal microbiome and metabolites as a result of different levels of milk production, Spearman correlation coefficient analysis was conducted in the present study. Clear correlations were identified between the perturbed ruminal microbiome and altered metabolite profiles (r>0.5 or <-0.5, p<0.05). Figure 6 shows that several definite metabolites were greatly correlated with specific ruminal bacteria, which established the functional correlation between the microbiome and metabolites in rumen. It was found that Firmicutes, Actinobacteria, and Synergistetes phyla were highly correlated with most metabolites. Additionally, Bacteroidetes phylum was significantly correlated with only D-xylose, Saccharibacteria was significantly correlated with only apigenin, and Fibrobacteres was significantly correlated

DISCUSSION
In the present study, LC/MS metabolomics was used to evaluate the metabolites in the ruminal fluid of HY and LY dairy cows. We also wanted to investigate whether the observed changes in metabolites could provide further insight into specific ruminal microflora-related changes and aimed to provide a new perspective regarding perturbations of ruminal metabolism to reveal the mechanism of ruminal digestion in dairy cows with different milk yields.
Our results showed a clear separation of the ruminal metabolites in HY and LY dairy cows, indicating significant differences in metabolic composition. In the present study, L-isoleucine, L-valine, L-tyrosine, L-lysine, and dimethylglycine were significantly increased in LY group than HY group. Allison [13] previously reported that valine is generated from the main substrates pyruvate and isobutyrate. In addition, pyruvate can be converted to lactate by lactate dehydrogenase [14] or to acetyl coenzyme-A and formate by pyruvate formate-lyase [15], then acetyl-CoA is converted either to acetate or ethanol. Moreover, Andries et al [16] reported that branched chain fatty acids are also synthesized from deamination of the amino acids valine, isoleucine, leucine, and proline. Our previous study showed that acetate was remarkably increased in LY group compared with HY group, which further suggested increased substrate for amino acid synthesis. Furthermore, the enriched pathways observed in the present study, including aminoacyl-tRNA biosynthesis, amino acid biosynthesis, cyanoamino acid metabolism, glucosinolate biosynthesis, pantothenate and CoA biosynthesis as well as protein digestion and absorption, also confirmed the above hypothesis. Thus, we can reasonably speculate that high levels of acetate could promote the biosynthesis of amino acids in the rumen, and thereby this metabolite could reflect the condition of digestion in the rumen.
Moreover, compared with LY group, myristic acid, palmitaldehyde, linoleic acid, and alpha-linolenic acid were significantly increased in HY group. The concentrations of ruminal long-chain fatty acids are indicative of active lipolysis, biohydrogenation and microbial fatty acid synthesis in the rumen [17]. Additionally, before these unsaturated fatty acids can be rapidly hydrogenated by microbes into saturated Table 2. Metabolic pathways of the significantly enriched metabolites between high-yield and low-yield dairy cows end products, lipases, galactosidases, and phospholipases produced by ruminal microbes remove nonesterified fatty acids [18]. Furthermore, palmitic acid, linoleic acid and alphalinolenic acid were reported could be as sensitive and specific candidate biomarkers to distinguish low-and high-quality steers [19]. Therefore, it could be concluded that these fatty acids were detected at relatively high levels due to the microbial fermentation in the rumen of HY group; microbial fermentation and the efficient absorption of violate fatty acids may explain the mechanism of different milk yields. One remarkable alteration observed in this study were the increased levels of D-maltose and D-glucose in HY group compared with LY group. The nutrients in the diet could be degraded to maltose by amylase and then to glucose by malt- Figure 5. Pathways associated with the metabolites identified in both the high-and low-milk-yield dairy cows. The x-axis represents the pathway richness factor, and the y-axis represents the pathway name. Large sizes and dark colours represent a large number of metabolites and high pathway impact values, respectively. ase or maltose phosphorylase in the rumen. Previous research has also indicated that high levels of glucose promote the production of pyruvate in the rumen [20,21]. Glucose is the major monosaccharide liberated during the degradation of starch, and it can be converted into a number of different polyols and amino acids via the glycolytic pathway and its branches [22]. Microorganisms convert carbohydrates to pyruvate and acetyl-CoA by glycolytic pathway and pentose phosphate pathway [14]. Furthermore, previous studies demonstrated that carbohydrate metabolism caused significant increases in lactate and propionate concentration as well as decline in the acetate concentration [8,23], which was in accordance with our previous results that acetate levels were lower and propionate levels were greater in HY dairy cows relative to the LY dairy cows. It has also been reported that diets rich in readily available carbohydrates are associated with alterations in the rumen microbiota, which are followed by significant changes in the metabolic pathways of the rumen [8,24]. Interestingly, our correlation results revealed that Dglucose was significantly positively correlated with Firmicutes and that D-maltose was significantly positively correlated with Actinobacteria. Based on our integrated pathway analysis, it is reasonable to conclude that D-maltose and D-glucose metabolites were increased in HY group compared with LY group.
One of the most interesting observations from this study was that phyla of the ruminal microbiome were specifically correlated with the different metabolites detected in HY and LY dairy cows. It is well known that ruminal metabolites do not independently exist, and a comprehensive understanding of the variation in other metabolites is warranted. In the present study, Firmicutes was mostly correlated with fatty acyls. Furthermore, it has been reported that Firmicutes plays an important role in milk production [25]. Thus, it could be suggested that the fatty acyl metabolites contributed to the increased abundance of these bacteria in the rumen fluid. Therefore, these correlations of ruminal metabolites with the ruminal microbiome can reflect metabolic processes. Additionally, Bacteroidetes was significantly negatively correlated with only D-xylose, and D-xylose level was significantly increased in HY group compared with LY group. Accumulated evidence strongly suggests that the metabolic alterations associated with microbiome perturbations are important biomarkers that indicate the health [19], nutrition [26], or dietary changes in response to the physiology [27] of dairy cows. Taken together, these results indicate that rumen microbes promote protein degradation in HY dairy cows, which provides a better understanding of the difference in the milk proteins of the two groups. Therefore, correlations between ruminal metabolites and ruminal microbiome in the present study may provide new information for advanced understanding of complex rumen metabolism.

CONCLUSION
These data indicated that not only the metabolites in the ruminal fluids were significantly different between HY and LY cows, but also these metabolites associated with the ruminal microbiota affected the metabolic function. These highly correlated metabolites may be potential biomarkers for the measurement of digestive and rumen function. In this study, the integration of the high-throughput sequencing of the ruminal microbiome and metabolomics should provide new insight to improve the understanding of the physiological and metabolic mechanisms involved in milk production. Moreover, these findings suggest that due to host-microbiota interactions, dairy cows with different milk yields have different metabolite content, which are most likely correspond to definite differences between ruminal fermentation and milk component.

CONFLICT OF INTEREST
We certify that there is no conflict of interest with any financial organization regarding the material discussed in the manuscript. Tong J is an employee of Beijing Bei Nong Enterprise Management Co., Ltd..