High‐production dairy cattle exhibit different rumen and fecal bacterial community and rumen metabolite profile than low‐production cattle

Abstract Our aim was to simultaneously investigate the gut bacteria typical characteristic and conduct rumen metabolites profiling of high production dairy cows when compared to low‐production dairy cows. The bacterial differences in rumen fluid and feces were identified by 16S rDNA gene sequencing. The metabolite differences were identified by metabolomics profiling with liquid chromatography mass spectrometry (LC‐MS). The results indicated that the high‐production dairy cows presented a lower rumen bacterial richness and species evenness when compared to low‐production dairy cows. At the phylum level, the high‐production cows increased the abundance of Proteobacteria and decreased the abundance of Bacteroidetes, SR1, Verrucomicrobia, Euryarchaeota, Planctomycetes, Synergistetes, and Chloroflexi significantly (p < 0.05). At the genus level, the rumen fluid of the high‐production group was significantly enriched for Butyrivibrio, Lachnospira, and Dialister (p < 0.05). Meanwhile, rumen fluid of high‐production group was depleted for Prevotella, Succiniclasticum, Ruminococcu, Coprococcus, YRC22, CF231, 02d06, Anaeroplasma, Selenomonas, and Ruminobacter significantly (p < 0.05). A total of 92 discriminant metabolites were identified between high‐production cows and low‐production cows. Compared to rumen fluid of low‐production dairy cows, 10 differential metabolites were found up‐regulated in rumen fluid of high‐production dairy cows, including 6alpha‐Fluoropregn‐4‐ene‐3,20‐dione, 3‐Octaprenyl‐4‐hydroxybenzoate, disopyramide, compound III(S), 1,2‐Dimyristyl‐sn‐glycerol, 7,10,13,16‐Docosatetraenoic acid, ferrous lactate, 6‐Deoxyerythronolide B, vitamin D2, L‐Olivosyl‐oleandolide. The remaining differential metabolites were found down‐regulated obviously in high‐production cows. Metabolic pathway analyses indicated that most increased abundances of rumen fluid metabolites of high‐yield cows were related to metabolic pathways involving biosynthesis of unsaturated fatty acids, steroid biosynthesis, ubiquinone and other terpenoid‐quinone biosynthesis. Most down‐regulated metabolic pathways were relevant to nucleotide metabolism, energy metabolism, lipid metabolism and biosynthesis of some antibiotics.


| INTRODUC TI ON
Many studies showed that the gut microbiota play an important role on health, metabolism and immunity of the host recently (Amato et al., 2014;Ridaura et al., 2013;Trompette et al., 2014). The bovine rumen harbors a symbiotic microbiota capable of converting indigestible plant mass into energy (Flint, Bayer, Rincon, Lamed, & White, 2008), which is of vital importance for production of milk and meat. Among this complex microbial community, 95% of the microorganisms are bacteria (Brulc et al., 2009). Jami, White, and Mizrahi (2014); recently identified the connections between the ratio of Firmicutes to Bacteroides and daily milk-fat yield. Variation in the rumen microbiome of dairy cattle has also been correlated with methane emission levels (Ross, Moate, Marett, Cocks, & Hayes, 2013b), and metagenomic profiling of the rumen microbiome can actually be used to predict phenotypes relating to enteric methane gas production (Ross, Moate, Marett, Cocks, & Hayes, 2013a). These studies revealed connections between the gut microbial community and certain physiological host parameters, which can be applied on improving animal production by manipulation of relevant beneficial bacterial flora.
Some studies suggested correlations between rumen microbial groups and bovine feed efficiency traits (Carberry, Kenny, Han, McCabe, & Waters, 2012;Guan, Nkrumah, Basarab, & Moore, 2008;Hernandez-Sanabria et al., 2012). The relative proportion of sequences that were assigned to Prevotella appeared to be positively associated with high residual feed intake (RFI) in bulls, whereas an unidentified group within the order Bacteroidales was positively associated with low RFI in bulls (McCann et al., 2014). Lima et al. (2015) have characterized the rumen fluid microbiomes of prepartum and postpartum high-producing Holstein cows and revealed that some bacteria have strong correlations with milk production. Moreover, they built a multivariable regression model using bacterial taxa significantly associated with average milk yield in the first 150 days postpartum to predict the weekly milk production; this microbiomepredicted milk yield was significantly correlated with the actual weekly averages of milk production.
Additionally, an early study on mice observed that obese mice (ob/ob) exhibited a different ratio of the phyla Firmicutes to Bacteroidetes when compared with lean littermates (Ley et al., 2005). This difference is not totally attributable to differences in food consumption, as a runted ob/ob mouse weighed less than the ob/ob littermates owing to reduced chow consumption, but still demonstrated a significantly greater percent body fat and ratio of Firmicutes to Bacteroidetes. Moreover, analogous differences have been observed in the distal gut microbiota of obese versus lean humans as well; the relative abundance of Bacteroidetes increases as obese individuals lose weight on either a fat-or a carbohydraterestricted low calorie diet. And the increase in Bacteroidetes was markedly relevant to weight loss but not to total caloric intake. (Ley, Turnbaugh, Klein, & Gordon, 2006). A later study by the same team demonstrated that the obesity-associated gut microbiome have increased capacity for energy harvest by transplantation of lean and obese cecal microbiotas into germ-free wild-type mouse recipients  However, there is still a lack of information investigating the gut bacterial profile of high-production dairy cows. In this work, we explored the bacterial differences in rumen fluid and feces from high-yield and low-yield dairy cows under the same diet, region and surroundings. Meanwhile, we investigated the metabolites differences in rumen fluid between these two groups as well. The aim of the study was to explore the rumen and fecal bacteria, and the rumen metabolite profiles of high-production dairy cows. Bacteria which are associated with high-production dairy cattle have the potential to be cultured and applied as probiotics to improve the performance of low-production dairy cattle.

| Experimental animals and sample collection
We selected 22 Holstein dairy cows, with 11 high-yield and 11 low-yield cows; 11 of each were matched as pairs, and each pair of cows was reared under the same diet regimes, feeding environment, and were paired for similar age, parity and nearing lactation days (Table 1) on the aote dairy farm in Qingdao. Cows were fed twice (05:00 hr and 17:00 hr) and milked twice daily; all cows had free access to clean water.
Representative rumen fluid samples were obtained from all Holstein cows via the cow's mouth with the oro-ruminal sampling device within two hours before morning feeding. Fecal samples were acquired by the rectum pick dung method. All samples were immediately placed into liquid nitrogen, and were transferred to laboratory for −80°C storage.
The ethics committee of Shandong agriculture university approved the study (SDAU2015-18).

| Experimental procedures for 16S rDNA sequencing
The samples were slowly thawed at 4°C. Total DNA was extracted from the rumen fluid and fecal samples using the Stool DNA Isolation Kit (Tiangen, Beijing, China). DNA samples were quantified, using a Nanodrop spectrophotometer (Nyxor Biotech, Paris, France), and then transferred to BGI Genomics for V4 region of the 16S rDNA K E Y W O R D S bacterial diversity, high-production dairy cows, low-production dairy cows, metabolites, rumen fluid gene sequencing with PE250 Miseq. The PCR primer used for 16S rDNA amplicon libraries was 515F-806R.

| Bioinformatics analysis for 16S rDNA sequencing
The raw data were filtered to eliminate the adapter pollution and low quality to obtain clean reads (Douglas et al., 2014). Truncation of sequence reads not having an average quality of 20 over a 30 bp sliding window based on the Phred algorithm, and trimmed reads having less than 75% of their original length, as well as its paired read, were removed. Then paired-end reads with overlap were merged into tags by FLASH (Magoc & Salzberg, 2011) (Fast Length Adjustment of Short reads, v1.2.11). Tags were clustered to OTU at 97% sequence similarity by scripts of software USEARCH (v7.0.1090) (Edgar, 2013). Chimeras were filtered out by using UCHIME (v4.2.40) (Edgar, Haas, Clemente, Q uince, & Knight, 2011). OTU representative sequences were taxonomically classi- Wilcoxon Rank-Sum Test was used for comparison of two groups, using the alpha diversity indices. Beta diversity analysis was done by software QIIME(v1.80) (Caporaso et al., 2010). PCoA (Principal coordinate analysis) is used to exhibit the differences between the samples according to the matrix of beta diversity distance.
Unweighted Pair Group Method with Arithmetic mean (UPGMA) is a type of hierarchical clustering method using average linkage and was used to interpreting the distance matrix produced by beta diversity. To measure the robustness of this result to sequencing effort, we perform a jackknifing analysis, wherein 75% of the smallest sample sequences from each sample are chosen at random, and the F I G U R E 1 PCA analyses of bacteria in rumen fluid and fecal samples. Note: Unweighted UniFrac was used to create the PCA. X-axis, 1st principal component; Y-axis, 2nd principal component. The number in brackets represents the contributions of principal components to differences among samples. A dot represents each sample, and different colors represent different groups. HighS and LowS represent rumen fluid from high-yield and low-yield cows, respectively. HighF and LowF represent groups of high-yield and low-yield cow feces, respectively resulting UPGMA tree from this subset of data is compared with the tree representing the entire available data set by QIIME(v1.80). This process is repeated with 100 random subsets of data, and the tree nodes which prove more consistent across jackknifed datasets are deemed more robust. And the figure is drawn by software R(v3.0.3).

TA B L E 1 Fundamental information of cows in experiment
The abundance differences in microbial communities between samples were obtained using statistical methods, and FDR (false discovery rate) was adopted to assess the significance of differences.

| Metabolites extraction
Metabolites extraction is achieved by organic solvent to precipitate protein (Sarafian et al., 2014). A quantity of 100 μl liquid sample with micropipette was plunged into 96-well plate. 300 μl isopropanol (−20°C) was added to the well above and vortex mixed and then stored overnight at −20°C. Centrifuged at 14,000g for 20 min at 4°C and collected supernatants to new tube until LC-MS analysis.

| Bioinformatics analysis for LC-MS
For metabolomics profiling, we used Xevo G2-XS QTOF (Waters, U.K.) to detect metabolites in the samples. Progenesis QI software (Waters, U.K. http://www.nonlinear.com/progenesis/qi/) was used for the preprocessing and identification steps. Each detected metabolic feature was normalized to the QC sample using LOESS Signal Correction (LSC). The correction effect is evaluated by RSD analysis and PCA analysis. The stability of the analytical system is evaluated by intensity distribution of QC samples. Ions with no signal(intensity equals 0) in more than 80% samples of any group are discarded.
The above ions with RSD <30% are included for further analysis.
All filtered metabolomics data were searched against the KEGG database. For statistics analysis, we used multivariate analysis (PCA/ PLS-DA) and univariate analysis (FDR/Fold change) to gain differential ions. Cluster analysis was generated using the 'pheatmap' package in R (v3.0.3). Pathway analysis was performed through the KEGG database.  Figure S1).

| OTU classification and statistics
Based on the OTU abundances, OTUs of each group were listed. The number of OTUs found in HighS, LowS, HighF, and LowF were 2,329, 2,592, 2,265, and 2,574, respectively. A total of 1,929 OTUs were The taxonomic composition distribution in rumen fluid and fecal samples at the phylum level. Note. Abscissa represent samples of rumen fluid or feces and each serial number represent a different cow. S represents rumen fluid samples. F represents fecal samples. The vertical axis represent relative abundance of each phylum shared between HighS and LowS; similarly, 1851 OTUs were shared between HighF and LowF (Supporting Information Figure S2).
The OTU composition was distinctly different between HighS and LowS, while there was no significant difference between HighF and LowF (Figure 1).
For the rumen fluid samples, a total of 20 phyla were identified, and the predominant phyla the species of whose abundance was more than 1% were Bacteroidetes, Proteobacteria, Firmicutes, Spirochaetes, and Cyanobacteria ( Figure 2). As for the fecal samples, a total of 22 phyla were identified, and the predominant phyla Note. HighS represents for rumen fluid of high-production cows. LowS represents for rumen fluid of low-production cows. Wilcoxon Rank-Sum Test is used for two group comparation. If p value is less than 0.05, there is significant difference in alpha diversity between the two groups.
TA B L E 2 Alpha diversity indices of rumen fluid between high-production dairy cows and low-production dairy cows Note. HighF represents for fecal samples of high-production cows. LowF represents for fecal samples of low-production cows. Wilcoxon Rank-Sum Test is used for a two group comparison. If p value is less than 0.05, there is significant difference in alpha diversity between the two groups.
TA B L E 3 Alpha diversity indices of fecal samples between high-production dairy cows and low-production dairy cows F I G U R E 3 Weighted unifrac PCoA analyses of bacteria in rumen fluid and fecal samples. Note. Green triangles represent rumen fluid of high-sproduction cows. Blue dots represent rumen fluid of low-production cows. Red triangles represent fecal samples of high-production cows. Orange squares represent fecal samples of low-production cows were Bacteroidetes, Firmicutes, Spirochaetes, Proteobacteria, Euryarchaeota and Tenericutes (Figure 2 Figure S3).
A total of 95 and 114 genera were detected in the HighF and LowF groups, respectively. Genera with abundances greater than 1% in the HighF and the LowF were almost identical, including 5-7N15, CF231,

| Diversity analysis within and among samples
Alpha diversity (Patrick et al., 2009)  There were significant differences in both species richness and species evenness between HighS and LowS. LowS samples had a higher richness and species evenness than HighS samples (Table 2).
But there were no significant differences in richness or evenness between HighF and LowF (Table 3).
Beta diversity analysis was used to evaluate sample differences in species complexity. PCoA (Principal coordinate analysis) was used to exhibit the differences between the samples according to the matrix of beta diversity distances. PCoA analysis and the clustering results showed that bacterial communities in the rumen fluid were separated from those in the feces (Figures 3 and 4). There were marked differences between HighS and LowS but almost no differences between HighF and LowF samples (Figures 3 and 4). Note. HighS represents for rumen fluid of high-production cows. LowS represents for rumen fluid of low-production cows. Metastats is used for a two group comparison study. If p value is less than 0.05, there is significant difference in alpha diversity between the two groups.

| Significant differences analysis between groups of samples
F I G U R E 5 (a) The difference of Prevotella in rumen fluid between highproduction cows and low-production cows. (b) The taxonomic distribution of genera differing significantly in abundance among rumen fluid samples. Note. HighS represents for rumen fluid of highproduction cows. LowS represents for rumen fluid of low-production cows Dialister in Firmicutes (p < 0.05, Figure 5a,b).

| Metabolic profiling
The PLS-DA model showed a clear separation of samples between high-yield and low-yield dairy cows (Figure 7). The most discriminant metabolites were selected by filtering for fold changes of >1.2 or <0.8, simultaneous with q-value of <0.05 and vip of >1.0.

| D ISCUSS I ON
In agreement with other previous studies, the three dominant phyla observed in all rumen fluid samples were Bacteroidetes, Firmicutes, and Proteobacteria. Contrasted with the higher ratio of Firmicutes to Bacteroidetes in feces, the abundance of Firmicutes in rumen fluid was far less than that of Bacteroidetes, which was consistent with other studies (Jami, Israel, Kotser, & Mizrahi, 2013). It was known that the fiber content in the rumen was far higher than that in the hindgut; thus, we inferred the extra Bacteroidetes present in rumen fluid may be enriched for cellulolytic bacteria. Analogous differences were observed in a recent study on goat (Do et al., 2018).
The study indicated that increasing the members of Bacteroidetes to keep low ratio of Firmicutes versus Bacteroidetes in goat rumen resulted an increased lignocellulose digestion. More interestingly, the high-production cows showed a significant increase in phylum Proteobacteria compared to low-production cows. In fact, the F I G U R E 6 Distribution of taxonomic compositions for significantly different genera in fecal samples. Note. HighF represents for fecal samples of high-production cows. LowF represents for fecal samples of low-production cows The of genera F I G U R E 7 PLS-DA score plot of rumen metabolites between groups. Note. The abscissa represents the first principal component PC1, the ordinate represents the second principal component PC2, and the model parameter R2 is above the graph. Each point in the plot corresponds to an observation. The groups are shown in different colors. Group 1 represents high-production cows. Group 2 represents low-production cows abundance of Proteobacteria in the high-production group was even greater than that of Firmicutes, demonstrating a reversed result with the low-production group and with most previous studies (Jami et al., 2014;Jewell, McCormick, Odt, Weimer, & Suen, 2015).
However, when we further analyzed at the genus level, any genus within Proteobacteria that could account for the marked differences between the two groups was detected, which could be explained by the high unclassified ratio (31%-45%). Much more work still needs to be done to investigate the genus-level difference in phylum Proteobacteria between the high and low-production dairy cows.
Among the identified genera, Prevotella represented the highest percentage in spite of the milk production. The abundance of Prevotella in the high-yield group (37.85%) was lower than that in the low-yield group significantly (47.29%). Prevotella was found negatively associated with RFI in dairy cows (Jami et al., 2014), and the same study also suggested there was a strong negative correlation (Pearson R = −0.69, p = 5 × 10 −3 ) between Prevotella and milk-fat yield.
Moreover, a study on Korean Adolescents showed that Prevotella was associated with triglycerides (TG) and total cholesterol positively, and ultimately induced obesity (Hu et al., 2015). In our study, we did not measure the milk fat ratio of the cows. But the low production cows were fatter than the high production cows generally.
The HighS group was significantly enriched for the genera  TA B L E 6 Discriminant metabolites with down-regulated in the rumen fluid of high production dairy cows and their corresponding metabolic pathways suggested unsaturated fatty acids (UFA) treatments supplemented at 2% of diet DM as either soybean FA distillate or soybean oil increased milk yield, but did not effectively reduce milk fat yield.
Ubiquinone had been suggested to play an important role in the mitochondrial generation of hydrogen peroxide (Boveris, Cadenas, & Stoppani, 1976). However, the metabolic pathways of reduced abundance metabolites in high-production dairy cows were mainly relevant to nucleotide metabolism, energy metabolism, lipid metabolism and biosynthesis of some antibiotics.
The microbiome interacted with the host immune system to regulate metabolism by various mechanisms: direct physical contact, production of metabolites, and shedding of structural components (Zmora at al., 2017). These affected metabolic homeostasis by local mucosal immune modulation and by remote alteration of metabolic organs, such as adipose tissue, muscle, and the liver. It was a pity that we did not detect immune indicators in this study. So we are planning to explore the differences in the blood immunity indices between high-production and low-production dairy cows in the following study.

| CON CLUS ION
In this study, significant bacterial differences were presented between high-yield and low-yield dairy cows, which were mainly reflected by the relative abundances of some special bacteria. Furthermore, there existed significant metabolic differences including biosynthesis of unsaturated fatty acids, steroid biosynthesis, energy metabolism, fatty acid metabolism, amino acid metabolism, biosynthesis of some antibiotics, etc. between the two groups. However, much more work still needs to be done to identify the detailed differences in bacterial abundances between high-yield and low-yield dairy cows.
Accordingly, we can isolate specific beneficial dominant strains in high production cows sequentially to provide material for carrying out microorganism mediated nutritional regulation.

ACK N OWLED G M ENTS
The study was financially supported by the earmarked fund for

CO N FLI C T O F I NTE R E S T
The authors declare no conflict of interest.