Changes in Metabolically Active Bacterial Community during Rumen Development, and Their Alteration by Rhubarb Root Powder Revealed by 16S rRNA Amplicon Sequencing

The objective of this present study was to explore the initial establishment of metabolically active bacteria and subsequent evolution in four fractions: rumen solid-phase (RS), liquid-phase (RL), protozoa-associated (RP), and epithelium-associated (RE) through early weaning and supplementing rhubarb root powder in 7 different age groups (1, 10, 20, 38, 41, 50, and 60 d) during rumen development. Results of the 16S rRNA sequencing based on RNA isolated from the four fractions revealed that the potentially active bacterial microbiota in four fractions were dominated by the phyla Proteobacteria, Firmicutes, and Bacteroidetes regardless of different ages. An age-dependent increment of Chao 1 richness was observed in the fractions of RL and RE. The principal coordinate analysis (PCoA) indicated that samples in four fractions all clustered based on different age groups, and the structure of the bacterial community in RE was distinct from those in other three fractions. The abundances of Proteobacteria decreased significantly (P < 0.05) with age, while increases in the abundances of Firmicutes and Bacteroidetes were noted. At the genus level, the abundance of the predominant genus Mannheimia in the Proteobacteria phylum decreased significantly (P < 0.05) after 1 d, while the genera Quinella, Prevotella, Fretibacterium, Ruminococcus, Lachnospiraceae NK3A20 group, and Atopobium underwent different manners of increases and dominated the bacterial microbiota across four fractions. Variations of the distributions of some specific bacterial genera across fractions were observed, and supplementation of rhubarb affected the relative abundance of various genera of bacteria.

The objective of this present study was to explore the initial establishment of metabolically active bacteria and subsequent evolution in four fractions: rumen solid-phase (RS), liquid-phase (RL), protozoa-associated (RP), and epithelium-associated (RE) through early weaning and supplementing rhubarb root powder in 7 different age groups (1,10,20,38,41,50, and 60 d) during rumen development. Results of the 16S rRNA sequencing based on RNA isolated from the four fractions revealed that the potentially active bacterial microbiota in four fractions were dominated by the phyla Proteobacteria, Firmicutes, and Bacteroidetes regardless of different ages. An age-dependent increment of Chao 1 richness was observed in the fractions of RL and RE. The principal coordinate analysis (PCoA) indicated that samples in four fractions all clustered based on different age groups, and the structure of the bacterial community in RE was distinct from those in other three fractions. The abundances of Proteobacteria decreased significantly (P

INTRODUCTION
The complex rumen ecosystem harbors various prokaryotic (bacteria and archaea) and eukaryotic (protozoa and fungi) microorganisms which symbiotically convert feed stuffs into microbial biomass and fermentation end products that can be utilized by the host (Dehority, 2003;Kong et al., 2010). The ruminal microbiome is characterized by its high population density, wide diversity, and complexity of interactions, and it was suggested that the abundance of various microbial genotypes within the rumen can be significantly bound up with host feed efficiency and diet (Zhou et al., 2010;Carberry et al., 2014). In young ruminants, ingested microbes colonize and establish in a defined and progressive sequence shortly after birth during rumen development (Fonty et al., 1987(Fonty et al., , 2007Abecia et al., 2013). During the last few decades, extensive efforts have been taken to explore the relationship between the microbial colonization and rumen development applying different methods (Fonty et al., 1987;Skillman et al., 2004;Jami et al., 2013;Rey et al., 2014;Guzman et al., 2015;Jiao et al., 2015b), as the composition of ruminal microflora directly influences the digestive and metabolic performance of the host animal, and therefore the developing rumen in the newborn ruminants provides a unique opportunity to manipulate the diverse commensal microflora (Yáñez-Ruiz et al., 2015). It has been reported that early dietary experiences of animal have a greater and more lasting effect than those occurring later in life (Distel et al., 1994), which implies the feasibility to regulate the rumen microbial community at the early period of rumen development, i.e., microbial programming. Li et al. (2012) found that the microbiome in the developing rumen of 14 days old calves was responsive to dietary modifications as well as physiological changes in the host. Further studies (Yáñez-Ruiz et al., 2010;Abecia et al., 2014) suggested that it would be possible to promote different microbial populations inhabiting the rumen of the young animal by controlling the feeding management in early life that persisted in later life. However, more research is essential to gain insight into the development of the rumen and its microbiome, and the method (e.g., the alteration of diets and the supplementation of specific additives) and timing to manipulate the ruminal microflora in early life.
Rhubarb (Rheum spp.) is a commonly used herb in traditional Chinese medicine to cure indigestion, constipation, and other diseases since ancient times, and it has been characterized with antimicrobial and antitumor properties in previous research (Kang et al., 2008;Kosikowska et al., 2010). Previous in vitro and in vivo studies reported that rhubarb could modify rumen fermentation by lowering methane production and the acetate: propionate ratio without inhibiting feed intake or the degradation of substrates (Bodas et al., 2008;García-González et al., 2010Kim et al., 2016). As all these investigations were focused on mature ruminants and rumen fermentation, it would be promising to explore the effects of rhubarb as feed additive in shaping the rumen bacteria community in early life during rumen development.
Within the rumen ecosystem, bacterial species are considered to have a more important role compared to protozoa and fungi in determining the extent and rate of feed degradation and utilization for the production of microbial protein and VFA, which significantly contributes to the maintenance and to the production of meat and milk of the host ruminant (Miron et al., 2001). Hence, improving the understanding of rumen bacterial ecology and acquiring methods to shape the bacterial community may help increase feed efficiency and enhance production performance of animals. The bacterial community in the rumen can been divided into four classes according to their different spatial locations: (1) free-living bacteria associated with the rumen liquid phase; (2) bacteria associated with feed particles; (3) bacteria associated with rumen epithelium; and (4) bacteria attached to the surface of protozoa or contained as endosymbionts inside the protozoal cells (Koike et al., 2003;Yu and Forster, 2005). Amongst the four bacterial groups, the bacteria associated with feed particles, i.e., the solid-phase bacteria, are the most predominant and occupy up to 75% of the total microbial population, and are estimated to produce up to 90% of the endoglucanase and xylanase activities in the rumen (Wang et al., 2011). In contrast, bacteria associated with the liquid-phase (RL) take up 20 to 30% of the total microbes on high-forage rations (Miron et al., 2001). To our knowledge, the majority of the studies on rumen bacteria had been focused on the solid-phase and/or (RL) bacteria (Skillman et al., 2004;Kong et al., 2010;Wang et al., 2011;Jami et al., 2013;Guzman et al., 2015;Jiao et al., 2015b), while research aimed at the bacterial associated with ruminal protozoa and epithelium was limited despite the fact that the presence/absence of protozoa has been found to relate to the structure of different bacterial communities and different rumen fermentation pattern (Yáñez-Ruiz et al., 2007;Belanche et al., 2015). Moreover, the epimural bacterial community which significantly differs from that of rumen contents may influence the extent of development of the rumen epithelium and the immune function (Malmuthuge et al., 2012(Malmuthuge et al., , 2014. Therefore, it would be of great significance to comprehensively explore the simultaneous evolution of the bacterial populations associated with these four fractions during rumen development. Up to now, most of the investigations on the rumen microbiota using 16S rRNA gene sequencing have been conducted based on DNA-derived amplicons (Li et al., 2014;Guzman et al., 2015;Liu et al., 2016;Wang et al., 2016), which could indicate the comprehensive diversity of both the living and inactive microorganisms. Nevertheless, DNA-based studies do not reflect the potential biological activity of the rumen microbial community in real time (Hugoni et al., 2013;Salter et al., 2015). By contrast, RNA-based techniques could help gain insights into the metabolic state of microbes and thus could be used to indicate the most active rumen microorganisms and their metabolic potential (Kang et al., 2013;Li et al., 2016). Lettat and Benchaar (2013) reported that despite the minor disparities in the abundance of different taxa, the RNA-based analysis was more discriminative than the DNA-based counterpart in identifying diet-induced alterations within the microbial community.
In the present study, we aim to deepen the understanding of bacterial colonization during rumen development and offer theoretical foundations for the manipulation of rumen microbiome and fermentation in early life. We used RNAbased 16S rRNA amplicon sequencing to investigate the initial colonization and diversity of metabolically active bacteria in four fractions (i.e., the ruminal solid-phase, RL, protozoa, and epithelium) and the subsequent evolution from 1 to 60 d after birth, and the influences of early weaning and the supplementation of rhubarb on the rumen bacterial population.

MATERIALS AND METHODS
The experiment was approved by the Animal Care Committee, Institute of Subtropical Agriculture, the Chinese Academy of Sciences, Changsha, China.

Animals, Diets, and Management
Forty-five newborn Xiangdong (native breed) black goats (Capra hircus) were used in this study and housed in a wellventilated room with controlled temperature and humidity. The experimental start for each goat was staggered to accommodate differing birth dates. After birth, the goats were left with their dams until weaning. On 1, 10, and 20 d, 8, 7, and 6 goats were slaughtered, respectively. The remaining goats were gradually weaned off goat milk and supplied with ad libitum a mixture of fresh grass (Miscanthus sinensis, 40% of total dry matter [DM]) and starter concentrate (60% of total DM) from 15 d until they were weaned at 40 d. Four goats were further slaughtered at 38 and 41 d, respectively. Sixteen goats were randomly assigned to two diet treatments: the control diet and the diet supplemented with rhubarb root powder, and then reared separately from the dams after weaning. The rhubarb root powder was purchased from a local herbalist retailer and consisted of the dried and milled rhizomes of Rheum offcinale Baill.
The control diet (per kg DM) contained 400 g fresh grass (in DM) and 600 g starter concentrate (in DM), and every 600 g starter concentrate was composed of the following components: 193 g extruded soybean, 69 g whey powder, 100 g maize flour, 109 g fat powder, 80 g soybean meal, 6 g CaCO 3 , 15 g CaHPO 4 , 8 g NaCl, and 20 g premix. Every 1 kg premix contained 2.5 g FeSO 4 •7H 2 O, 0.8 g CuSO 4 •5H 2 O, 3 g MnSO 4 •H 2 O, 5 g ZnSO 4 •H 2 O, 10 mg Na 2 SeO 3 , 40 mg KI, 30 mg CoCl 2 •6H 2 O, 95000 IU vitamin A, 17500 IU vitamin D, and 18000 IU vitamin E. In the control treatment, goats were fed 150 g control diet twice per day at 08.00 and 17.00 h, and four goats were slaughtered before morning feeding separately at 50 and 60 d. In the rhubarb supplemented group, goats were gradually accustomed to the supplementation of rhubarb from 1 week before weaning. Two goats were removed for the reason irrelevant to the experiment, the remaining six goats received 150 g control diet plus 2 g rhubarb root powder per meal, and three goats were slaughtered at 50 and 60 d, respectively. The management of goats and sampling is further illustrated in Figure S1, and the increase of body weight is shown in Figure S2.

Sample Fractionation
After the goats were slaughtered, the rumen was immediately removed for sampling of the four fractions, i.e., the ruminal solid-phase samples (RS), the ruminal (RL) samples (RL), the ruminal protozoa samples (RP), and the ruminal epithelium samples (RE). RS and RL samples were collected and separated using a French press filter (Bodum Inc., Triengen, Switzerland) according to the method described by Kong et al. (2010). To obtain the RP samples a modified method of Leng (1982) was used, briefly 10 mL of rumen fluid was centrifuged at 500 g for 1 min and the protozoal pellet was then rinsed with sterile anaerobic saline solution and collected by centrifugation (500 × g). The rinsing was repeated three times. Three RP samples from each goat were pooled for analysis. For the RE samples, 3 pieces of 2 g (approximately 4 cm 2 ) epithelium samples were excised at different sites of the same rumen and washed with sterile saline solution and then combined. As the rumen was underdeveloped and the contents were limited, no RS and RP sample was collected on 1 and 10 d, and only 4 RL samples were collected on 1 d. 5 RS samples and 4 RP samples were collected on 20 d. All the samples were immediately flash-frozen in liquid nitrogen and then stored at −80 • C for subsequent use.

RNA Extraction and First-Strand cDNA Synthesis
To isolate total RNA a modification of the method described by Wang et al. (2011) was used. Briefly, samples were first manually ground into crude powder in liquid nitrogen using a mortar and pestle, and then 2 g of crude powder was, respectively, weighed and further ground for 5 min in liquid nitrogen using a Retsch RM100 grinder (Retsch GmbH, Haan, Germany). After grinding, 0.3 g frozen fine powder was weighed into each 50-mL tube and mixed with 3 mL of Ambion TRIzol reagent (Life Technologies, Carlsbad, USA). Subsequent procedures were conducted in accordance with the method of Wang et al. (2011). After the extraction, an Ambion MEGAclear kit (Life Technologies, Carlsbad, USA) was used to purify the isolated RNA. The RNA concentration and integrity were estimated using an Agilent 2100 bioanalyzer and RNA 6000 Nano kit (Agilent Technologies, Santa Clara, USA). The prokaryotic total RNA nano assay protocol was used, as prokaryotes account for the majority of RNA in rumen contents (Yu and Forster, 2005).
Five hundred nanogram of isolated total RNA from each sample was used to synthesize the first-strand cDNA using an Invitrogen SuperScript III RT kit (Life Technologies, Carlsbad, USA), and the cDNA synthesis reactions were stored at −20 • C until further analysis was performed.

PCR Amplification and 16S rRNA Amplicon Sequencing
The PCR amplification of bacterial 16S rRNA genes was conducted on a Dyad Peltier Thermal Cycler (AL056543, Bio-Rad Laboratories, Hercules, USA) using primers Bact_341F (5 ′ -TATGGTAATTGTACTCCTACGGGNGGCWGCAG-3 ′ ) and Bact_806R (5 ′ -AGTCAGTCAGCCGGACTACHVGGGTATCTA AT-3 ′ ) (Herlemann et al., 2011;Caporaso et al., 2012). A dual barcode assay adapted for the Illumina MiSeq sequencer (Illumina Inc., San Diego, USA) was used (Table S1). Each primer contained the Illumina adapter sequence, unique barcode, spacer and forward or reverse primer. For each cDNA sample, 20 µL of reaction mix was prepared containing 1 µL cDNA, 1 µL of each barcoded primer (1 µM), 7 µL molecular biology grade H 2 O, and 10 µL KAPA2G Robust Hotstart ReadyMix (Kapa Biosystems, Wilmington, USA). The PCR procedures were as follows: initial denaturation at 95 • C for 5 min; 20 cycles of denaturation (95 • C, 20 s), annealing (55 • C, 15 s) and elongation (72 • C, 5 min); and a final 10-min extension at 72 • C. Each cDNA sample was amplified in duplicates, and three wells per run served as a negative control for the master mix. After amplification, duplicate PCR products were pooled, and the correct sizes of PCR products and the absence of signal from negative controls were further verified through agarose gel electrophoresis. Quantitation of amplicons was performed in a Synergy HTX Multi-Mode Microplate Reader (model SIAFRM, Bio-Tek Instruments Inc., Winooski, USA) using a Quant-iT dsDNA Assay Kit (Thermo Fisher Scientific, Waltham, USA). The amplicons were pooled in equimolar concentrations and purified using Agencourt AMPure XP beads (Beckman Coulter Inc., Brea, USA) and then further quantified as described above. The amplicon library was combined with 10% PhiX control library and sequenced in the Illumina Miseq (Illumina Inc., San Diego, USA) using a v3 600 cycle (300 cycles per read) paired end kit.

Bioinformatic Analysis
The quality of the raw fastq files were checked with the FastQC program (http://www.bioinformatics.bbsrc.ac.uk/projects/ fastqc). Trimmomatic v0.33 (Bolger et al., 2014) was used to trim the raw reads, to remove ambiguous and low quality reads. Reads with average quality score <20 over a 4 bp sliding window and reads with lengths shorter than 36 bp were removed. Merging of the paired-end reads was effected with PEAR v0.9.8 using default options . Reads which did not get assembled were discarded. High quality sequence reads from the various samples were then combined into a single dataset and subsequent analysis was carried out using the open-source software package, QIIME V1.8.0 (Caporaso et al., 2010b). This involved picking Operational Taxonomic Units (OTUs), assigning taxonomy, inferring phylogeny, creating OTU tables and computing microbial community diversity indices. The sequences were clustered into OTUs using the de novo OTU picking workflow with a 97% similarity threshold. Taxonomic assignment of OTUs was performed by comparing the most abundant "representative sequences" within each OTU to the SILVA v119 database (Yilmaz et al., 2014). To enable calculation of Unifrac distances (Lozupone et al., 2011) and to facilitate downstream diversity analysis the picked OTUs were aligned by PyNAST (Caporaso et al., 2010a) against the core alignment template of SILVA v119, and a phylogenetic tree was built using FastTree (Price et al., 2009). To differentiate the conserved from the non-conserved regions of the alignment and remove sections comprised of only gaps (useful in phylogenetic tree construction) a lanemask file was applied. This was constructed from the SILVA v119 core alignment file using a python script. The alpha (within sample) diversity of the samples was estimated using the Chao1, Shannon and observed_otus indices. The Chao1 index was used to further compare the alpha diversity of the samples. Beta (between sample) diversity of the samples was also computed and visualized with three dimensional PCoA plots generated using the Bray-Curtis dissimilarity index (Bray and Curtis, 1957) and the unweighted UniFrac distances. Information on the summary of sequencing data is displayed in the Table S2. All the sequences in the present study were deposited to the sequence read archive (SRA) of the NCBI database using files generated by Mothur V1.33.3 (Schloss et al., 2009), under the accession number SRP081114.

Statistical Analysis
Before the analysis of relative abundance data at phylum and genus level, the compliance of data with the assumptions of normality and homogeneity of variances were first examined visually through residual plots created by the UNIVARIATE and PLOT procedures (SAS Institute Inc., 2001), and variables that were deemed non-normal were then arcsine transformed. Data of bacterial relative abundances at phylum and genus level were analyzed as a completely randomized design using the PROC MIXED procedure of SAS (SAS Institute Inc., 2001) to test the effect of fractions of samples, the model included fraction, age, and fraction × age as the fixed effects, with individual animal as the experimental unit. To test the effect of age on the relative abundances of bacteria at phylum and genus level, the PROC MIXED procedure of SAS (SAS Institute Inc., 2001) was used, with animal nested within age as the random effect and individual animal as the experimental unit. Linear, quadratic, and cubic effects of age were analyzed using orthogonal polynomial contrasts. To compare the relative abundances between the control diet group and the rhubarb group, the PROC MIXED procedure of SAS (SAS Institute Inc., 2001) was used with a model which included the fixed effects of diet, age and diet × age interaction, with individual animal as the experimental unit. Least squares means are reported throughout the text, and statistical significance was declared at P < 0.05.

Rumen Bacterial Community Composition in Four Fractions during Rumen Development and with Supplementation of Rhubarb amongst Individuals
In total, 32 bacterial phyla were identified throughout the four fractions during rumen development. Proteobacteria, Firmicutes, and Bacteroidetes were detected as the three dominant phyla, which together took up from 37.2 to 99.7% of the bacterial microbiota across different fractions in individuals ( Figure  S3). By comparison, the phyla Fibrobacteres, Spirochaetae, Actinobacteria, and Synergistetes were less predominant and their total proportion in individuals ranged from 0.2 to 61.2%. The abundances of Proteobacteria in all the fractions except RS decreased dramatically as the age increased, while the phylum Synergistetes became more abundant on 50 and 60 d than it was before.
In the control diet group at 50 d, Proteobacteria accounted for up to 70.3% of the bacterial community in RE, while its abundance was reduced to as low as 13.1% in RE of the rhubarb treatment group ( Figure S4). No noticeable difference of the bacterial community composition was observed between two treatments on 60 d.
According to taxonomic assignment, a total of 30 bacterial genera with >0.5% of relative abundance amongst four fractions throughout rumen development were observed ( Figure S5). The proportion of the unclassified bacteria and the genera with abundance <0.5 was considerable and ranged from to 4.6 to 90.5% across different fractions in individuals. In the control diet treatment at 50 d, up to 63.0% of the bacterial community in RE was represented by the Suttonella spp. (Figure S6), however its abundance ranged from approximately 0.0-25.1% in RE of the rhubarb treatment. Moreover, in the control diet treatment on 60 d, the abundance of the genus Quinella in RL of individuals ranged from 23.0 to 60.0%, while it decreased to as low as 7.0% in RL of the rhubarb treatment.

Diversity of Rumen Bacteria in Four Fractions during Rumen Development and with Rhubarb Treatment
The Chao 1 index of Alpha diversity was used to measure and compare the bacterial diversity within different days in each fraction (Figure 1). An age-dependent rise of Chao 1 was positively observed in the fractions of RL and RE, but no similar pattern was shown in RS and RP. Furthermore, the comparison of Chao 1 within four fractions on different days revealed that at 50 d the Chao 1 richness in RE was the lowest compared to other three fractions, and at 60 d the richness in RP was less than the other three fractions ( Figure S7). There was no difference in the Chao 1 indices between the control diet and rhubarb treatment on 50 or 60 d (Figure 2).
The beta diversities of bacteria communities within different ages for each fraction were calculated and visualized through three dimensional PCoA analysis using the unweighted UniFrac distances (Figure 3). The samples in all the four fractions clustered based on different age groups. Moreover, the bacterial communities on 38 d and 41 d were separate in RS and RL, while no distinction between these two age groups was found in RP or RE. Difference between the bacterial communities in RE and RL was observed on 1 d and 10 d ( Figure S8), and the community in RE was always different from the communities in other three fractions which clustered together. On 50 d, a clear distinction between the bacterial communities in the control diet treatment and the rhubarb treatment was noted (Figure 4). However, there was no similar treatment-dependent clustering pattern between two treatments on 60 d.

Relative Abundance of Bacteria in Four Fractions during Rumen Development and with Rhubarb Treatment
After the arcsine transformation and subsequent statistical analysis, all the data at phylum level was converted back into the original percent relative abundance and is presented in Table 1. Result of the statistical analysis indicated that the relative abundances of all the main bacterial phyla, except Fibrobacteres, were significantly (P < 0.01) affected by the fraction of samples, and a significant (P < 0.05) interaction between fraction and age was observed for the relative abundances of the phyla Bacteroidetes, Spirochaetae, and Synergistetes. For the phylum Proteobacteria, its abundance in RE was always the greatest compared to other three fractions, and different types of decline in abundances over time were, respectively, noted in RS (cubic, P < 0.05), RL (cubic, P < 0.05), RP (linear, P < 0.01), and RE (quadratic, P < 0.01). Further, the abundances of Proteobacteria in four fractions on 41 d were significantly (P < 0.05) or numerically higher than those on 38 d. The relative abundance of Firmicutes in RE was generally less than those in the other three fractions, cubic (P < 0.05), linear (P < 0.01), and quadratic (P < 0.01) increments with age were shown in RS, RL, and RE, respectively. However, significant (P < 0.05) or numerical decreases of abundances in four fractions from 38 to 41 d were also observed. Age exerted quadratic effects on the abundances of Bacteroidetes in RL (P < 0.01), RP (P < 0.05), and RE (P < 0.05). The relative abundances of Fibrobacteres in RP (P < 0.05) and RE (P < 0.01) both linearly decreased as the age increased. For the phylum Spirochaetae, its abundance in RS was always the highest from 20 to 60 d, and different increases of abundances were found in RL (quadratic, P < 0.05), RP (cubic, P < 0.01), and RE (linear, P < 0.01), respectively. The abundances of Actinobacteria in RS and RP on 41 d were significantly (P < 0.05) improved compared with those on 38 d, and the abundances in RP and RE both increased quadratically (P < 0.05) with age, while a cubic (P < 0.05) effect of age was observed in RL. An increasing quadratic (P < 0.01) effect of age on the relative abundance of Synergistetes was presented in RS, RP, and RE, while the abundance of Synergistetes in RL rose cubically (P < 0.01) with age.
The supplementation of rhubarb root powder significantly (P < 0.05) affected the relative abundance of Proteobacteria throughout four fractions, and a significant (P < 0.05) interaction between diet and age was revealed on the relative abundances of Proteobacteria and Firmicutes ( Table 2). In the rhubarb treatment at 50 d, the abundance of Proteobacteria in RE was significantly (P < 0.05) reduced compared to that of the control diet, while the abundance of Firmicutes in RE and the abundance of Bacteroidetes in RL were significantly (P < 0.05) higher than the control. Moreover, the supplementation of rhubarb significantly (P < 0.05) decreased the abundance of Spirochaetae in RE when compared with the control on 60 d. It is also noted that the abundances of Synergistetes in all the four fractions of the rhubarb treatment on 50 d were numerically less than those of the control diet.
Data at genus level was converted back into the original percent relative abundance after the statistical analysis and is displayed in Table 3. Except the genera Arcobacter, Bibersteinia, Lactobacillus, Ruminiclostridium, Streptococcus, and Fibrobacter, the distributions of all the remaining 24 genera with relative abundances >0.5% were significantly (P < 0.05 or 0.01) influenced by the sample fraction. Furthermore, significant (P < 0.05 or 0.01) interactions between fraction and age were, respectively, found on the relative proportions of the genera Actinobacillus, Alysiella, Arcobacter, Bibersteinia, Mannheimia, Moraxella, Suttonella, Butyrivibrio, Christensenellaceae R-7 group, Quinella, Ruminococcaceae NK4A214 group, Ruminococcus, Streptococcus, Bacteroides, Porphyromonas, Prevotella, Prevotellaceae UCG-001, Fretibacterium, and Treponema. In accordance with the fact that Proteobacteria abundance decreased significantly (P < 0.05) with age, declines to different extents in the relative abundances within four fractions of the genera Actinobacillus, Alysiella, Arcobacter, Bibersteinia, Mannheimia, and Moraxella were observed, respectively. Nevertheless, the proportions of the Ruminobacter spp. in four fractions significantly (P < 0.05) or numerically increased, and the abundances of the Campylobacter (cubic, P < 0.05) and Suttonella (quadratic, P < 0.05) phylotypes in RE rose as the age increased. It is also observed that the fraction of RE accommodated the greatest proportion of Campylobacter spp. compared to other three fractions from 20 to 60 d, and the abundance of Suttonella was significantly (P < 0.05) or numerically higher than those in the other three fractions for most of the days. Within the phylum Firmicutes, the relative abundances of the Lactobacillus spp. and Streptococcus spp. in four fractions generally underwent decreases in different patterns. On the contrary, age had different increasing effects on the relative proportions of the genera Acetitomaculum, Blautia, Lachnospiraceae NK3A20 group, Quinella, Ruminiclostridium, Ruminococcaceae NK4A214 group, and Ruminococcus in four fractions. When compared to 38 d, the abundances of the genera Lachnospiraceae NK3A20 group and Quinella in four fractions on 41 d were significantly (P < 0.05) reduced, while the proportion of Ruminiclostridium spp. was significantly (P < 0.05) enhanced. For the Christensenellaceae R-7 group phylotypes, a cubic (P < 0.01) effect of age on its abundances in RL and RE was presented, while its proportions in RS (P < 0.05) and RP (P < 0.01) reduced quadratically with age. The abundance of Succiniclasticum spp. in RS declined linearly (P < 0.01) with age, while age exerted quadratic (P < 0.01) and linear (P < 0.01) influences on the proportions in RL and RE, respectively. Moreover, RE harbored the most Butyrivibrio spp. but basically the least Quinella spp. compared with the other three fractions, while the abundance of the Blautia phylotypes in RS was the highest from 38 to 60 d. In the Bacteroidetes phylum, the relative proportions of the genera Prevotella and Prevotellaceae UCG-001 in four fractions increased to different extents with age, while FIGURE 4 | Principal coordinate analysis (PCoA) of bacterial community structure using unweighted Unifrac matrix between two treatments on 50 and 60 d.
the decrease in the abundance of Porphyromonas spp. was noted. The proportions of Prevotella in four fractions and Prevotellaceae UCG-001 in RL and RE on 41 d were significantly (P < 0.05) or numerically lower than those on 38 d. For the genus Bacteroides, its abundances in RS (P < 0.05) and RP (P < 0.01) underwent linear decline with age, while a cubic (P < 0.01) effect of age on the abundances in RL and RE was also noticed. For the most of the days during rumen development, the abundance of the Prevotellaceae UCG-001 phylotypes in RE was significantly (P < 0.05) or numerically higher, while the ratio of Prevotella in RE was significantly (P < 0.05) or numerically fewer than those in the other three fractions. In the rest 4 minor phyla, the abundances of Fibrobacter spp. in RP (P < 0.05) and RE (P < 0.01) decreased linearly with age, while different increasing effects of age on the proportions of Treponema in RL, RP, and RE were, respectively, observed. Different manners of increases in the relative abundances of the Atopobium spp. in four fractions were noted, and the abundances in RS, RP, and RE reached their peak on 41 d. In contrast, the proportions of Fretibacterium in all the four fractions quadratically (P < 0.05 or 0.01) increased with age, achieving their maxima on 60 d.
On the genus level, the addition of rhubarb significantly (P < 0.05 or 0.01) influenced the relative proportions throughout fractions of the genera Actinobacillus, Suttonella, Butyrivibrio, Christensenellaceae R-7 group, Lachnospiraceae NK3A20 group, Quinella, and Ruminococcaceae NK4A214 group (Table 4), and significant interaction between diet and age was observed on the abundances of Suttonella (P < 0.05), Lachnospiraceae NK3A20 group (P < 0.05), Ruminiclostridium (P < 0.01), and Selenomonas (P < 0.05). On 50 d in the rhubarb treatment, the proportions of the Campylobacter spp. and Suttonella spp. in RE, and Selenomonas spp. in RL were significantly (P < 0.05) reduced compared with those of the control diet, while the relative abundances of Butyrivibrio spp. and Prevotellaceae UCG-001 in RE, Christensenellaceae R-7 group in RS and RE, Lachnospiraceae NK3A20 group in RS, RP, and RE, and Ruminococcaceae NK4A214 group in RS and RP were greater (P < 0.05) than those of the control. On 60 d, the relative proportions of Actinobacillus, Christensenellaceae R-7 group, and Ruminococcaceae NK4A214 group in RL, and Bibersteinia and Butyrivibrio in RE were significantly (P < 0.05) improved by the supplementation of rhubarb, while the abundances of Quinella in RL and Treponema in RE of the rhubarb treatment were lower (P < 0.05) than those of the control diet. Moreover, the addition of rhubarb significantly (P < 0.05) enhanced the abundance of Ruminiclostridium spp. in RL on 50 d, however a reverse effect on its proportions in RS and RE was observed as well.

DISCUSSION
In the past few decades, research has been performed to explore the initial establishment of the bacterial microbiota in the rumen. By using a culture-dependent technique, Fonty et al. (1987) reported that the rumen of lamb was quickly colonized by an ample microbiota and dominated by the strictly anaerobes on day 2 after birth. The bacterial community was verified to colonize the rumen of goat on 0 d after birth (Jiao et al., 2015b), and more specifically, the existence of Proteobacteria (Geobacter spp.), and fibrolytic bacteria (Fibrobacter succinogenes, Ruminococcus flavefaciens, and Prevotella ruminicola) were detected in the rumen fluid and tissue 0-20 min after birth through realtime qPCR (Guzman et al., 2015). With the aid of amplicon sequencing, it was found that Proteobacteria, Firmicutes, and Bacteroidetes were the three predominant phyla in the rumen fluid in 1-day-old calves (Jami et al., 2013), 2-day-old dairy cattle (Rey et al., 2014), and 7-day-old goats (Wang et al., 2016). One significant similarity amongst those previous studies which employed culture-independent methods is that they were all conducted based on DNA extracted from the rumen contents. In contrast, the current research is the first to investigate the initial colonization and the subsequent evolution of the metabolically active bacterial microbiota, using RNA isolated across four fractions within the rumen of black goats. The taxonomic classification of bacteria in this study reveals that on the first day after birth, the potentially active bacterial populations in RL and RE were both dominated by the Proteobacteria, Firmicutes, and Bacteroidetes, which is in agreement with the previous findings. One of the major findings of the Global Rumen Census (Henderson et al., 2015) was that diet is the main driver of rumen microbial diversity and other studies have noted that individual animals have unique bacterial patterns (Jami and Mizrahi, 2012;Rey et al., 2014;Jiao et al., 2015a).
In the present study, individual variation existed in initial establishment as well as the subsequent evolution of the bacterial microbiome across four fractions. The initial colonization of microbiota in the rumen could stem from the microorganisms in the maternal vagina and skin, colostrum, or the surrounding environment (Dominguez-Bello et al., 2010;Hunt et al., 2011;Di Mauro et al., 2013). These sources and their actual effects which might vary amongst the individuals, and the host genetics could result in the individual variations in bacterial establishment and development (Mayer et al., 2012). Welkie et al. (2010) observed that Holstein cows fed the same ration  possessed substantially different ruminal bacterial community compositions, but displayed similar ruminal pH, VFA profile, and performance (milk production and composition). Using a network bipartition approach for linking the microbial network to a bovine metabolic network, Taxis et al. (2015) reported that two bovine ruminal microbiotas which are rather dissimilar at the taxonomic level were considerably more similar at the metabolic network level. Therefore, it can be assumed that despite the individual variation in bacterial community composition under the same management, it is still possible for the goats to have similar metabolic and production performance in this study.
In the current study, the age of the goats increased along with the alterations of diets. From after birth to 15 d, the dam's milk was the unique source of feed for the goats. From 16 to 39 d, apart from the milk, the goats also had free access to solid starter. After goats were weaned off on 40 d, they were fed only solid starter. Therefore, as was illuminated in previous studies (Rey et al., 2014;Wang et al., 2016), the effect of age in this research actually involved the influence of diet. The Chao 1 index underwent an age-dependent rise both in RL and RE, which is consistent with the report of Wang et al. (2016) that the Chao1 richness of the methanogenic population in rumen fluid increased from 7 days to 2 years. Similarly, an age-dependent increment in Chao 1 of the bacteria community in RE was observed in our previous study (Wang et al., submitted). However, this pattern was not noted in RS and RP and it could be explained by the fact that the Chao 1 indices in these two fractions on 20 d were already relatively high, suggesting that the microbiota from 20 d on were more diverse than it was before. When compared with the microorganisms in RS and RL, knowledge on the microbial community associated with the RE are few. In this study, the Chao 1 richness of the bacterial population in RE was significantly lower than that of the other fractions on 50 d. An explanation for this phenomenon could be that RE locates at the border of host tissues and thus its interaction with different feed, diverse microbes, and complex microscale activities is limited compared to the microbiota in the other fractions, and the environmental heterogeneity could influence the microbial density, diversity, and composition (Horner-Devine et al., 2004;Liu et al., 2016). Moreover, as is reflected by the PCoA, samples in all the four fractions clustered according to different age groups, which is in line with the previous reports on the structure of bacterial population across different ages in rumen fluid (Jami et al., 2013;Rey et al., 2014) and rumen epithelium (Jiao et al., 2015a). This indicates that the structure of bacterial community progressively changed along with the shifts in diets over time, and the microbiota in each age group is distinct. As a strategy to accustom the suckling ruminants to a diet composed of forage and concentrates and lower the cost of production, early weaning and the related strategies have been studied intensively and regarded as an effective approach in regulating the microbial community and improving rumen fermentation (Windeyer et al., 2014;Guzman et al., 2016). In addition to the shift in diet structure and components, weaning could also cause psychological as well as physiological stress since the kids are no longer bred with their mothers (Lay et al., 1998;Loberg et al., 2008). In the present study, the bacterial communities in RS and RE on 38 d were separate from those on 41 d, implying that weaning altered the structure of the bacterial microflora in these two fractions. As is manifested in the PCoA, the structure of bacterial community in RE was always different from those in the other three fractions during rumen development, being supported by those researches in which discrepancies between the epithelial tissue-associated and the rumen contents-associated bacterial populations were observed (Sadet-Bourgeteau et al., 2010;Petri et al., 2013;Liu et al., 2016). Previous studies performed have reported that Proteobacteria, Firmicutes, and Bacteroidetes were the three predominant phyla throughout different age groups (Jami et al., 2013;Rey et al., 2014;Jiao et al., 2015a;Wang et al., 2016). In contrast to the DNA-based 16S rRNA gene sequencing used in those above investigations, the RNA-derived amplicon sequencing could reflect the metabolically active microbes and their potential activities in rumen fermentation (Hugoni et al., 2013;Li et al., 2016). Similarly, results of the current study revealed that the potentially active bacterial microbiota in four fractions was dominated by the phyla Proteobacteria, Firmicutes, and Bacteroidetes regardless of different days. More specifically, the relative proportions of Proteobacteria in RE was consistently the highest compared to other three fractions, which is in agreement with the findings of Chen et al. (2011) andLiu et al. (2016). Rey et al. (2012) reported that the ruminal redox potential value in dairy calves at 1 d after birth was positive but then decreased dramatically and remained negative as the calves aged, and it was also noticed that the ruminal reducing power in lambs rose from 2 to 10 d after birth (Chaucheyras-Durand and Fonty, 2002), which together indicating the foundation and maintenance of a reducing environment in the rumen contents during rumen development. In the present study, the explanation for a greater abundance of Proteobacteria in RE might be that unlike the reducing conditions in the rumen contents, there is trace amount of oxygen in RE, and the phylum Proteobacteria includes many microaerophilic or facultative anaerobic bacteria which are capable of thriving in the environment with trace amount of oxygen (Liu et al., 2016). The consumption of oxygen by the epithelial microflora could contribute to the colonization and activities of the obligate anaerobes in the rumen (Jost et al., 2012). Moreover, RS, RL, and RP accommodated more Bacteroidetes than RE, being supported by previous report (Liu et al., 2016). It is believed that some species within the phylum Bacteroidetes could degrade soluble polysaccharides in plant cell wall, a great proportion of Bacteroidetes in the rumen contents could hence benefit the rumen fermentation of feed matter (Salyers et al., 1995;Power et al., 2014).
In this study, Proteobacteria was the most predominant phylum which occupied 70.25 and 89.05% of the bacterial community in RL and RE on 1 d, respectively. Nevertheless, substantial declines of its proportions in all the four fractions were subsequently observed, while the relative abundances of the phyla Firmicutes and Bacteroidetes increased in inverse proportion to that of the Proteobacteria, being the most and second dominant, respectively. This is partly in line with the findings of a few parallel reports (Jami et al., 2013;Rey et al., 2014;Jiao et al., 2015a). Furthermore, since most of the genera detected within the phyla Firmicutes and Bacteroidetes in current research are anaerobes, the above phenomenon could be defined as a transition from an aerobic or facultative anaerobic to an exclusively anaerobic environmental niche, reflecting the rapid change of the bacterial microflora during normal development of rumen. Meale et al. (2016) observed the increases in the abundances of Proteobacteria and Firmicutes, but the decrease in that of Bacteroidetes in the rumen fluid of dairy calves post weaning vs. before weaning. By contrast, result of the   present study shows that on 41 d the ratios of Proteobacteria in four fractions increased, however most of the proportions of Firmicutes and Bacteroidetes were reduced compared to those on 38 d. This inconsistency could result from the differences in animal species or the methods of weaning used. In contrast to the phyla, the composition and its dynamics of the bacterial microbiota at genus level in the present sequencingbased study were much more diverse and complex. On 1 d after birth, the most prevalent phylum Proteobacteria was mainly dominated by the facultative anaerobic genus Mannheimia which underwent dramatical decline afterwards. Conversely, different types of increases in the relative proportions of the genera Ruminobacter and Suttonella within the phylum Proteobacteria and most of the genera in Firmicutes and Bacteroides were observed. As the age increased and the goats were successively supplemented with solid feed and weaned, the bacterial populations across the four fractions became mainly predominated by the genera Quinella, Prevotella, Fretibacterium, Ruminococcus, Lachnospiraceae NK3A20 group, and Atopobium. Some members of the genus Quinella were found to possess a metabolism of utilizing lactate similar to Selenomonas ruminantium (Krumholz et al., 1993;Belenguer et al., 2010), while the colonization of the genus Fretibacterium were previously reported in human oral cavity (Szafranski et al., 2015). The functional roles of these two genera in rumen fermentation during rumen development are yet to be explored and examined. It has been observed that some ruminal species of genus Prevotella have the capability to utilize starches, simple sugars, and other non-cellulosic polysaccharides as energy sources (Purushe et al., 2010). Within the genus Ruminococcus, R. flavefaciens and R. albus are two major cellulolytic species frequently found in the adult rumen (Flint and Bayer, 2008).
In the current study, the relative abundances of Ruminococcus in four fractions remained low on 1 and 10 d but started to increase from 20 d on, which is in accordance of the initial intake of solid feed on 15 d. Since it was reported that strains of the family Lachnospiraceae detected in the rumen and are potentially able to biohydrogenate fatty acids (Boeckaert et al., 2009), the Lachnospiraceae NK3A20 group found in this research might possess similar potential. Members of the genus Atopobium are Gram-positive and anaerobic bacteria which together consistute a distinct group (Coriobacteria) within the Actinomycetes (Mao et al., 2013), Harmsen et al. (2000) found that the growth of the Atopobium was elevated by sugars. Similarly in present study, the abundances of Atopobium in four fractions were enhanced as the starter concentrate was introduced. Furthermore, the proportions of the genera Selenomonas, Streptococcus, and Fibrobacter across fractions stayed at relatively low level during rumen development in this study, being supported by the finding that these three major genera were identified with relative abundances <1% in the rumen of adult cows (Zened et al., 2013). Henderson et al. (2015) suggested that the pivotal function of Fibrobacter could still be achieved despite of its low relative abundance. The prevalence of the genera Butyrivibrio and Campylobacter in RE was reported in previous researches (Jiao et al., 2015a;Liu et al., 2016), higher proportions of these two genera were also noted in the current study. It is hypothesized that, as butyrate producers, these genera could release butyrate close to the epithelium and thus improve butyrate bioavailability for the host, which might be beneficial to the proliferation of rumen and reticulum epithelium (Siavoshian et al., 2000). Jiao et al. (2015a) proposed that Butyrivibrio might participate in the metabolism of ammonia and volatile fatty acid, and the anatomic development of the rumen. Moreover, it was reported that specific symbiotic microbes could influence the host immune function by regulating host immune responses through proinflammatory and anti-inflammatory pathways, and the epithelial-cell-mediated signals as well (Klaenhammer et al., 2012;Malmuthuge et al., 2014).
In the present study we found that after the introduction of solid starter, the relative abundance of Campylobacter in RE was constantly higher than the other three fractions. The human foodborne pathogen Campylobacter has been previously detected in rumen contents and at levels of up to 30.00% in RE samples (Jiao et al., 2015a). The presence of increased levels of Campylobacter in RE samples is in agreement with its microaerophilic metabolism, and may be an opportunistic, transient member of the RE microbial community. In the present study, when compared to the other three fractions, the abundances of Suttonella and Prevotellaceae UCG-001 in RE were higher, while less relative ratios of Quinella and Prevotella were presented in RE. It is also remarked that RS harbored the most Blautia than other three fractions. These discrepancies amongst fractions necessitate further investigations to explore the mechanisms in the distributions of bacterial populations across different fractions, and deepen the understandings on the commensal interactions between the epimural microflora and the host.
Effects of rhubarb treatment on rumen fermentation have been examined in a few investigations in vitro and in vivo, indicating that rhubarb might manipulate rumen fermentation by reducing methane production and the acetate: propionate ratio without negative side effects, such as suppressing feed intake or digestibility (Bodas et al., 2008;García-González et al., 2010. Recently, it was found that the relative abundances of Prevotella and Lactobacillus were improved after the treatment of rhubarb (Kim et al., 2016). By contrast, result of the present study shows no difference in the abundances of Lactobacillus between the control diet and the rhubarb treatment, and there were only numerical increases in the proportions of Prevotella in rhubarb treatment compared to the control. Further, the addition of rhubarb enhanced the relative abundances of the genera Christensenellaceae R-7 group, Lachnospiraceae NK3A20 group, and Ruminococcaceae NK4A214 group to different extents, but decreased the proportions of Quinella spp. in four fractions. Further investigations are warranted to explore the mechanisms of rhubarb in manipulating the rumen bacterial population and its connection with the methanogen community and methanogenesis. Moreover, in comparison with the 50 d, the effects of rhubarb on the bacterial community structure as well as the relative abundances of bacterial genera were less significant on 60 d in the present study. This phenomenon was also noted in our previous study on the methanogen community during rumen development (Wang et al., submitted), and it could be probably ascribed to the potential rumen microbial adaptation (Hart et al., 2008;García-González et al., 2010.
This study observed individual variations in the metabolically active bacterial community composition, and the rapidly changing bacterial microbiota across fractions in response to the change of diet and age during rumen development. The diversities and structures of bacterial population, and the distributions of diverse bacteria differed across the four fractions during the normal development of rumen, implying that the disparity amongst these four fractions should be taken into account when exploring the overall bacterial community and formulating strategies for microbial programming. The interactions within the complex bacterial microbiota during rumen development and the relevant factors and mechanisms require further investigations. This study provides information on the development of the rumen bacterial community and the related modulation. Future research should try to explain the interactions of the anatomical, functional, and microbial development, as well as the impact of manipulation during early life on ruminant production in the long term.