Transcriptomic and proteomic analyses of ovarian follicles reveal the role of VLDLR in chicken follicle selection

Abstract

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 that are less than 3 mm in diameter, large white follicles that are 3-5 mm in diameter, small yellow (SY) follicles that are 6-8 mm in diameter, large yellow follicles that are 9-12 mm in diameter and five to six hierarchical follicles of increased sizes, i.e., F6 to F1 [1]. Follicle selection in reproductively active domestic hens refers to the daily collection of one follicle from a pool of 6-8 mm small yellow follicles, which becomes a hierarchical follicle [2] and continues to develop rapidly from the F6 follicle stage to the F1 follicle stage until ovulation. In the process of chicken follicle selection, granulosa cells rapidly proliferate and differentiate to produce high levels of progesterone [3, 4], and vitellin synthesized by the liver enters oocytes via the very low-density lipoprotein receptor (VLDLR) [5, 6].

Changes in the transcripts involved in steroidogenesis, paracrine signaling and transcription during the early stage of follicular growth and development were identified by transcriptome analysis [7]. Comparison among the transcriptomes of small white, F1 and postovulatory chicken follicles identified differentially expressed genes that are involved in the adhesion, apoptosis and steroid biosynthesis pathways [8]. Candidate genes, including ANXA2, Wnt4 and transforming growth factor genes, were shown to play several roles in chicken follicle growth [9–14]. However, the dynamics of the transcriptome during chicken follicle development from the SY follicle to the F6 follicle are unclear.

In addition, the mRNA abundance may not accurately predict the quantities of the corresponding functional proteins, while a proteomic approach can provide a systemic overview of protein levels [15]; therefore, a proteomic approach has certain advantages over mRNA expression profiling [16]. Proteomic analyses of ovarian function [17] and maturation of oocytes [18], such as polycystic ovarian syndrome (PCOS) and cancer [19], in human ovarian diseases and early embryonic development [20–22] in mammals were reported, while in chicken, 2,889 proteins were identified in the white yolk and ovarian stroma of small white follicles in Bovan's white laying hen [23]. However, the temporal changes in the proteome during chicken follicle selection are unknown. In this study, we compared the proteomes and transcriptomes of 6-8 mm SY follicles and the smallest hierarchical follicles (F6) in laying hens and found several differentially expressed genes/proteins (DEGs/DEPs) that might play certain roles in chicken follicle selection.

Results

Transcriptomic analysis

RNA-seq was used to compare the transcriptomes of three SY follicles and three of the smallest hierarchical follicles (F6), which are referred to here as S1, S2, and S3 and F1, F2 and F3, respectively. High-throughput RNA-seq generated 61.66 Gb of clean data from the six samples of chicken follicles, and 91.07%-93.42% of the reads could be mapped to the chicken genome. For all six samples, at least 93.55% of the reads were equal to or exceeded 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) according to the significance criteria of |log₂ (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-/downregulated DEGs is shown in Table 2. The entire list of DEGs is shown in Table S1.

The DEGs were then assessed by GO and KEGG pathway analyses. The GO functional analysis revealed that most of the DEGs were involved in circulatory system processes, cell differentiation and transition metal ion binding (Fig. 1C). KEGG pathway analysis of the DEGs showed that the most enriched pathways were those involved in TGF-β signaling, tyrosine metabolism and cytokine-cytokine receptor interactions (Fig. 1D).

To validate the RNA-seq data, seven DEGs (Table S2), including very low density lipoprotein receptor 1 (VLDLR1), nerve growth factor receptor (NGFR), WNT inhibitory factor 1 (WIF1), anti-Mullerian hormone (AMH), bone morphogenetic protein 15 (BMP15), growth differentiation factor 6 (GDF6) and matrix metallopeptidase 13 (MMP13), were chosen and quantified by quantitative real-time PCR (qRT-PCR). The results showed that the mRNA levels of these genes were similar to those revealed by the sequencing data, suggesting that the RNA-seq results were reliable (Fig. 2).

Table 1 Summary of the RNA-seq metrics for chicken follicles.

Table 2 The top 10 up- and downregulated genes of chicken F6 vs SY follicles.

Proteomics analysis

The proteins from three SY follicles and three of the smallest hierarchical follicles (F6) that were used for the above transcriptome analysis, i.e., the S1, S2, and 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 the MS data. The distribution of the mass error was close to zero, and most of the absolute values were less than 5 ppm, which meant that the mass accuracy of the MS data was compliant with the requirements (Fig. 3A). The length of most peptides was between eight and 16 amino acids, which was in agreement with the general characteristics of tryptic peptides (Fig. 3B). In this study, a total of 5,883 proteins were identified in the samples, and 5,236 proteins were quantified. According to the relative levels, the quantified proteins were divided into two categories: proteins with a quantitative ratio over 1.5 were considered upregulated, and proteins with a quantitative ratio less than 1/1.5 were considered downregulated (P < 0.05) (Fig. 3C). In the F6 follicles, the levels of 175 and 84 proteins were significantly increased and decreased, respectively, compared with those in the SY follicles. A detailed analysis of the top 10 up-/downregulated differentially expressed proteins (DEPs) is shown in Table 3. The entire list of DEPs is shown in Table S3.

The pathways of the DEPs were constructed using KEGG software. Several important pathways were enriched in the F6 follicles compared with the SY follicles (Fig. 3D), including pathways involved in the ribosome, neuroactive ligand-receptor interactions and cytokine-cytokine receptor interactions.

Table 3. Top 10 up- and downregulated proteins in chicken F6 vs SY follicles.

Nine DEPs were randomly selected for parallel reaction monitoring (PRM) analysis to verify the accuracy of the proteome analysis by LC-MS/MS, including apovitellenin-1 (APO1), apolipoprotein B (APOB), prostate stem cell antigen (PSCA), coagulation factor X (F10), vitellogenin-1 (VTG1) and vitellogenin-3 (VTG3); all of these proteins were significantly increased in F6 follicles, and zona pellucida sperm-binding protein 2 (ZP2), zona pellucida sperm-binding protein 3 (ZP3) and very low-density lipoprotein receptor (VLDLR) were significantly decreased in F6 follicles (Table 4). The PRM results (Fig. 4) showed that the relative abundances of the peptides from the nine selected individual proteins were consistent with the proteome data.

Table 4. Nine proteins selected for the parallel reaction monitoring analysis of the chicken follicle proteome data

Transcriptome and proteome association analysis

The association analysis of the proteomic and transcriptomic data of the F6 and SY follicles revealed a weak relationship between protein and mRNA expression with a Pearson’s correlation coefficient of 0.23 (Fig. 5A). The number of items for which “Transcript up and Protein up” and “Transcript down and Protein down” was 14 (1.3%) and 26 (2.4%), respectively (Fig. 5B). To further understand the relationship between transcripts and proteins, we compared the intersection between DEGs and DEPs (Fig. 6). Most genes were significantly expressed at the mRNA level but not at the protein level. At both the protein and mRNA levels, fourteen and twenty-six genes were revealed as significantly up- and downregulated in SY follicles compared with F6 follicles, respectively. In addition, the expression of two genes were inconsistent in terms of changes in the mRNA levels and protein levels. Table 5 shows the specific regulation information for several genes at the mRNA and protein levels. The specific comparative analysis results are shown in Table S4.

Table 5 Differentially expressed genes at both the mRNA and protein levels between SY and F6 follicles in chicken

Dynamics and regulation of VLDLR mRNA by FSH

Both transcriptomic and proteomic analyses indicated that VLDLR expression was significantly downregulated in F6 follicles compared with SY follicles; therefore, we further analyzed the expression of VLDLR mRNA in chicken tissues and found that it was predominantly expressed in the ovary (Fig. 6A), and VLDLR expression in the prehierarchical follicles was significantly higher than that in the hierarchical follicles (P < 0.01) (Fig. 6B). In both the hierarchical and prehierarchical follicles, the VLDLR mRNA expression was significantly higher in the granulosa cells (GCs) than in the thecal cells (TCs) (P < 0.05) (Fig. 6C). Follicle-stimulating hormone (FSH) treatment stimulated the expression of VLDLR in the GCs of prehierarchical follicles in a dose-dependent manner (P < 0.01) (Fig. 6D), while in hierarchical follicles, the effect was not significant at concentrations of ≤10 ng/ml (Fig. 6E).

Discussion

Follicle selection is an important stage of follicle development that affects many egg production traits in the poultry industry. The mechanism of follicle selection in chickens is becoming a hot topic of research in poultry reproduction biology. Previous studies revealed the effect of Wnt4 [9], bone morphogenetic protein 4 [11], BMP15 [12], AMH [24], vasoactive intestinal peptide [25] and parathyroid hormone-like hormone [26] on chicken follicle selection and the changes in the RNA N6-methyladenosine methylation profile [27] during chicken follicle selection; however, high-throughput screening of functional genes at the protein level involved in chicken follicle selection is lacking. Therefore, in this study, by combined analysis of changes in the transcriptomic and proteomic profiles of follicles prior to and post selection in chickens, we revealed several DEGs and DEPs, including VLDLR, that may play important roles in chicken follicle selection.

At both the mRNA and protein levels, the expression of VLDLR, NGFR and WIF1 were significantly downregulated in chicken F6 follicles compared with that in SY follicles. During chicken follicle selection, VLDLR plays a pivotal role in the absorption of vitelin by oocytes, and without VLDLR, oocytes are unable to enter the rapid growth stage of follicle development [28]. During the development of the small white follicle, VLDLR migrates to the follicular wall, enabling the endocytosis of vitellogenin into the yolk, followed by follicular differentiation [7]. Studies have revealed that the expressed variant of VLDLR in chicken granulosa cells differs from the variant expressed in the oocyte, which contains an O-linked sugar domain (VLDLR 1) [29, 30], and the reduced level of VLDLR in granulosa cells is suggested to allow more VLDLR to reach the oocytes by passing through intercellular gaps rather than via receptor-mediated endocytosis into granulosa cells [31]. In this study, we provided further evidence that after follicle selection, the expression of VLDLR was significantly decreased. Moreover, we found that VLDLR is mainly expressed in chicken ovaries, SW and SY follicles and the GCs of prehierarchical follicles, and its expression was stimulated by FSH in the GCs of prehierarchical follicles. These data collectively suggest that, after follicle selection, the decreased expression of VLDLR in GCs might allow more VLDLR 1 to be expressed on the oocyte membrane, thus promoting the rapid growth of hierarchical follicles. Similarly, in geese, the expression of VLDLR mRNA is decreased concomitant with an increase in the follicular diameter [32]. For NGFR, in human mural and cumulus granulosa cells, the nerve growth factor receptor tropomyosin-related kinase A (TrkA) mRNA level was strongly correlated with the number of oocytes retrieved, and the number of oocytes retrieved was greater among women with a low p75(NTR)/TrkA ratio [33]. The role of WIF1 has not been reported in ovarian function or follicle growth. The function of NGFR and WIF1 in chicken follicle selection requires further investigation.

Hundreds of DEGs were identified by transcriptome comparison between SY and F6 follicles during chicken follicle selection. Among these DEGs, the expression of AMH and BMP15 was greatly decreased in F6 follicles compared with that in SY follicles. Similar dynamics for AMH and BMP15 expression were previously reported; AMH is mainly expressed in granulosa cells of 1-5 mm follicles in the early stage of follicular development [1] and is decreased in follicles after selection [24], while BMP15 may promote follicle selection and affect granulosa cell proliferation and steroidogenesis in hens [12].

Among the 259 DEPs identified by global proteome analysis of chicken SY and F6 follicles, the following proteins may play certain roles in chicken follicle selection. APO1 is one of five apoproteins that forms low-density lipoprotein (LDL) and is expressed in the egg yolk and vitelline membrane [34, 35]. APOB is the major apolipoprotein in very low-density lipoprotein (VLDL) and LDL and is a ligand for LDL receptors. The increased expression of APOV1 and APOB in F6 follicles is consistent with higher yolk incorporation after follicle selection. Vitellogenins are yolk precursor proteins produced by the liver and are essential for the growth of chicken ovarian follicles, and they circulate in the bloodstream until a follicle enters the stage of vitellogenesis, which triggers endocytosis of vitellogenins and transport into the yolk [6]. Vitellogenins 1, 2, and 3 were found in chicken egg yolk plasma and granules [35]. The increased expression of VTG1 and VTG3 in F6 follicles found in this study suggests that yolk deposition is active in F6 follicles. The zona pellucida (ZP) is a specialized extracellular matrix that surrounds the oocyte and early embryo and is composed of three or four glycoproteins, including zona pellucida sperm-binding proteins 1-4 (ZP1-4), with various functions in oogenesis and fertilization [36]. ZP3 has been shown to promote oocyte maturation in pigs [37]. In this study, we found that the expression levels of ZP2 and ZP3 were decreased to approximately half of the level found in F6 follicles; the mechanism of this process in follicle selection requires further study.

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 in line with the results of other studies showing that most genes can exhibit differential expression at the transcript and protein levels [38, 39], and this suggests that it is necessary to analyze the functions of genes at both levels.

Conclusions

By comparing 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 several DEGs and DEPs, including VLDLR, was further discussed. These data may contribute to the identification of the functional genes and proteins involved in chicken follicular development and selection.

Methods

Tissue collection

Fifteen randomly sampled Hy-line brown hens from the same batch were divided into three biological groups, with five hens in each group. These hens had been laying regularly for at least one month (28 weeks old, with a mean body weight of 2.1± 0.12 kg) and were housed under standard conditions with free access to food and water. Vaccination was performed according to the recommendations from Hy-line International. The laying time of the sampled hens were recorded and the hens were killed by cervical dislocation after approximately ten hours after laying. From each hen, small yellow follicles 6-8 mm in diameter and hierarchical follicles 12-15 mm in diameter were separately collected, and the egg yolk was carefully squeezed out with tweezers, washed with phosphate-buffered saline (Thermo Fischer Scientific, MA, USA), immediately frozen in liquid nitrogen and used for transcriptomic and proteomic analyses. The Institutional Animal Care and Use Ethics Committee of Shandong Agricultural University reviewed and approved all procedures described in this study. This study was performed according to the Guidelines for Experimental Animals of the Ministry of Science and Technology of China. Three biological replicates were prepared for the transcriptomic and proteomic analyses of total RNA and proteins.

Transcriptome analysis

Total RNA was isolated with TRIzol reagent (Invitrogen, CA, USA) in accordance with the manufacturer’s instructions. Then, the 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 the total RNA using poly-T oligo-attached magnetic beads. Fragmentation was carried out using divalent cations at an elevated temperature in NEBNext First Strand Synthesis Reaction Buffer (5X). The first-strand cDNA was synthesized using random hexamer primers and M-MuLV Reverse Transcriptase (RNase H-). Second-strand cDNA synthesis was subsequently performed using DNA Polymerase I and RNase H. The remaining overhangs were converted into blunt ends via exonuclease/polymerase activities. After adenylation of the 3’ ends of the DNA fragments, NEBNext Adaptors with hairpin loop structures were ligated to prepare for hybridization. To select the cDNA fragments that were 250~300 bp in length, the library fragments were purified with the AMPure XP system (Beckman Coulter, Beverly, USA). The library quality was assessed on the Agilent Bioanalyzer 2100 system (Agilent, CA, USA). The library preparations were sequenced on an Illumina Novaseq platform, and 150-bp paired-end reads were generated; the clean data (clean reads) were obtained by removing reads containing adapters and poly-N and low-quality reads from the raw data. The index of the reference genome was built and the paired-end clean reads were aligned to the reference genome using Hisat2 v2.0.5 [40]. The gene expression levels were quantified with featureCounts v1.5.0-p3 and FPKM [41].

Genes with an adjusted P-value <0.05 and a |log₂FoldChange| > 1 according to DESeq2[42] were considered differentially expressed. Gene ontology (GO) enrichment analysis and Kyoto Encyclopedia of Genes and Genomes database (KEGG) pathway enrichment analysis of the differentially expressed genes were implemented by the cluster Profiler R package, and genes with P-values less than 0.05 were considered significantly enriched as differentially expressed genes [43].

Sample processing and liquid chromatography coupled with tandem mass spectrometry (LC-MS/MS)

The proteins in the tissues were extracted with lysis buffer containing 8 M urea (Sigma Aldrich, MO, USA) and 1% protease inhibitor cocktail (Merck Millipore, MA, USA), and the protein concentration was determined with a BCA kit (Beyotime, Shanghai, China) according to the manufacturer’s instructions. Protein enzymolysis was performed using trypsin (Promega, WI, USA). After trypsin digestion, the peptides were desalted on a Strata X C18 SPE column (Phenomenex, CA, USA) and vacuum-dried. The peptides were reconstituted in 0.5 M TEAB and processed according to the manufacturer’s protocol for the TMT kit (Thermo Fischer Scientific, MA, USA). The tryptic peptides were fractionated by high-pH reverse-phase HPLC using an Agilent 300 Extend C 18 column (5 μm particle size, 4.6 mm ID, 250 mm length). The tryptic peptides were dissolved in 0.1% formic acid (Sigma Aldrich, MO, USA) (solvent A) and directly loaded onto an in-house-packed reversed-phase analytical column (150 mm length, 75 μm i.d.). The gradient included the following steps: from 6% to 23% solvent B (0.1% formic acid in 98% acetonitrile) over 26 min, from 23% to 35% solvent B over 8 min and from 35% to 80% over 3 min; then, the column was eluted with 80% solvent B for the last 3 min of the run. The flow rate was constant at 400 nL/min, and an EASY-nLC 1000 UPLC system was used. The peptides were exposed to an NSI source followed by tandem mass spectrometry (MS/MS) with a Q Exactive^TM Plus spectrometer (Thermo Fischer Scientific, MA, USA) coupled online to the UPLC system. The electrospray voltage was 2.0 kV. The m/z scan range was 350 to 1,800 for the full scan, and the intact peptides were detected using an Orbitrap at a resolution of 70,000. Peptides were then selected for MS/MS using an NCE set at 28, and the fragments were detected using an Orbitrap at a resolution of 17,500. A data-dependent procedure alternated between an MS scan and 20 MS/MS scans with 15.0 s dynamic exclusion. Automatic gain control (AGC) was set at 5E4. The fixed first mass was set as 100 m/z.

The resulting MS/MS data were processed using the Maxquant search engine (v.1.5.2.8). The tandem mass spectra were searched against the Gallus gallus database (http://www.uniprot.org/proteomes/UP000000539) concatenated with the reverse decoy database. Trypsin/P was specified as the cleavage enzyme, and up to 2 missing cleavages were allowed. The mass tolerance for the precursor ions was set to 20 ppm for the first search and to 5 ppm for the main search, and the mass tolerance for the fragment ions was set to 0.02 Da. Carbamidomethyl was specified as a fixed modification of Cys, and the oxidation of Met was specified as a variable modification. The FDR was < 1%, and the minimum score for the peptides was > 40.

The GO proteome was derived from the UniProt-GOA database (www.http://www.ebi. ac.uk/GOA/). The Kyoto Encyclopedia of Genes and Genomes (KEGG) database was used to identify the enriched pathways. Proteins with a threshold of P<0.05 and a fold change of > 1.5 or <1/1.5 were identified as differentially expressed proteins (DEPs) between SY and F6 follicles.

Parallel reaction monitoring

Protein extraction and trypsin digestion were performed as described above. Then, the tryptic peptides were dissolved in 0.1% formic acid (solvent A) and directly loaded onto an in-house-packed reversed-phase analytical column (150 mm length, 75 μm i.d.). The gradient included the following steps: from 6% to 23% solvent B (0.1% formic acid in 98% acetonitrile) over 38 min, from 23% to 35% solvent B over 14 min, from 35% to 80% solvent B over 4 min, and 80% solvent B for 4 min; a constant flow rate of 400 nL/min was maintained using an EASY-nLC 1000 UPLC system. The peptides were exposed to an NSI source followed by tandem mass spectrometry (MS/MS) using a Q ExactiveTM Plus spectrometer (Thermo Fischer Scientific, MA, USA) connected to the UPLC system. The electrospray voltage was 2.0 kV. The m/z scan range was 350 to 1,000 for the full scan, and intact peptides were detected using an Orbitrap at a resolution of 35,000. Then, the peptides were selected for MS/MS using an NCE set at 27, and the fragments were detected using an Orbitrap at a resolution of 17,500. The data-independent procedure alternated between an MS scan followed by 20 MS/MS scans. The automatic gain control (AGC) was set at 3E6 for the full MS and at 1E5 for MS/MS. The maximal IT was set at 20 ms for the full MS and was switched to automatic mode for MS/MS. The isolation window for MS/MS was set to 2.0 m/z.

The resulting MS/MS data were processed using Skyline (v.3.6). The following peptide settings were used: the enzyme was set as trypsin [KR/P]; the maximal number of missed cleavages was set as 2; the peptide length was set as 8-25; variable modifications were set as carbamidomethylation of Cys and oxidation of Met; and the maximal variable modification was set as 3. The following transition settings were used: The precursor charges were set as 2, 3; the ion charges were set as 1, 2; and the 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.

Cell culture and cell treatment

All of the small yellow follicles and the hierarchical follicles were isolated from hens and placed in PBS. The small yellow follicles were treated with 0.1% collagenase II (MP Biomedicals, Santa Ana, CA, USA) at 38°C for 15 min to obtain the GCs. The yolk of the follicles was removed carefully with ophthalmic forceps, and the GCs were isolated from the hierarchical follicles and then dispersed by treatment with 0.25% trypsin-EDTA (Gibco-BRL, NY, USA) at 38°C for 15 min with gentle oscillation in a centrifuge tube. After centrifugation, the GCs were suspended in culture medium that contained M199 (Gibco-BRL, NY, USA), 10% fetal bovine serum (Biological Industries, Kibbutz Beit Haemek, Israel) and 1% penicillin/streptomycin (Solarbio, Beijing, China) and subsequently seeded in 24-well culture plates at a density of 1 × 10⁵ cells/well. The number of viable cells (90%) was estimated using Trypan blue staining. The cells were cultured at 38°C in a water-saturated atmosphere of 95% air and 5% CO₂. After 24 h, the cultured cells were subsequently treated with different concentrations (0, 5, 10, and 50 ng/ml) of recombinant FSH (R&D Systems, MN, USA). All the treated cells were collected after another 24 h for RNA extraction and qRT-PCR analysis.

Real-time quantitative PCR

Total RNA was extracted from the S and F ovarian follicles that were used for proteome analysis using TRIzol reagent (InvivoGen, CA, USA). Synthesis of the cDNA was performed using 1 μg of the RNA by using a PrimeScript RT reagent kit with gDNA Eraser (TaKaRa, Dalian, China) according to the manufacturer’s protocol. Real-time quantitative PCR was performed using a SYBR Premix Ex Taq^TM II kit (TaKaRa, Dalian, China) on a Light Cycler 480 real-time PCR system (Roche, Basel, Switzerland) as follows: 95°C for 30 s, followed by 40 cycles of denaturation at 95°C for 10 s and annealing and extension at 58°C for 20 s. The melting curves were obtained, and quantitative analysis of the data was performed using the 2^−ΔΔCT relative quantification method [44]. The mRNAs of VLDLR, VLDLR1, WIF1, NGFR, AMH, BMP15, GDF6 and MMP13 were quantitatively analyzed by this method. Quantification was performed by standardizing the reaction results versus those of β-actin. All primers are listed in Table S5.

Statistical analysis

Differences in mRNA expression between follicles 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 in SAS (version 9.2). For one experiment, each treatment was repeated four times, and at least three independent experiments were performed. All data are presented as the mean ± SEM (n = 4). The GO and KEGG analyses were performed by using Fisher’s t-test. P<0.05 was considered to indicate a significant difference.

Abbreviations

BP: Biological process; CC: Cellular component; DEGs: Differentially expressed genes; DEPs: Differentially expressed proteins; GO: Gene ontology; KEGG: Kyoto Encyclopedia of Genes and Genomes; GCs: granulosa cells; TCs: thecal cells.

Declarations

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.

Availability of data and materials

For transcriptome sequencing, the raw reads were submitted to the Sequence Read Archive (SRA, https://www.ncbi.nlm.nih.gov/sra/docs/) repository, https://www. ncbi.nlm.nih.gov/sra/?term=SRP236909. The MS proteomics data have been deposited in the ProteomeXchange Consortium (http://www.proteomex change.org/) via the PRIDE partner repository with the dataset identifier PXD011470 (Username: [email protected], Password: AfwLJZRb).

Competing interests

The authors declare that they have no competing interests.

Funding

This research was financially supported by grants from the National Natural Science Foundation of China (NSFC 31672414, 31272435, 31772588), the Funds of Shandong“Double Tops” Program (SYL2017YSTD12) and Shandong Agricultural Breed Project (2019LZGC019). The funding body had no role in the design of the study and collection, analysis, and interpretation of data and in writing the manuscript.

Authors’ contributions

QC, LK and YJ designed and drafted the manuscript. QC, YW and ZL carried out animal care, prepared samples and performed the experiments. QC and XG 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.

Acknowledgments

The authors thank PTM-Biolabs Co., Ltd. (Hangzhou, China) for mass spectrometry analysis and help in preparing the manuscript. The authors would also like to thank Novogene Bioinformatics Technology Co., Ltd. (Beijing, china) for their bioinformatics supports.

References

1. Johnson AL. Ovarian follicle selection and granulosa cell differentiation. Poult Sci, 2014;94(4):781–785.
2. Onagbesan O, Bruggeman V, Decuypere E. Intra-ovarian growth factors regulating ovarian function in avian species: a review. Anim Reprod Sci. 2009;111:121–140.
3. Tilly JL, Kowalski KI and Johnson AL. Stage of ovarian follicular development associated with the initiation of steroidogenic competence in avian granulosa cells. Biol Reprod, 1991;44(2):305–314.
4. Johnson AL, Woods DC. Ovarian dynamics and follicle development. In: Jamieson, B.G.M. (Ed.), Reproductive Biology and Phylogeny of Aves. 2007;Science Publishers, Inc., pp: 243–277 (Chapter 6).
5. Stifani S, Barber DL, Nimpf J, Schneider W. A single chicken oocyte plasma membrane protein mediates uptake of very low density lipoprotein and vitellogenin. Proc Natl Acad Sci. 1990;87:1955–1959.
6. Schneider W,. Receptor-mediated mechanisms in ovarian follicle and oocyte development. Gen Comp Endocrinol. 2009;163:18–23.
7. Diaz FJ, Anthony K, Halfhill AN. Early avian follicular development is characterized by changes in transcripts involved in steroidogenesis, paracrine signaling and transcription. Mol Reprod Dev. 2011;78:212–223.
8. Zhu G, Yong M, Zhou W, Jiang Y. Dynamic changes in the follicular transcriptome and promoter DNA methylation pattern of steroidogenic genes in chicken follicles throughout the ovulation cycle. PLoS One. 2015;10:e0146028.
9. Wang Y, Chen Q, Liu Z, Guo X, Du Y, Yuan Z, Guo M, Sun Y, Kang L, Jiang Y. Transcriptome analysis on single small yellow follicles reveals that WNT4 is involved in chicken follicle selection. Front Endocrinol. 2017;8:317.
10. Ocón-Grove OM, Poole DH, Johnson AL. Bone morphogenetic protein 6 promotes FSH receptor and anti-müllerian hormone mRNA expression in granulosa cells from hen prehierarchal follicles. Reproduction. 2012;14:3825.
11. Kim D, Ocón-Grove O, Johnson AL. Bone morphogenetic protein 4 supports the initial differentiation of hen (Gallus gallus) granulosa cells. Biol Reprod. 2013;88: 161.
12. Stephens CS, Johnson PA. Bone morphogenetic protein 15 may promote follicle selection in the hen. Gen Comp Endocrinol. 2016;235:170–176.
13. Knight PG, Glister C. TGF-beta superfamily members and ovarian follicle development. Reproduction. 2016;132:191–206.
14. Zhu G, Chen X, Mao Y, Kang L, Ma X, Jiang Y. Characterization of annexin A2 in chicken follicle development: Evidence for its involvement in angiogenesis. Anim Reprod Sci. 2015;161:104–111.
15. Ewen K, Baker M, Wilhelm D, Aitken RJ, Koopman P. Global survey of protein expression during gonadal sex determination in mice. Mol Cell Proteomics. 2009;8:2624.
16. Ma M, Guo X, Wang F, Zhao C, Liu Z, Shi Z, Wang Y, Zhang P, Zhang K, Wang N. et al. Protein expression profile of the mouse metaphase-ii oocyte. J Proteome Res. 2008;7:4821–4830.
17. Wang L, Zhu YF, Guo XJ, Huo R, Ma X, Lin M, Zhou ZM, Sha JH. A two-dimensional electrophoresis reference map of human ovary. J Mol Med. 2005;83:812–821.
18. Pfeiffer MJ, Siatkowski M, Paudel Y, Balbach ST, Baeumer N, Crosetto N, Drexler HC, Fuellen G, Boiani M. Proteomic analysis of mouse oocytes reveals 28 candidate factors of the "reprogrammome". J Proteome Res. 2011;10:2140–2153.
19. Upadhyay RD, Balasinor NH, Kumar AV, Sachdeva G, Parte P, Dumasia K. Proteomics in reproductive biology: beacon for unraveling the molecular complexities. Biochim Biophys Acta. 2013;1834:8–15.
20. Massicotte L, Coenen K, Mourot M, Sirard MA. Maternal housekeeping proteins translated during bovine oocyte maturation and early embryo development. Proteomics. 2006;6:3811–3820.
21. Memili E, Peddinti D, Shack LA, Nanduri B, McCarthy F, Sagirkaya H, Burgess SC. Bovine germinal vesicle oocyte and cumulus cell proteomics. Reproduction. 2007;133: 1107–1120.
22. Peddinti D, Memili E, Burgess SC. Proteomics-based systems biology modeling of bovine germinal vesicle stage oocyte and cumulus cell interaction. PLoS One. 2010;5:e11240.
23. Nepomuceno AI, Muddiman DC, Petitte JN. Global proteomic analysis of functional compartments in immature avian follicles using laser microdissection coupled to LC-MS/MS. J Proteome Res. 2015;14:3912–3923.
24. Johnson PA. Follicle selection in the avian ovary. Reprod Domest Anim. 2012;47 Suppl 4:283–287.
25. Kim D, Johnson AL. Vasoactive intestinal peptide promotes differentiation and clock gene expression in granulosa cells from prehierarchal follicles. Mol Reprod Dev. 2016;83(5):455–463.
26. Guo X, Wang Y, Chen Q, Yuan Z, Chen Y, Guo M, Kang L, Sun Y, Jiang Y. The role of PTHLH in ovarian follicle selection, its transcriptional regulation and genetic effects on egg laying traits in hens. Front Genet. 2019;10:430.
27. Fan Y, Zhang C, Zhu G. Profiling of RNA N6-methyladenosine methylation during follicle selection in chicken ovary. Poult Sci. 2019;98(11):6117–6124.
28. Nimpf J, Radosavljevic MJ, Schneider WJ. Oocytes from the mutant restricted ovulator hen lack receptor for very low density lipoprotein. J Biol Chem. 1989;264:1393–1398.
29. Hayashi K, Nimpf J, Schneider WJ. Chicken oocytes and fibroblasts express different apolipoprotein-B-specific receptors. J Biol Chem. 1989;264:3131–3139.
30. Eresheim C, Leeb C, Buchegger P, Nimpf J. Signaling by the extracellular matrix protein reelin promotes granulosa cell proliferation in the chicken follicle. J Biol Chem. 2014;289:10182–10191.
31. Lin X, Ma Y, Qian T, Yao J, Mi Y, Zhang C.Basic fibroblast growth factor promotes prehierarchical follicle growth and yolk deposition in the chicken. Theriogenology. 2019;139:90–97.
32. Hu S, Liu H, Pan Z, Xia L, Dong X, Li L, Xu F, He H, Wang J. Molecular cloning, expression profile and transcriptional modulation of two splice variants of very low density lipoprotein receptor during ovarian follicle development in geese (Anser cygnoide). Anim Reprod Sci. 2014;149:281–296.
33. .Buyuk E, Santoro N, Cohen HW, Charron MJ, Jindal S. Reduced neurotrophin receptor tropomyosin-related kinase A expression in human granulosa cells: a novel marker of diminishing ovarian reserve. Fertil Steril. 2011;96(2):474-478.
34. Mann K. Proteomic analysis of the chicken egg vitelline membrane. Proteomics. 2008;8:2322–2332.
35. Mann K, Mann M. The chicken egg yolk plasma and granule proteomes. Proteomics. 2010;8:178–191.
36. Wassarman PM, Jovine L, Qi H, Williams Z, Darie C, Litscher ES. Recent aspects of mammalian fertilization research. Mol Cell Endocrinol. 2005;234:95–103.
37. Sakaguchi Y, Uzuhashi R, Iwata H, Monji Y, Kuwayama T. Changes in the sperm-zona pellucida binding properties during porcine oocyte maturation. Journal of Mammalian Ova Research. 2010;27:130–135.
38. Song J, Wang C. Transcriptomic and proteomic analyses of genetic factors influencing adductor muscle coloration in QN Orange scallops. BMC Genomics. 2019;20:363.
39. Abdulghani M, Song G, Kaur H, Walley JW, Tuteja G. (2019). Comparative analysis of the transcriptome and proteome during mouse placental development. J Proteome Res. 2019.
40. Kim D, Langmead B, Salzberg SL. HISAT: a fast spliced aligner with low memory requirements. Nat Methods. 2015;12(4):357–360.
41. Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28(5):511–515.
42. Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010;11:R16.
43. Young MD, Wakefield MJ, Smyth GK, Oshlack A. Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biol. 2010;11:R14.
44. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. Methods. 2001;25:402–408.