Transcriptomic and proteomic analyses between small yellow follicles and the smallest hierarchical follicles reveal a role of VLDLR in chicken follicle selection

Background Follicle selection in chicken refers to the process of selecting one from a group of small yellow follicles (SY, 6-8mm in diameter) to enter the 12-15 mm hierarchical follicles (usually F6 follicles), which is a an important process affecting laying performance in the poultry industry. Although transcripromic analysis on chicken ovarian follicles was reported, integrated analysis on chicken follicles around selection by using both transcripromic and proteomic approaches was still lacking. In this study, we compared the proteomes and transcriptomes of SY and F6 follicles of laying hens and found some genes involved in chicken follicle selection. Results Transcriptomic analysis revealed 855 differentially expressed genes (DEGs) between SY follicles and F6 follicles of laying hens, among which 202 were upregulated and 653 were downregulated. Proteomic analysis revealed 259 differentially expressed proteins (DEPs), including 175 upregulated and 84 downregulated proteins. Among the identified DEGs and DEPs, the expression changes of seven genes including VLDLR1,WIF1, NGFR, AMH, BMP15, GDF6 and MMP13 , and nine proteins including VLDLR, VTG1, VTG3, PSCA, APOB, APOV1, F10, ZP2 and ZP3L2 were validated. In addition VLDLR expression was significantly downregulated in F6 follicles compared with SY follicles, was signifcantly higher in the GCs than in the TCs and was stimulated by FSH in GCs of both hierarchical and prehierarchical follicles. comparing transcriptomes SY and F6 laying hens, we identified some differentially expressed proteins/genes that might play certain roles in chicken selection. data may identification of the functional genes and proteins

development and selection.

Background
The ovary is a dynamic organ and a pivotal component of the reproductive system in hens. In the abdomen of laying hens, ovarian follicles of various sizes exist, including small white follicles less than 3 mm in diameter, large white follicles 3-5 mm in diameter, small yellow (SY) follicles 6-8 mm in diameter, large yellow follicles 9-12 mm in diameter and five to six hierarchical follicles with increased sizes, i.e. F6 to F1 [1]. In a reproductively active domestic hen, the hierarchical follicles are recruited daily from a pool of 6-8 mm prehierarchical small yellow follicles, which is called follicle selection [2]. The hierarchical follicles then continue to develop rapidly from the F6 follicle to F1 follicle until ovulation [3]. In the process of follicle selection, the SY follicle being selected absorbs egg yolk at accelerating rate, its granulosa cells rapidly proliferate and differentiate to produce high levels of progesterone [4,5]. The absorption of the egg yolk by the SY follicles is mediated by receptors, and vitellin synthesized by the liver enters the oocytes via the very low-density lipoprotein receptor (VLDLR) [6,7].
Transcriptome analysis of chicken prehierarchical follicles revealed changes in the transcripts involved in steroidogenesis, paracrine signaling and transcription during the early stage of follicular growth and development [8]. Comparison of chicken small white, F1 and post-ovulatory follicles transcriptomes identified differentially expressed genes involved in the adhesion, apoptosis and steroid biosynthesis pathways [9]. Candidate genes including ANXA2, Wnt4 and transforming growth factors play some roles in chicken follicle selection [10][11][12][13][14][15]. However, the dynamics of the transcriptome during chicken follicular selection are unclear. In addition, the mRNA abundance may not accurately predict the quantity of the corresponding functional proteins, while the proteomic approach can provide a systemic overview of protein levels [16]; therefore, has certain advantages over mRNA expression profiling [17].
Proteomic analyses on the ovarian function [18], maturation of oocytes [19], and early embryonic development [20][21][22] in mammals and in human ovarian diseases, such as polycystic ovarian syndrome (PCOS) and cancer [23] were reported, while in chicken, 2889 proteins were identified in the white yolk and ovarian stroma of small white follicles in Bovan's white laying hen [24]. However, the temporal changes in the proteome during chicken follicular selection are unknown. In this study, we compared the proteomes and transcriptomes of 6-8 mm SY follicles and the smallest hierarchical follicles 12-15 mm in diameter in laying hens and found some differentially expressed proteins/genes (DEPs/DEGs) that might play certain roles in chicken follicle selection.

Transcriptomic analysis
RNA-seq was used to compare the transcriptomes of three SY follicles and three smallest heirachical follicles, which are referred to here as S1, S2, S3, and F1, F2 and F3, respectively. High-throughout RNA-seq generated 61.66 Gb clean data for six samples of chicken follicles, and 91.07%-93.42% reads could be mapped to the chicken genome. All six samples had at least 93.55% reads equal to or exceeding Q30 (Table 1).
A total of 855 DEGs, including 202 upregulated and 653 downregulated genes, were identified between the SY follicles (S) and F6 follicles (F) at the significant creteria of (|log2 (FoldChange)| >1 and padj <0.05) (Fig. 1a). A hierarchical clustered map of DEGs was then constructed and is shown in Fig 1b. Detailed analysis of the top 10 up-/down-regulated DEGs is shown in Table 2. The entire list of DEGs is shown in Additional file 1: Table S1.
DEGs were then assessed by GO and KEGG pathway analyses. Go functional analysis revealed that most of the DEGs were assigned to circulatory system process, cell differentiation and transition metal ion binding (Fig.1c). KEGG pathway analysis of DEGs showed that the most enriched pathways are TGF-β signaling pathway, tyrosine metabolism and cytokine-cytokine receptor interaction (Fig. 1d).
To validate the RNA-seq data, seven DEGs (Table 3)

Proteomics analysis
Proteins from three SY follicles and three smallest heirachical follicles, which are used for above transcriptome analysis, i.e. S1, S2, S3, and F1, F2 and F3 follicles, were used for TMT labeling and HPLC fractionation followed by LC-MS/MS analysis.
The first step was to validate MS data. The distribution of mass error was close to zero, and most of the absolute values were less than 5 ppm, meaning that the mass accuracy of the MS data is compliant with the requirements (Fig. 3a). The length of most peptides was between eight and 16 amino acids, in agreement with the general characteristics of tryptic peptides (Fig. 3b). In this study, a total of 5883 proteins were identified in the samples and 5236 proteins were quantified.
According to relative levels, the quantified proteins were divided into two categories: a quantitative ratio over 1.5 was considered up-regulation, and a quantitative ratio less than 1/1.5 was considered down-regulation (P < 0.05) (Fig.   3c). In the F6 follicles, the levels of 175 and 84 proteins significantly increased and decreased, respectively, compared with SY follicles. Detailed analysis of the top 10 up-/down-regulated differentially expressed proteins (DEPs) is shown in Table 4. The entire list of DEPs is shown in Additional file 2: Table S2.
The pathways of DEPs were constructed using the KEGG software. Several important pathways are enriched in the F6 follicles compared with the SY follicles (Fig. 3d).
The DEPs were mainly enriched in the ribosome, neuroactive ligand-receptor interaction pathway and cytokine-cytokine receptor interaction.    Table 5). The PRM results ( Fig. 4) showed that the relative abundance of the peptides from above selected nine individual proteins is consistent with the proteome data.  Fig. 4 The histogram of nine significantly abundant proteins in F6 follicles (F) vs SY follicles (S) using PRM (P < 0.05).

Transcriptome and proteome association analysis
The association analysis of the protepme and transcriptome data of F6 and SY follicles revealed a weak relationship between protein and mRNA expression with a Pearson's correlation coefficient of 0.23 (Fig. 5). To further understand the relationship between transcripts and proteins, we compared the intersection between differentially expressed genes (DEGs) and differentially expressed proteins(DEPs) (Fig. 6). Most genes were significantly expressed at mRNA level but not at protein level. At both the protein and mRNA levels, 14 and 26 genes were revealed as significantly up-and down-regulated in SY follicles as compared with that in F6 follicles, respectively. In addition, the expression of two genes is inconsistent at changes in mRNA levels and protein levels. Table 6 shows specific regulation information of some genes at mRNA and protein levels. The specific comparative analysis results are shown in Additional file 3: Table S3.

Dynamics and regulation of VLDLR mRNA by FSH
Both transcriptomic and proteomic analyses indicated that VLDLR expression was significantly down-regulated in F6 follicles compared with SY follicles, therefore, we further analyzed the expression of VLDLR mRNA in chicken tissues and found that chcken VLDLR is predominantly expressed in the ovary (Fig. 7a), and its expression level in the prehierarchical follicles was signifcantly higher than that in the hierarchical follicles (P < 0.01) (Fig. 7b). In both the hierarchical and prehierarchical follicles, the VLDLR were signifcantly higher in the granulosa cells (GCs) than in the thecal cells (TCs) (P < 0.05) (Fig. 7c). Follicle-stimulating hormone treatment increased the expression of VLDLR in the GCs of both hierarchical and prehierarchical follicles (P < 0.01) (Fig. 7d, 7e). Through bioinformatics analysis and qRT-PCR validation on the DEGs identified by RNA-seq, hundreds of DEGs were revealed to be involved in follicular development, some of which have been shown to be involved in follicular development and selection (Fig. 2). Among these DEGs, anti-mullerian hormone (AMH) is a member of the TGF-β superfamily and is mainly expressed in granulosa cells of 1-5 mm follicles in the early stage of follicular development [1]. MMP13, also known as collagenase-3, initiates the breakdown of thefibrillar collagens that form a key structural element of membranes. Our previous study also revealed that MMP13 is mainly expressed in chicken F5 and POF1 follicle and polymorphisms in the promoter region are associated with age at the first laying trait in laying hens [25]. Bone morphogenetic proteins 15 (BMP15) may promote follicle selection and effects on granulosa cells (GCs) proliferation and steroidogenesis in hens [13,26]. In this study, we also found that, compared with SY follicles, the expression of AMH and BMP15 were significantly decreased in F6 follicles, while that of MMP13 was significantly increased (Table 3). We also detected significant changes in the expression of tumor necrosis factor receptor superfamily member 16 (NGFR), Wnt inhibitory factor 1 (WIF1) and growth differentiation factor 6 (GDF6), the function of which in chicken follicle selection requires further investigation. follicles compared with SY follicles. Apovitellenin-1 is one of five apoproteins that forms low-density lipoproteins (LDL) and is expressed in the egg yolk and vitelline membrane [27,28]. Apolipoprotein B is the major apolipoprotein of very low-density lipoproteins (VLDL) and LDL and is a ligand for LDL receptors. The increased expression of apovitellenin-1 and apolipoprotein B in 12-15 mm chicken hierarchical follicles is consistent with higher yolk incorporation after follicle selection. Prostate stem cell antigen is reported to be involved in some functions of stem cells, such as self-renewal and proliferation [29]. Coagulation factor X is a vitamin K-dependent serine protease that plays a key role in the common coagulation pathway [30,31].
Vitellogenins are the yolk precursor proteins produced by the liver and are very important for growth of the chicken ovarian follicles. Vitellogenins circulate in the bloodstream until a follicle enters a stage of vitellogenesis that triggers endocytosis of vitellogenins and transport into the yolk [7]. Vitellogenins 1, 2, and 3 were found in the chicken egg yolk plasma and granules [28]. The increased expression of vitellogenins 1 and 3 in F6 follicles from this study is predictable because yolk deposition is active in F6 follicles. The zona pellucida (ZP) is a specialized extracellular matrix that surrounds the oocyte and early embryo. It is composed of three or four glycoproteins including zona pellucida sperm-binding proteins 1-4 (ZP1-4) with various functions in oogenesis and fertilization [32]. zona pellucida sperm-binding protein 3 can promote oocyte maturation in pigs [33]. In this study, we found that the expression levels of zona pellucida sperm-binding proteins 2 and 3 are grdually decreased.
Very low-density lipoprotein receptor (VLDLR) plays a pivotal role in chicken reproduction, and without VLDLR, oocytes are unable to enter the rapid growth stage of follicle development [34]. During the development of the small white follicle, the VLDLR migrates to the follicular wall, enabling endocytosis of vitellogenin into the yolk, followed by follicular differentiation [7]. In chicken granulosa cells, the expressed variant of VLDLR differs from the variant expressed in the oocyte which contains an O-linked sugar domain (VLDLR 1) [36,37]. The lower level of VLDLR in granulosa cells is suggested to allow more VLDLR to pass through the intercellular gaps to reach the oocytes rather than receptor mediated endocytosis into the granulosa cells [35]. In this study, by RNA-seq and qRT-PCR, we found that VLDLR1 mRNA is significantly down-regulated, while proteomic analysis and PRM validation both indicated that VLDLR protein mRNA and protein is significantly down-regulated. Further analysis on VLDLR expression indicated that it is mainly expressed in chicken ovary, in SW and SY follicles and in the GCs of prehierachical follicles (Fig. 7). Moreover, its expression in the GCs of both prehierachical and hierachical follicles was stimulated by FSH (Fig. 7d, 7e). These data collectively suggest an important role of VLDLR in chicken follicle selection.
Similarly, in geese, the expression of VLDLR mRNA is decreased concomitant to an increase in follicular diameter [38].
In this study, by association analysis between the transcriptome and the proteome, we found that the expression changes of most genes are inconsistent at the mRNA and protein levels, which is consistent with other studies showing that most genes can exhibit differential expression at the transcript and protein level [39][40][41][42][43], suggesting that it is necessary to analyze the function of genes at both levels.

Conclusions
By comparison between the transcriptomes and proteomes of SY follicles and F6 follicles in laying hens, this study revealed 855 DEGs and 259 DEPs. The possible functional significance of some DEGs and DEPs including VLDLR was further discussed. These data may contribute to identification of the functional genes and proteins involved in chicken follicular development and selection.

Tissue collection
Fifteen randomly sampled Hy-line brown hens were randomly divided into three

Transcriptome analysis
Total RNA was isolated with TRIzol reagent(Invitrogen, UK) in accordance with the manufacturer's instructions. Then RNA purity and integrity were checked using the NanoPhotometer® spectrophotometer(IMPLEN, CA, USA) and the RNA Nano 6000 Assay Kit of the Bioanalyzer 2100 system (Agilent Technologies, CA, USA), respectively. Subsequently, mRNA was purified from total RNA using poly-T oligoattached magnetic beads. Fragmentation was carried out using divalent cations under elevated temperature in NEBNext First Strand Synthesis Reaction Buffer(5X). Genes with an adjusted P-value <0.05 and |log 2 FoldChange| > 1 found by DESeq2 were assigned as differentially expressed. Gene Ontology (GO) enrichment analysis and Kyoto Encyclopedia of Genes and Genomes database (KEGG) pathway enrichment analysis of differentially expressed genes were implemented by the clusterProfiler R package,and the P-value less than 0.05 were considered significantly enriched by differential expressed genes.

Spectrometry (LC-MS/MS)
Proteins in the tissues were extracted with lysis buffer containing 8 M urea (Sigma) and 1% protease inhibitor cocktail (Calbiochem) and the protein concentration was determined with a BCA kit (Beyotime) according to the manufacturer's instructions.
Protein enzymolysis was performed using trypsin (Promega ion charges were set as 1, 2; and ion types were set as b, y, p. The product ions were derived from ion 3 to the last ion, and the ion match tolerance was set to 0.02 Da.

Real-time quantitative PCR
The total RNA was extracted from the S and F ovarian follicles that were used for proteome analysis using TRIzol reagent (Invitrogen  Table S4.

Statistical analysis
Differences in mRNA expression between fillicles and follicular cells were evaluated by one-way ANOVA followed by Duncan's multiple range test (P < 0.05) using the General Linear Model procedure of SAS (version 9.2). In one experiment, each treatment was repeated four times, and at least three independent experiments were performed. All data were presented as the mean ± SEM (n = 4). The GO and

Availability of data and materials
For transcriptome sequencing, the raw reads have been submission to the Sequence performed the data processing and biological information analysis. YJ, LK and YS conceived the study and the experimental design and helped draft the manuscript.
All authors read and approved the final manuscript.

Ethics approval and consent to participate
The animal experiments were carried out in accordance with the protocols of the 'Guidelines for Experimental Animals' of the Ministry of Science and Technology (Beijing, China) and all efforts were made to minimize suffering. The animal experiments were approved by the Institutional Animal Care and Use Ethics Committee of Shandong Agricultural University (Permit Number: 2,007,005). We have obtained the informed consents of a written form from Mr. Xiaoyan Yang who was the owner of the Hy-line brown hens.

Consent for publication
Not applicable.

Figure 2
The mRNA expression levels of genes examined by qRT-PCR. All data are presented as the m Figure 3 TMT analysis of the differentially expressed proteins (DEPs) data in chicken F6 and SY follicle The histogram of nine significantly abundant proteins in F6 follicles (F) vs SY follicles (S) usin Association analysis of proteome and transcriptome differences between F6 and Sy follicles i Figure 6 Venn diagram of differentially expressed genes/proteins from TMT and DEGs analyses.