Abstract
For two libraries (PDF1 and ODF1) using Illumina sequencing 44,082,301 and 43,708,132 clean reads were obtained, respectively. After being mapped to the bovine RefSeq database, 15,533 genes were identified to be expressed in both types of follicles (cut-off RPKM > 0.5), of which 719 were highly expressed in bovine follicles (cut-off RPKM > 100). Furthermore, 83 genes were identified as being differentially expressed in ODF1 versus PDF1, where 42 genes were upregulated and 41 genes were downregulated. KEGG pathway analysis revealed two upregulated genes in ODF1 versus PDF1, CYP11A1, and CYP19A1, which are important genes in the steroid hormone biosynthesis pathway. This study represents the first investigation of transcriptome of bovine follicles at predeviation and onset of deviation stages and provides a foundation for future investigation of the regulatory mechanisms involved in follicular development in cattle.
1. Introduction
The ovarian follicle is an essential component of the reproductive process. It plays an important role in controlling the estrous cycle, determining estrous behaviour, ensuring oocyte competency and subsequent embryo survival rate, and determining both postovulation corpus luteum function and progesterone synthesis [1]. In a number of species, follicular growth is characterized by a wave-like pattern, with two or three waves occurring during the normal course of estrous cycles in cattle [2]. During each wave of follicular development, a cohort of antral follicles are induced to begin accelerated growth [3]. After a period of concurrent growth, a species specific number of follicles will then be selected to become dominant, while the remaining follicles will be lost through a process known as atresia. Diameter deviation is defined as the divergence in growth rates between the two largest follicles in a follicular wave [3]. The onset of diameter deviation occurs when the largest follicle reaches 8.5 mm in dairy cattle and marks initiation of divergence in growth rate and estradiol producing capacity between the F1 or largest (future dominant follicle) and F2 or second largest (future subordinate) growing follicles culminating the process of dominant follicle selection. While the exact mechanisms of dominant follicle selection are not completely understood, there have been many studies on the hormones and factors involved in follicular development. Antral follicles are dependent upon FSH for growth and each follicular wave is preceded by a transient rise in FSH concentrations [4]. Many growth factors linked to regulation of follicular development, such as inhibins, activins, and insulin-like growth factors and their binding proteins, have been identified in follicular fluid of individual bovine follicles [5]. These molecules can regulate follicular cell survival, proliferation, or death. Recent studies have attempted to understand the molecular regulation of follicular development in cattle [6, 7]. However, the molecular mechanisms governing the wave-like pattern of follicular development are incompletely described, particularly at the onset of diameter deviation which is the first morphological indication of follicular dominance.
Traditionally, gene expression studies in the field of follicular development focus on the study of expression of candidate genes of interest. With the development of next-generation sequencing technologies, transcriptome profiling has become a powerful approach for identification of genes globally expressed in various tissues including ovarian follicles [8]. In the present study, we performed RNA-Seq of granulosa cell RNA from bovine ovarian follicles at predeviation (PD) and onset of deviation (OD) stages of a follicular wave in cattle to catalog the transcriptome and identify potential differentially expressed genes associated with these key stages of a follicular wave. This study provides a comprehensive sequence resource for future studies on follicular development in cattle.
2. Materials and Methods
2.1. Materials
All materials were obtained from Sigma-Aldrich (St. Louis, MO) unless otherwise stated.
2.2. Animal Model and Sample Collection
All animal procedures were approved by the Institutional Animal Care and Use Committee at Michigan State University. Estrus was synchronized in nonlactating Holstein dairy cows with two injections of prostaglandin F2α (PGF2α; Prostamate; IVX Animal Health, St. Joseph, MO) administered 14 days apart, and follicular growth was observed and recorded by daily ultrasonography.
Ovaries were removed from cows at the following stages of the first follicular wave: predeviation (PD; approximately Day 3 after estrus; 1.5 days after emergence [emergence is the first scan where a new follicle at least 4 mm is detected by ovarian ultrasonography]) and onset of deviation (OD; first scan where growth of the F1 [largest; future dominant] follicle to >8.5 mm was detected by ovarian ultrasonography). The F1 follicles were isolated from the PD and OD groups. Granulosa cells were isolated from the two types of follicles (PDF1 and ODF1), lysed, and stored at −80°C immediately.
2.3. RNA Isolation
Total RNA was isolated from the lysed granulosa cells using the RNeasy mini kit (Qiagen) and DNase treated on column according to the manufacturer’s protocol. The RNA integrity was evaluated by Agilient Bioanalyser and the RNA concentration was measured using a Nanodrop-1000 spectrophotometer. RNA samples with a RNA integrity number greater than 8 were selected for deep sequencing.
2.4. Library Preparation and Illumina Sequencing
RNA sequencing was performed by the WM Keck Center for Comparative and Functional Genomics at the University of Illinois at Urbana-Champaign. RNA samples from four follicles were pooled within each group (ODF1 or PDF1). RNA-Seq libraries were prepared with a TruSeq RNA Sample Preparation kit (Illumina) according to the manufacturer’s instructions. The cDNA libraries were sequenced on one lane for 100 cycles using Illumina HiSeq™ 2000 by a TruSeq SBS kit v5 (Illumina) and analyzed with pipeline version 1.8.
2.5. Identification of Differentially Expressed Genes and Pathway Analysis
The CLC Genomics Workbench (CLC bio, Aarhus, Denmark) was used to map the sequence reads to the bovine RefSeq database. The reads per kilobase per million reads (RPKM) values were calculated as the normalized transcript expression values [9]. A -test [10] was used to identify differentially expressed genes between ODF1 and PDF1 (FDR corrected value , RPKM cut-off > 0.5, and RPKM fold change >1.5) using the CLC genomics workbench. DAVID software (https://david.ncifcrf.gov/gene2gene.jsp) was used to perform GO annotations and KEGG pathway analysis for highly expressed (RPKM > 100) and differentially expressed genes.
3. Results
3.1. Illumina Sequencing
To identify differentially expressed genes involved in bovine follicular development, Illumina sequencing was used on two libraries constructed from RNA isolated from ODF1 and PDF1 follicles. After filtering, a total of 44,082,301 and 43,708,132 clean reads were obtained from PDF1 and ODF1 libraries, respectively. The clean reads were mapped to the bovine RefSeq database (containing 35,325 annotated transcripts). Using a cut-off value of RPKM > 0.5, a total of 15,533 genes were identified in both types of follicles (Additional file 1: Table S1 in Supplementary Material available online at http://dx.doi.org/10.1155/2016/3472748), among which 719 are considered to be highly expressed (RPKM cut-off > 100) in bovine follicles (Additional file 2: Table S2).
3.2. GO Functional Classification and KEGG Pathway Analysis of Highly Expressed Genes
The top 30 highly expressed genes in granulosa cells of bovine follicles at the predeviation and onset of deviation stages are shown in Table 1. Many of them are known to be important for follicular growth and development, such as Serpin peptidase inhibitor clade E member 2 (SERPINE2), Inhibin alpha (INHA), Inhibin beta A (INHBA), and Follistatin (FST). GO functional classification of these highly expressed genes was performed using DAVID software. All 719 highly expressed genes can be assigned into 22 groups under three categories (biological process, 39%; cellular component, 44%; and molecular function, 17%) based on their putative functions (Figure 1). Many of the highly expressed genes are involved in metabolic process, multicellular organismal process, and binding. KEGG pathway analysis showed that the highly expressed genes are involved in 12 major pathways (Figure 2), of which the most significantly enriched genes are involved in ribosome pathway.
3.3. Differentially Expressed Genes in ODF1 versus PDF1
Using RPKM cut-off > 0.5 and fold change cut-off > 1.5 at FDR corrected value , a total of 83 differentially expressed genes were identified, with 41 downregulated genes (Table 2) and 42 upregulated genes (Table 3) in ODF1 versus PDF1. To understand the functions of these differentially expressed genes, GO analysis was performed. The upregulated genes were categorized into 14 functional groups under 3 major GO classifications: biological process (35%), cellular component (30%), and molecular function (35%) (Table 4). Many of the differentially expressed genes are known to play a role in ovarian follicular development (Table 5). For example, serine protease 23 (PRSS23) is expressed in granulosa cells and may play a crucial role in follicular atresia, whereas serine protease 35 (PRSS35) is also expressed in granulosa cells and may be involved in ovulation and CL formation and regression. KEGG pathway analysis of the upregulated genes demonstrated that two important genes (CYP11A1 and CYP19A1) in the steroid hormone biosynthesis pathway are upregulated in ODF1 versus PDF1.
4. Discussion
Follicular growth occurs in a characteristic wave-like pattern in monotocous species such as cattle [3, 5, 11]. A transient increase in FSH triggers initiation of each follicular wave [5, 11, 12]. Emergence is defined as the first day a new follicle >4 mm in diameter is detected and is the first chronological event marking a new follicular wave that is detectable by ultrasonography. After emergence, follicles in the cohort initially grow at a similar rate (common growth phase) prior to deviation [3]. However, the molecular mechanisms involved regulating the onset of deviation are not well understood, in order to characterize the differences in gene expression that associated with follicular development in different follicles sized in diameter, which the previous studies examined using microarray technology [13–17]. To further investigate the bovine granulosa cell transcriptome and molecular alterations associated with onset of deviation, we examined the transcriptome at specific stages of the estrous cycle.
Illumina sequencing technology was used to determine gene expression levels in ODF1 and PDF1 follicles. A total of 15,533 genes were identified in both types of follicles and 83 of them were identified as differentially expressed between ODF1 and PDF1. Our study provided novel information on the bovine granulosa cell transcriptome and identified specific transcripts highly expressed in granulosa cells of bovine follicles prior to and at onset of deviation, including transcripts encoding for several housekeeping genes (e.g., ribosomal proteins L18a, S27a, and L4) and genes with well-established roles in regulation of ovarian function (e.g., INHBA, INHBB, and FST). Of particular interest was SERPINE2, which is abundantly expressed in granulosa cells of follicles collected at both the predeviation and onset of deviation stages of a follicular wave, illustrating its potential importance in bovine ovarian follicular development. Estradiol and SERPINE2 secretion are positively correlated, but estradiol treatment cannot alter the expression of SERPINE2. FSH and growth factors can directly regulate the expression and secretion of SERPINE2 in granulosa cells, and SERPINE2 is an antiapoptotic factor, which may regulate atresia in bovine follicles [18]. Eleven SERPINE genes are expressed in bovine follicles, but only SERPINE2, SERPINE5, and SERPINE6 are expressed in the granulosa cells [19].
KEGG analysis revealed upregulated genes associated with onset of deviation (CYP11A1 and CYP19A1) involved in the steroid hormone biosynthesis pathways that play an essential role during follicular development. Proteins encoded by CYP11A1 and CYP19A1 genes are members of the cytochrome P450 superfamily, which are monooxygenases that catalyze many reactions involved in steroidogenesis. Previous studies suggested that CYP19A1 was regulated by multiple pathways, including estrogen receptors and cAMP/protein kinase A which are activated by FSHR in granulosa cells, and these regulatory mechanisms are likely critical for acquisition of follicular dominance in cattle [20]. Our transcriptome sequencing data is consistent with these results. The increase in transcript abundance for CYP19A1 in ODF1 versus PDF1 follicles is consistent with the increase in estradiol producing capacity associated with diameter deviation [3].
It is acknowledged that study design was not optimal due to limited biological replication because single pooled samples ( per group) were used in Illumina sequencing analysis. Despite such limitations, results have significantly enhanced understanding of bovine follicle transcriptome composition and potential differences in gene expression associated with follicular development that are foundational to further study in the future; several interesting candidates were revealed for future investigation, particularly genes linked to regulation of cell proliferation and survival. For example, results of present studies suggest that PPM1K, a Mn2+/Mg2+-dependent protein phosphatase of PPM family, is potentially upregulated with the onset of deviation in the granulosa cell layer of bovine follicles. This protein is critical for cell survival and embryonic development and can regulate the mitochondrial membrane permeability transition pore opening [21]. Potential upregulation of granulosa cell BEX2 and GREB1 transcript abundance was also noted in association with onset of dominance and may be associated with enhanced granulosa cell survival. BEX2 can downregulate apoptosis and activate the JNK (Jun NH2-terminal kinase) pathway, and these effects can be abolished by administration of a JNK specific inhibitor [22]. GREB1 is an estrogen receptor and coactivator linked to cell proliferation and GREB1 expression is estrogen dependent. It is possible that increased expression of PPM1K, BEX2, and GREB1 may be associated with granulosa cell proliferation and survival during the onset of deviation [23]. TNFAIP6/TSG6 is tumor necrosis factor and alpha-induced protein 6; it is suggested that TSG-6 plays a role in cell-cell or cell-cell matrix interactions during inflammation and tumorigenesis. High LH/hCGR gene expression intensity was associated with TNFAIP6/TSG6 gene expression which has a pivotal importance in the mucification of the COC during the preovulatory period [24]. It is suggested that the expression levels of TNFAIP6/TSG6 were nearly 280-fold in granulosa cells of large follicles than that of small follicles [14]. In our study, TNFAIP6/TSG6 was also differentially expressed in ODF1 and PDF1 with a 4.26-fold change. Altogether, these characteristics suggest that TNFAIP6/TSG6 plays a crucial role in accelerating follicle growth during follicular waves in cattle.
5. Conclusions
The present study characterized the granulosa cell transcriptome of bovine follicles at specific stages of follicular development and identified 83 differentially expressed genes potentially associated with onset of deviation, many of which are linked to regulation of follicular development. The study provides a foundation for future studies to investigate regulation of granulosa cell expressed genes and the regulatory mechanisms controlling antral follicle development during follicular waves in cattle.
Competing Interests
The authors declare that they have no competing interests.
Acknowledgments
This study was supported by the Chinese Natural Science Foundation (Grant no. 31172211), Introduction of Advanced International Agricultural Science and Technology Plan Project Fund (Grant no. 2010-Z43) to Lihua Lyu, Shanxi Province Science and Technology Research Project (Grant no. 20130311027-2) and Shanxi Agriculture University Introduction of Doctor Scientific Research Startup Fund (Grant no. 2014ZZ04) to Pengfei Li. Agriculture, and Food Research Initiative Competitive Grant no. 2009-65203-05700 from the USDA National Institute of Food and Agriculture to George W. Smith and by Michigan AgBioResearch (George W. Smith). The authors thank Liying Qiao and Wenge Hu, technicians of the College of Animal Science and Technology, Shanxi Agricultural University, for their technical assistance.
Supplementary Materials
Table S1: Shows all the genes (15,533) expressed in both types of follicles (ODF1 and ODF2) with a cut-off RPKM of 0.5.
Table S2: Shows the highly expressed genes (719) in PDF1 and ODF1 follicles with a cut-off RPKM > 100.