Characterization of the Intestinal Microbiota and Metabolites of Broiler Breeders with Different Egg Laying Rate

Background The aim of the present study was to determine differences in the intestinal microbiota and metabolites of broiler breeders with different egg laying rate; and its possible relationships between intestinal microbiota and intestinal metabolites were also explored. Methods A total of 200 AA+ parent broiler breeders (41-week-old) were separated into two groups according to their different egg laying rate: average reproductive group (AR: 77.26%±0.88%) and high reproductive group (HR: 86.67%±0.75%), with 10 replicates and 10 birds each 42-days study. Results Feed conversion ratio (FCR), ovary cell apoptosis rate (ApoCR) and abdominal fat pad weight percentage were lower (p < 0.05), while the hatchability of qualied egg rate (HQR) was higher (p < 0.05) in HR group than that in AR group. At phylum level, Firmicutes and Proteobacteria were the dominant bacteria in small intestines, while Bacteroidetes was dominant bacteria in cecum. At genus level, Lactobacillus was dominant in small intestines while Bacteroides were dominant in cecum. In HR Phascolarctobacterium were lower (p < 0.05) in ileum. A positive correlation was observed between Actinobacteria in the small intestines and ApoCR (r = 0.34 & p < 0.05). A positive correlation was observed between FCR and Spirochaetes in cecum (r = 0.54 & p < 0.05). The metabolomics analysis on ileal digesta content revealed a higher concentration of Octadecanoic acid, 9,12-Octadecadienoic acid, 11-trans-Octadecenoic acid, d-Ribose and Inositol in HR broiler breeders. Moreover, inositol was enriched in "Phosphatidylinositol signaling system", "Inositol phosphate metabolism", "Ascorbate and aldarate metabolism", "Galactose metabolism", "ABC transporters", and Octadecanoic acid was enriched in "Fatty acid biosynthesis", of unsaturated fatty acids". results suggest that glycerol, inositol modulate energy partitioning in fat metabolism to inuence reproductive performances. Microbiota, such as Actinobacteria (phylum), Spirochaetes (phylum), and Phascolarctobacterium (genus) may have negative relationship with broiler breeders’ reproductive performances. The present study obtained a total of 5,328,583 high quality 16S rRNA gene sequences from duodenum (1,253,578), jejunum (1,320,806), ileum (1,347,555), cecum (1,406,644), with an average of 83,259 sequences per sample and a range of 67,290 to 95,053. A total of 36,679 OTUs were identied at 97% sequence identity. The Chao1 and ACE indexes were lower in duodenum and ileum of HR (p < 0.05) (Fig 2 a,b) and Shannon index was lower in ileum of HR (p < 0.05) (Fig 2 c). The Venn and Flower diagrams were used to describe common and unique OTUs between goups. In the duodenum, jejunum, ileum and cecum of the AR and HR birds, we observed, respectively, 713, 595 756 and 1046 common OTUs (Fig 3). inositol play a signicant role in fat metabolism though by modulating energy partitioning which can subsequently inuence reproductive performances. Also, higher diversity and richness of microbial communities in the small intestines could induce nutrient competition with the host, negatively impacting reproduction. What’s more, the higher relative abundance of Phascolarctobacterium (genus), Actinobacteria (phylum) and Spirochaetes (phylum) in the AR birds seems to be correlated to their lower reproduction performances.


Sample Collection and Measurements
Egg production, broiler breeder mortality, quali ed hatching egg production and feed consumption were recorded daily. Hatching performance was recorded at day 42. At end of the study (day 42), birds were sacri ced by cervical dislocation and the intestines were removed for the measurement of duodenum (Duo), jejunum (Jej), ileum (Ile), cecum (Cec) weight and length. The content of the duodenum, jejunum, ileum and cecum was immediately transferred into 1.5ml sterile centrifuge tubes for microbiota and metabolomics analysis. The ovaries were collected and placed into 4% paraformaldehyde (pH = 7.2) for cell apoptosis analysis.

Follicle Apoptosis by TUNEL Method
Sections were histochemically stained using the TUNEL technique by an in situ apoptosis detection kit (Roche, Switzerland; [24]. The follicle section of ovary was aseptically removed and placed immediately into a solution of 3-ethoxy-3-aminopropyl silane (Zhongshan Golden Bridge Biotechnology Co., Ltd. Beijing, China). Set into oven to 60 °C for 60 min. Dewaxing and hydration. After an incubation with proteinase K (Boehringer Mannheim, Germany), rinse thrice with PBS. With 50 μL TUNEL in the reaction mixture specimens, with slide in the dark wet box in response to 37 °C for 1 h. Rinse thrice with PBS. Add 50 μL (converter-POD) in liquid samples, or sealing membrane reaction 37 °C in the dark wet box with glass for 30 min. Rinse thrice with PBS. 50-100 μL DAB chromogenic agent, reaction of 25 °C, 10 min. Rinse thrice with PBS Hematoxylin staining, rinse water after a few seconds; alcohol dehydration, xylene transparent, neutral balata. Using BA200Digital (Mike Audi Industrial Group Co., Ltd.) to image acquisition. Apoptotic color is brown yellow, and negative expression is blue with white background.

DNA Extraction and Microbiota Analysis
Total microbial genomic DNA was extracted from the chyme of duodenum, jejunum, ileum and cecum, using the TIANamp Bacteria DNA isolation kit (DP302-02, Tiangen, Beijing, China) according to the manufacturer's instructions. 16S rRNA genes of distinct regions (16S V4) were ampli ed used the speci c primers 515F (5'-GTTTCGGT GCCAGCMGCCGCGGTAA-3') and 806R (5'-GCCAA TGGACTACHVGGGTWTCTAAT-3') with the barcode. All PCR reactions were carried out with Phusion® High-Fidelity PCR Master Mix (New England Biolabs). PCR products, a bright main strip between 400-450 bp, were mixed in equal density ratios for puri cation with Qiagen Gel Extraction Kit (Qiagen, Germany). Sequencing libraries were generated via the TruSeq®DNA PCR-Free Sample Preparation Kit (Illumina, USA) following manufacturer's recommendations, index codes were added. Finally, the library was sequenced on an IlluminaHiSeq2500 platform and 250 bp paired-end reads were generated. Splicing sequences were using FLASH (V1.2.7) to assemble [25], and then quality ltration was using QIIME (V1.7.0) [26] to nish. Chimera tags were removed using UCHIME algorithm [27]. Sequence analysis were performed via Uparse software (Uparse v7.0.1001) [28]. The original GC/MS data were converted into les in netCDF format [32] and subsequently preprocessed, including the identi cation, ltration, and alignment of peaks by XCMS software (www.bioconductor.org). Finally, a data matrix of information, including mass-to-charge ratios (m/z), retention times, and peak intensities, was obtained. Identi cation of metabolites was performed using the AMIDS system (Automated Mass Spectral Deconvolution and Identi cation System) by searching the National Institute of Standards and Technology library (http://srdata.nist.gov/gateway/) [33].
The resultant GC-MS data matrix were imported into the MetaX software package (http://metax.genomics.cn) for multivariate analysis, including principal component analysis (PCA) and partial least squares discriminant analysis (PLS-DA). PCA was utilized to visualize global clustering and display differences in metabolic pro les between samples. To achieve the maximum separation between samples, PLS-DA was used to identify differential metabolites that explained the separation between groups, missing values were imputed using a k-Nearest Neighbor (KNN) method [34]. The model quality was evaluated using the R2X, R2Y, and Q2 parameters. The R2X and R2Y parameters, which represent the fractions of explained X-variation and Y-variation, respectively, can be used to evaluate the quality of the model. The Q2 parameter represents the predictive ability of the model. In general, the model is acceptable when the values of R 2 X, R 2 Y, and Q 2 are greater than 0.5 [35]. In addition to cross-validation, the PLS-DA models were veri ed via 200 repeated permutation tests.
As described above, PLS-DA was employed to identify differential metabolites between groups. Furthermore, variable importance in the projection (VIP) values in the PLS-DA model were used to rank the metabolites based on their importance indiscriminating between groups. Metabolites with the highest VIP values are the most powerful group discriminators. Traditionally, VIP values > 1 are signi cant. In this study, metabolites with VIP values (VIP > 1.0) in the PLS-DA model, signi cance (p < 0.05) in Student's t-test (R package) and fold change < 0.83 or > 1.2 were selected as differential metabolites. Then, receiver operating characteristic (ROC) curves analysis was performed (R package pROC), and the areas under the curves (AUCs) were calculated to determine the diagnostic value of these differential metabolites [36]. The discriminatory power of the differential metabolites was ranked and visualized using heat maps [36].
Metabolic pathways were obtained for these identi ed differential metabolites from the KEGG (Kyoto Encyclopedia of Genes and Genomes) database (http://www.genome.jp/kegg).

Statistical Analysis
Data was analyzed by independent-sample T-test in the SPSS version 25.0 statistical software package (IBM® SPSS® Statistics, New York, USA), differences among treatments were considered signi cant at p < 0.05. Alpha indexes were analyzed by Wilcox rank sum test, differences among treatments were considered signi cant at p < 0.05. Adonis test was based on Bray-Curtis distance, differences among treatments were considered signi cant at p < 0.05.

Reproductive Performance
Average daily feed intake (ADFI) in both HR and AR birds; similarly, no differences in egg weight, quali ed and unquali ed egg rates were observed between the two groups. Feed conversion ratio (FCR), however, was lower in the HR birds (p < 0.05). As expected, the hatchability of quali ed egg rate (HQR) was higher in the HR group (p < 0.05) ( Table 1).

Digestive Organs
Compared to the AR groups, birds in the HR group presented a lower abdominal fat content (p < 0.05), however no differences were observed for chicken weight, relative weight of crop, proventriculus, gizzard, intestines (duodenum, jejunum, ileum, cecum) ( Table 2) and length of intestines (Supplementary Table  1, Additional File 2).

Ovary cell apoptosis
Ovary cell apoptosis rate (ApoCR) in the ovaries of the HR birds was lower than to that observed in the AR group (p < 0.05) (Fig 1).

Summary of High-Throughput Sequencing and Alpha Diversity
The present study obtained a total of 5,328,583 high quality 16S rRNA gene sequences from duodenum (1,253,578), jejunum (1,320,806), ileum (1,347,555), cecum (1,406,644), with an average of 83,259 sequences per sample and a range of 67,290 to 95,053. A total of 36,679 OTUs were identi ed at 97% sequence identity. The Chao1 and ACE indexes were lower in duodenum and ileum of HR (p < 0.05) (Fig 2 a,b) and Shannon index was lower in ileum of HR (p < 0.05) (Fig 2 c). The Venn and Flower diagrams were used to describe common and unique OTUs between goups. In the duodenum, jejunum, ileum and cecum of the AR and HR birds, we observed, respectively, 713, 595 756 and 1046 common OTUs (Fig 3).

Taxonomic Composition of GIT microbiota
A total top 10 phyla and top 10 genera were based on the identi ed OTUs. At phylum level, Firmicutes, Proteobacteria in duodenum, Firmicutes, Cyanobacteria in jejunum and ileum, Firmicutes, Bacteroidetes in cecum were the predominant phyla (Fig 4, a) ( The principal coordinate analysis (PCoA) based on weighted unifrac distance was applied to examine differences in taxonomic community composition and structure in the duodenum, jejunum, ileum and cecum. The PCoA analysis indicates that the bacterial community composition of the duodenum, jejunum, ileum and cecum of the HR birds can be separated from that of the AR group (Fig 5), however this separation did not reach signi cant difference (Table 4).
Spearman correlation analysis was employed to explore correlations between microbiota richness and performance parameters (EPR, FCR, HQR and ApoCR) measured in this study. At phylum level, Actinobacteria in the small intestines were positively correlated with ApoCR. Euryarchaeota in the small intestines were negatively correlated with EPR and HQR and positively correlated with FCR. At genus level, Phyllobacterium, Romboutsia, Megamonas, Coprococcus_1 in the small intestines were negatively correlated with EPR; Phyllobacterium in the small intestines were positively correlated with FCR (Table 5-1). At phylum level, Spirochaetes in the cecum was positive correlated with FCR while Planctomycetes in the cecum was negatively correlated with HQR. At genus level, Sphaerochaeta in the cecum was negatively correlated with EPR (Table 5-2).
A total of 71 metabolites were identi ed in both ileal and caecal samples. Ileal metabolites presented a higher level of expression than the caecal ones (Fig 6,  a). In the ileum Octadecanoic acid, 9,12-Octadecadienoic acid, 11-trans-Octadecenoic acid, d-Ribose and inositol were expressed at higher level in the HR birds, according to their VIP values in the PLS-DA model (VIP > 1) and p-value from Student's t-test (p < 0.05): (Fig 6, b; Table 6). In the cecum, we observed that inositol was lower in the HR birds (Table 6). A ROC analysis was carried out to determine the ability of the identi ed ileal metabolites to discriminate between the AR and HR groups. The heatmap (Fig 6, c) reports the ileal metabolites ranked in order of their AUC values, and indicated that the ve metabolites were effective in discriminating between the HR and AR birds. The KEGG pathway database was utilized to analyze the metabolic pathways relating to the ve metabolites identi ed in the ileum. We observed that inositol was enriched in the "Phosphatidylinositol signaling system", "Inositol phosphate metabolism", "Ascorbate and aldarate metabolism", "Galactose metabolism", "ABC transporters", while octadecanoic acid was enriched in the "Fatty acid biosynthesis", "Biosynthesis of unsaturated fatty acids" pathways (Fig 6, d).

Reproductive performances
The main function of the GIT is to digest feed and absorb nutrients. In the small intestine (duodenum, jejunum and ileum), the host absorbs mostly the nutrients of dietary origin [37], while in the cecum, the host absorbs nutrients that product by microbial fermentation [38,39]. It is well accepted that productive and reproductive performances are in uenced by absorption and utilization of nutrients in poultry [40,41]. In the present study, even though both AR and HR birds were fed the same diet, we observed that abdominal fat was lower in the HR birds, but both body and organ weight was not different from that of AR birds (Table 2). Moreover, FCR was lower in HR birds (Table 1), which suggests that energy and nutrients were used for production purposes rather than being diverted to abdominal fat depots as observed in the HR birds. It has been reported that changes in ovarian function are accompanied by increased abdominal adiposity [42][43][44][45][46] and that body fat is negatively correlated with fertility [47]. These observations support our ndings that ApoCR was higher in the AR birds (Fig 1).

Metabolites
In the present study, with the KEGG enrichment pathway analysis, we found that in "Fatty acid biosynthesis" and "Biosynthesis of unsaturated fatty acids" pathways, hexadecanoic acid (palmitic acid) was lower while octadecanoic acid (stearic acid) was higher in the HR birds ( Supplementary Fig 1, Fig 2, Additional File 2). In agreement with our observations are the nding of Xie et al (2012) where broiler hens had increased abdominal fat pad weight and reduced egg production. it has been observed that high concentrations of palmitic acid can be toxic to granulosa cells [48]. Palmitic acid is a saturated fatty acid that can induce cellular apoptosis and exert cytotoxic effect as it is a precursor of the complex lipid ceramide, which has been reported to induce apoptosis in granulosa cells of hens [49]. Palmitic acid can also induce cellular apoptosis through a ceramide-independent pathway [50]. Finally, palmitic acid and stearic acid could induce apoptosis in human ovarian granulosa cells [51]. In our study we observed higher stearic acid expression in the HR birds, however in light of the lower ApoCR and levels of palmitic acid observed in this group, it is likely that the higher stearic acid expression was not su cient to exert any negative effect on ovarian function in the HR broiler breeders.
We observed that in the "Inositol phosphate metabolism" and "Galactose metabolism" pathways, glycerol was lower while inositol (D-myo-inositol) was higher in the HR birds ( Supplementary Fig 3, Fig 4, Additional File 2). Glycerol can be converted into glycerol-3-phosphate and glyceric acid via glycolysis to supply energy [52] which probably was diverted to the abdominal fat depots in the AR birds. In the "Inositol phosphate metabolism" pathway, phosphatidyl-1D-myoinositol, an intermediate in inositol phosphate metabolism, was also indirectly linked to the "Glycerophospholipid metabolism", which plays an important role in fat metabolism [53]. Although in present study, we observed that inositol was lower in the cecum of the HR birds (Table 6), the overall expression of inositol in the GIT of HR birds was higher than that observed in the AR birds.
Overall, our ndings suggest that intestinal metabolites play a signi cant role in fat metabolism though by modulating energy partitioning which can subsequently in uence reproductive performances.

Microbiota
The nutritional requirements of the gastrointestinal microbiota are species-dependent as bacteria need different dietary substrates for their growth and metabolic activity [39]. Therefore it is obvious that in the gastrointestinal there is a resources competition between the microbiota and the host. In present study, we observed that the duodenum and jejunum microbiota of AR birds has more richness than that of the HR ones (Fig 2 a, b), and the ileum microbiota was more diverse in the AR birds than in the HR one (Fig 2 c). The above mentioned higher diversity and richness of bacteria in the AR birds could have resulted in an increased nutrients competition with the host and the microbiota in the small intestines; the higher FCR observed in the AR birds could support this hypothesis. In present study, Firmicutes and Proteobacteria were the dominant phyla in the small intestines while Bacteroidetes was the dominant bacteria in the cecum. Lactobacillus was thedominant genus in the small intestines while Bacteroides was the dominant genus in the cecum. Similar observations have been made in broiler chickens and laying hens [54,55]. Diet is a well known factor that can in uence the microbiota composition and metabolic activity [56,57]. For example, an altered abundance of Bacteriodetes and Firmicutes has been observed in animals fed a high fat diet [58][59][60], however as the same diet was to the fed AR and HR birds this could probably explain the lack of differences in bacterial communities between the two groups. Indeed, while the PCoA analysis suggested that there might be distinct bacteria communities along the different segments of the gastrointestinal tract of AR and HR birds (Fig 5), these were not signi cant (Table 4).
It is worth noting that we observed a higher relative abundance of Phascolarctobacterium in the ileum of AR birds (Fig 4, c). A higher abundance of Phascolarctobacterium has been correlated with elevated oxidative stress in sows during late gestation and at parturition [61]. Moreover the higher abundance of Phascolarctobacterium has been linked to intestinal in ammation and pathology [62]. Finally, Phascolarctobacterium seems to be correlated with increased body weight and fat mass, and altered metabolic pro le [63]. Therefore, in present study, the higher abundance of Phascolarctobacterium observed in the AR birdsmight be linked to the lower production and reproductive performances observed in the AR birds.
It has been reported that the microbiota can also metabolize steroids and produce hormones [64,65]. For example, Actinobacteria can regulate the sterols: sterol uptake, eliminate the aliphatic side chain at C17, oxidize the steroid nucleus [66] and synthetize testosterone [67,68]. In our study, Actinobacteria abundance in the small intestines was positively correlated with ApoCR (Table 5-1), suggesting that the increased apoptosis cell rate in ovaries could be ascribed to a modi cation in steroids metabolism induced by the increased abundance of Actinobacteria.
In our study we observed that Spirochaetes (phlum) abundance in the cecum was positively correlated with FCR (Table 5-2). Spirochaetes are responsible for infectious typhlitis in hens[69] which is characterized by decreased egg production and increased the feed consumption [70]. Broiler breeders infected by Spirochaetes had poorer feed conversion and the impact of the disease was also extended to their offspring, with chicks being weak, slow growing and with impaired gastrointestinal functionality [71]. In this study, the lower productive and reproductive performances observed in the AR birds could be partially ascribed to the higher relative abundance of Spirochaetes in the cecum of the AR birds (

Conclusions
Overall, our ndings suggest that intestinal metabolites: glycerol, inositol play a signi cant role in fat metabolism though by modulating energy partitioning which can subsequently in uence reproductive performances. Also, higher diversity and richness of microbial communities in the small intestines could induce nutrient competition with the host, negatively impacting reproduction. What's more, the higher relative abundance of Phascolarctobacterium (genus), Actinobacteria (phylum) and Spirochaetes (phylum) in the AR birds seems to be correlated to their lower reproduction performances.

Consent for publication
Not applicable.

Availability of data and materials
The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.

Competing interests
The authors declare that they have no competing interests.             Principal coordinate analysis (PCoA) of microbial comunities in the duodenum (a), jejunum (b), ileum (c), cecum (d), of AR and HR broiler breeders. Figure 6 a: Expression and distribution density of metabolites in ileum and cecum of AR and HR broiler breeders; the Y axis represented the peak-area (log 10), and the X axis represented each groups. b: volcano diagram showing different metabolites in the ileum of AR and HR birds; the Y axis represents the p-value (-log 10), and the X axis represents the fold change (log 2) of metabolites. The size of the spots represent the metabolite VIP value; black spots represent no-different metabolites; red spots represent up-different metabolites. c: heatmap clusering of different metabolites in the ileum of AR and HR broiler breeders. d: KEGG enrichment pathway analysis from different metabolites in ileum of AR and HR broiler breeders; the Y axis represents the enriched pathwat terms, and the X axis represents the ratio of individual metaboltes to the whole identi ed metabolites. The size of the spots represents the number of different metaboltes enriched in each pathway, while the color shade of the spot represented the p-value (-log10) of each pathway.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download.