# An overview of gene expression dynamics during early ovarian folliculogenesis: specificity of follicular compartments and bi-directional dialog

Agnès Bonnet, Cédric Cabau, Olivier Bouchez, Julien Sarry, Nathalie Marsaud, Sylvain Foissac, Florent Woloszyn, Philippe Mulsant, Béatrice Mandon-Pepin

<figure><img src="image_1.png" /></figure>

## ► To cite this version:

Agnès Bonnet, Cédric Cabau, Olivier Bouchez, Julien Sarry, Nathalie Marsaud, et al.. An overview of gene expression dynamics during early ovarian folliculogenesis: specificity of follicular compartments and bi-directional dialog. BMC Genomics, 2013, 14 (december), pp.1-19. 10.1186/1471-2164-14-904. hal-01019557

HAL Id: hal-01019557

https://hal.science/hal-01019557

Submitted on 23 May 2019

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>
# An overview of gene expression dynamics during early ovarian folliculogenesis: specificity of follicular compartments and bi-directional dialog

Agnes Bonnet¹,<fnref n="2" />*, Cedric Cabau³, Olivier Bouchez¹,<fnref n="2" />,⁴, Julien Sarry¹,<fnref n="2" />, Nathalie Marsaud⁴,⁵,⁶,⁷, Sylvain Foissac¹,<fnref n="2" />, Florent Woloszyn¹,<fnref n="2" />, Philippe Mulsant¹,<fnref n="2" /> and Beatrice Mandon-Pepin⁸

## Abstract

**Background:** Successful early folliculogenesis is crucial for female reproductive function. It requires appropriate gene specific expression of the different types of ovarian cells at different developmental stages. To date, most gene expression studies on the ovary were conducted in rodents and did not distinguish the type of cell. In mono-ovulating species, few studies have addressed gene expression profiles and mainly concerned human oocytes.

**Results:** We used a laser capture microdissection method combined with RNA-seq technology to explore the transcriptome in oocytes and granulosa cells (GCs) during development of the sheep ovarian follicle. We first documented the expression profile of 15 349 genes, then focused on the 5 129 genes showing differential expression between oocytes and GCs. Enriched functional categories such as oocyte meiotic arrest and GC steroid synthesis reflect two distinct cell fates. We identified the implication of GC signal transduction pathways such as SHH, WNT and RHO GTPase. In addition, signaling pathways (VEGF, NOTCH, IGF1, etc.) and GC transzonal projections suggest the existence of complex cell-cell interactions. Finally, we highlighted several transcription regulators and specifically expressed genes that likely play an important role in early folliculogenesis.

**Conclusions:** To our knowledge, this is the first comprehensive exploration of transcriptomes derived from *in vivo* oocytes and GCs at *key* stages in early follicular development in sheep. Collectively, our data advance our understanding of early folliculogenesis in mono-ovulating species and will be a valuable resource for unraveling human ovarian dysfunction such as premature ovarian failure (POF).

**Keywords:** RNA-seq, Early folliculogenesis, Molecular dialog, Transcriptome, Sheep, Oocyte, Granulosa cells, Microdissection

## Background

The major function of an ovary is to produce oocytes intended for fertilization as eggs and to create a favorable environment for the beginning of gestation. In adult mammals, the ovary is a heterogeneous organ containing follicles at various stages of development and *corpora lutea*, either active or at various stages of involution. The formation of primordial follicles (oocytes surrounded by flattened pre-granulosa cells) occurs during

fetal development in many species including human, sheep (from 75 days of gestation) [1], cattle, and goat, or after birth like in rodents. The primordial follicles represent a reserve of germ cells for the entire reproductive life of the female. They remain dormant until their recruitment and irreversible growth towards the primary, secondary and tertiary stages (with an antral cavity). The gradual exit of primordial follicles begins shortly after the formation of the primordial follicle pool and continues throughout the reproductive years [2]. Consequently, this early development is important as it regulates the size of the resting primordial follicle pool and the fate of the follicles, which, in turn, affects fertility and the reproductive life span.

<footnote>* Correspondence: agnes.bonnet@toulouse.inra.fr
¹INRA, UMR444 Génétique Cellulaire, Auzeville, BP52627, F-31326 Castanet-Tolosan, France</footnote>

<footnote n="2">²ENVT, UMR444, Laboratoire de Génétique Cellulaire, F-31326 Castanet-Tolosan, France</footnote>

<footnote>Full list of author information is available at the end of the article</footnote>
Early follicular development requires the appropriate expression of many genes at different developmental stages. Recent studies on natural mutations in sheep [3] and on mutant mouse models demonstrated that the expression of different oocyte-specific genes is essential during early folliculogenesis in a stage specific expression pattern [2]. Figla, Nobox, Pten, Foxo3A and nerve growth factor genes are involved in the formation of primordial follicles, whereas genes such as *βFGF*, *GDF 9* and *BMP 4* are involved in the transition from primordial to primary follicles. Other genes such as *BMP15* are not expressed in the oocyte until the primary follicle stage and are involved in the transition from primary to secondary follicles, as shown in sheep (for a review see [4]).

This process also requires orchestrated communication between oocytes and granulosa cells (GCs), which are the structural components of the early follicle. Oocytes and GCs regulate follicle growth in an autocrine and paracrine manner via secreted factors and direct gap junctional communications. On one hand, it is assumed that oocytes control the proliferation of GCs and later, their differentiation into steroid secreting cells and their metabolism [5]. On the other hand, GCs are indispensable to oocyte growth, meiosis, cytoplasmic maturation and control of transcriptional activity within the oocyte [6]. This dialog also takes place through GC cytoplasmic extensions that penetrate the *zona pellucida* (ZP) and form gap junctions with oocyte cell membrane [7].

Although some large scale expression studies have been conducted on rodent and human oocytes [8-11], expression profiling of follicular compartments is difficult to accomplish due to the very small size of the follicles and the mixture of preantral stages in the ovary. Two studies in particular identified the specific transcriptome of human oocytes in the quiescent state [10,11]. These studies highlighted enriched functions (transcription, RNA post transcriptional modifications) and molecular mechanisms (cell cycle, androgen signaling) in human oocytes, which were also identified in the only preliminary analysis conducted in sheep species to date [12]. But exactly which granulosa cell factors play a specialized role in early follicular development is not known. Our knowledge of the key mechanisms and the bidirectional communication underlying early folliculogenesis is consequently still very limited and is mainly derived from rodents (poly-ovulating species). Spatial and temporal information on gene expression patterns would facilitate the identification of the regulatory gene networks that underlie cell growth, differentiation, and the functional specification of each follicular compartment.

The purpose of this study was thus to 1- describe global gene expression during early ovarian folliculogenesis

in each follicular compartment (oocyte, GCs) in mono-ovulating species using sheep as the experimental model (whose fertility, litter size, ovulation rate and fetus development is similar to humans, unlike rodents); 2- identify differential and specific gene expression between these two follicular compartments; 3- investigate specific functions and pathways, and 4- explore bi-directional communication between oocytes and GCs. This was achieved using laser capture microdissection (LCM) to isolate oocytes and GCs during the earliest stages of folliculogenesis. Gene expression profiles were monitored at each follicle stage according to cell type using RNA-seq technology. We report the successful isolation and generation of global gene expression data sets of pure oocytes and GC populations at key stages of early follicular development. We describe the expression profile of 15 349 genes with 33.4% of genes displaying differential expression between the two compartments. In addition, we identified oocyte and GC-specific gene expression and the pathways involved in oocyte-GC communication.

## Results

We combined laser capture microdissection (LCM) [12] and RNA-sequencing technologies to detect separate GCs and oocyte global transcriptome for each stage of follicle development: primordial (PD), primary (PM), secondary (SC) follicles and the small antral stage (SA) as control. Three/four biological replicates were obtained per condition (8 conditions: cell type X stage). For the experimental protocol, see Additional file 1: Figure S1.

The RNA-seq experiment generated around 2.647 billion 100 bp reads with an average of 73.7 million per LCM-derived amplified RNA sample (LCM-aRNA). The result of the bioinformatics and annotation processing is summarized in Figure 1 and detailed in Additional file 2. Briefly, the genome assembly approach (mapping to sheep genome sequence, assembly, followed by read counting) generated a collection of 382 933 fragments that aggregated 47.5% of LCM-aRNA reads. After 3' UTR selection, the final dataset included 15 349 genes. 221 716 genomic fragments remained unannotated.

## Exploring global gene expression

The 15 349 genes correspond to 62% of the bovine annotated genes in Ensembl UMD3.1 release 65, and to 78.6% of the ovine annotated genes in CSIRO Oarv2.0 (18 656 genes). As expected, the number of expressed genes was uniformly distributed between chromosomes according to the number of expected genes except for chromosome 15 on which significantly fewer genes were expressed (Additional file 1: Figure S2A). The percentage of genes transcribed for this chromosome was 57%, compared to an average of 80% for the other chromosomes
<figure><img src="image_3.png" /><figcaption>**Figure 1** Summary of bioinformatics processes. The RNA-seq experiment produced a collection of 382,933 fragments that aggregated 47.5% of the LCM-aRNA reads. The annotation strategy assigned 73% of the mapped reads. A total of 221,716 genomic fragments remain unannotated. A total of 86.8% of the annotated reads are located in stop codon or 3'UTR regions, whereas only 5.5% are located in exons, 1.2% in start codon/ 5'UTR regions and 6.5% in introns. The final dataset conserved a single fragment per gene corresponding to the nearest 3'UTR region with the highest number of reads. This dataset aggregated 89.4% of the annotated LCM-aRNA reads (86.8% were located in 3' UTR regions and 2.6% were located in exons).</figcaption></figure>

(SD = 7.9%). The chromosome coverage of the expressed genes overlaps Oarv2.0 annotated genes (Additional file 1: Figure S2B) except for chromosome 15 on which we identified two regions with fewer expressed genes. The density of expressed genes in these two regions (45–50 MB and 77–80 Mb) was 75%-90% lower.

For each sample, the expression levels were measured in read counts per gene (the 3' UTR fragment being the most representative) and a normalized expression (RPM = sample gene read count/(sample total mapped read /1 000 000)) was calculated to describe the transcriptome of the tissue. First, we visualized a dynamic range of gene expression between 0.2 and 1 000 RPM (Additional file 1: Figure S3A-B). We estimated that 95% of the genes aggregated 56% of the reads. We then examined the highest expressed genes and detected a clear distinction between cell types (Additional file 1: Figures S3C-D). The highest oocyte-expressed genes contributed less to the total number of reads than GCs and multi-tissue samples (MT: 12 different tissues). For example, the ten highest oocyte-expressed genes contributed 4.2% to 4.7% of the total reads, while the highest expressed GC and MT genes contributed respectively 5% to 9.2% and 16.84% of the total reads. In oocytes, the most expressed genes included ZP2, ZP3 and GDF9, genes known to be oocyte specific, but also ACTG, PAIP1, PDK3. In GCs, the most expressed genes included mainly ribosomal genes, such as RPL4, RPL14, RPS23,

genes involved in cellular assembly and organization, such as ACTG, collagens, VIM and cell growth regulators, such as MGP.

Among these 15 349 expressed genes, we observed a significant decrease in oocyte-expressed genes during early follicular development, i.e., from 13 219 genes in primordial oocytes to 11 693 genes in small antrum oocytes (Figure 2A). Conjointly 9 499 genes were found to be expressed in all the samples (O + GCs + MT) and were considered as ubiquitously expressed genes. However, their relative expression levels varied depending on the type of cell, for example ribosomal genes (RPL4, RPL14), VIM, SPARC, and low expressed genes such as BSCL2.

In addition to this qualitative assessment, the transcriptome data set was analyzed using the R software package, DESeq [13]. Principal component analysis of the DESeq normalized data confirmed the relevance of the data (Figure 2B). The first axis explained 52% of expression variability and separated the two follicular compartments (O and GC). The second axis explained 15% of expression variability and separated the samples according to their follicular stages.

## Global differential gene expression

To identify differentially expressed genes, a generalized linear model was applied to the data after DEseq normalization (DESeq's GLM method using the gene
<figure><img src="image_4.png" /></figure>

data set). This model allowed statistical evaluation of the factor effects (O, GC) and their interactions. The following criteria including FDR < 0.5% to 1% (Benjamini-Hochberg procedure) and a fold change >2 were chosen to determine significant differentially expressed genes. Using these criteria, 33.4% of genes (FDR < 0.5%: 5 130 genes) were significantly differentially expressed between the two compartments; 19.6% of genes were significantly differentially expressed during early development (FDR < 1%, 3 015 genes). We observed differences in expression profiles during early folliculogenesis between the two compartments for 4.3% of the genes (interaction

effect FDR < 0.5%, 674 genes). For example, the expression of the *KCNN3* gene decreased in GCs during early development whereas its expression increased in oocytes. Differential expression is detailed in Table 1, which, interestingly, shows that genes were under expressed in oocytes compared to in GCs (2 833 genes).

## Focus on differential expression between compartments (O/GCs)

We focused our analysis on the 5 130 significant differentially expressed genes (DEG) between compartments (Additional file 3). The aim was to provide an overview
<table><caption>**Table 1 Summary of differential gene expression**</caption><thead><tr><th rowspan="2"></th><th rowspan="2"></th><th rowspan="2">Total</th><th colspan="2">Number of differential genes</th></tr><tr><th>Up</th><th>Down</th></tr></thead><tbody><tr><td>Compartment effect: O/GCs</td><td></td><td>5129</td><td>2297</td><td>2832</td></tr><tr><td rowspan="3">Stage effect: SA/PD</td><td>total</td><td>3015</td><td>1357</td><td>1658</td></tr><tr><td>oocyte</td><td>2173</td><td>1064</td><td>1109</td></tr><tr><td>granulosa cells</td><td>1192</td><td>408</td><td>784</td></tr><tr><td rowspan="3">Interaction effect</td><td>total</td><td>674</td><td></td><td></td></tr><tr><td>oocyte</td><td>522</td><td>250</td><td>272</td></tr><tr><td>granulosa cells</td><td>262</td><td>91</td><td>171</td></tr></tbody></table>

of specific gene expression, mechanisms and functions in these cell types, which may be involved in the dialog between the oocyte and GCs, and regulate follicle development.

First, we observed numerous genes on the DEG list already known to have cell-type specific expression. The expression of 39 genes known to be compartment specific was investigated and the expression of 36 of them was confirmed by the RNA-seq statistical analysis (Additional file 4).

Last, the expression of 33% of the genes (1694 genes) varied during early follicular development. Supervised hierarchical clustering was performed to identify groups of genes with similar expression profiles. This led to a clear separation between the two compartments and identified clusters of genes according to the stage of follicular development, suggesting a difference in the dynamics of transcription between oocytes and GCs (Figure 2C).

## Compartment specific expression

Previous studies that were mainly devoted to rodent species identified oocyte-specific genes with important functions related to oocyte development and folliculogenesis. However, few data are available for GCs. To gain insight into the specific expression of each compartment, we performed a pair-wise comparison of the transcriptome of each compartment and both MT and the other compartment transcriptomes. In this way, we identified 55 highest differentially expressed genes in GCs (FDR < 0.5%, FC > 5) and 161 highest differentially expressed genes in oocytes (FDR < 0.5%, FC > 10. In so doing, we confirmed the specific expression of the genes already known as WNT4, FSHR, PTEN in GCs or WEE2, DDX4, NRLP9, NLRP13, BMP15 in oocytes. We also identified new compartment-specific expression of genes such as IFNE, DEFB112, BMP1, IGFBP6 and TCF23 in GCs and ST8SIA3, MUSK, TECTA, BTG4 and ACCSL in oocytes. Among these genes, we also identified oocyte-

specific expression of 11 transcription regulators (BSX, PATL2, DMRTC2, SIX6, etc.) and the GC-restricted expression of four transcription regulators (HDAC9, LZTS1, HOXC4, ASB9). For a list of the highest differentially expressed genes, see Additional file 5.

## Validation of differential expression by quantitative real time RT-PCR

To validate the RNA-seq analysis, the expression profiles of a subset of genes of interest were monitored using qRT-PCR. These genes were selected as follows: (i) genes showing significant up or down expression in oocytes (11 genes: Table 2), (ii) genes involved in particular canonical pathways (2 genes: Table 2) or (iii) compartment-specifically expressed genes (17 genes: 8 oocyte-specific genes and 9 GC-specific genes; Table 3). Statistical analysis confirmed the differential expression between compartments observed in RNA-seq data for all the genes tested (Table 2) and the compartment-specific expression of 13 out of 17 genes (Table 3). The accuracy of our gene expression measurement was also confirmed by the similarity between RNA-seq and qRT-PCR fold change (Additional file 1: Figure S4). Figure 2D illustrates the new compartment-specific expressions checked by qRT-PCR in oocytes (BTG4, MCF2L, ST8SIA3) and in GCs (IFNE, ITIH5, DEFB112).

## Biological trends

### Biological functions

To investigate the biological functions associated with the two compartments, the DEG list was subjected to Ingenuity Pathway Analysis (IPA) and tested for function

<table><caption>**Table 2 qRT-PCR validation**</caption><thead><tr><th rowspan="2">Canonical pathways</th><th rowspan="2">Gene</th><th rowspan="2">RNA-seq analysis</th><th colspan="2">qPCR analysis</th></tr><tr><th>O/GCs Fold change</th><th>O/GCs Fold change</th><th>P value class</th></tr></thead><tbody><tr><td rowspan="3">IGF1</td><td>IGF1</td><td>0.284</td><td>0.52</td><td>*</td></tr><tr><td>INSL3</td><td>0.121</td><td>0.482</td><td>*</td></tr><tr><td>IGF1R</td><td>3.334</td><td>7.35</td><td>*</td></tr><tr><td rowspan="2">Gap junction</td><td>GJA4</td><td>20.84</td><td>4.247</td><td>*</td></tr><tr><td>GJA1</td><td>0.117</td><td>0.278</td><td>***</td></tr><tr><td>Notch</td><td>JAG1</td><td>2.053</td><td>2.909</td><td>***</td></tr><tr><td rowspan="3">PI3K</td><td>PIK3R3</td><td>P</td><td>P</td><td></td></tr><tr><td>KITLG</td><td>0.266</td><td>0.218</td><td>**</td></tr><tr><td>KIT</td><td>2.244</td><td>3.962</td><td>**</td></tr><tr><td>Paxilin</td><td>PARVA</td><td>P</td><td>3.5</td><td>*</td></tr><tr><td rowspan="3">Others</td><td>GATA2</td><td>7.028</td><td>10.906</td><td>*</td></tr><tr><td>AMHR2</td><td>+∞</td><td>0.64</td><td>***</td></tr><tr><td>MAEL</td><td>3.255</td><td>5.661</td><td>***</td></tr></tbody></table>

<footnote>*: pval < 0.05; **: pval <0.01; ***: pval < 0.005.</footnote>
<table><caption>Table 3 Validation of gene expression specificity</caption><thead><tr><th rowspan="2">HUGO gene name</th><th colspan="2">RNAseq data</th><th colspan="2">qPCR data</th></tr><tr><th>O/GCs FC</th><th>Pvalue</th><th>O/G FC</th><th>Pvalue</th></tr></thead><tbody><tr><td>WEE2</td><td>40.25</td><td>6.46E-15</td><td>56.15</td><td>*</td></tr><tr><td>MUSK</td><td>22.52</td><td>1.40E-10</td><td>108.64</td><td>*</td></tr><tr><td>SPO11</td><td>13.16</td><td>7.71E-06</td><td>16.45</td><td>***</td></tr><tr><td>GAPT</td><td>+ 999<sup>1</sup></td><td>2.5E-05</td><td>60.91</td><td>*</td></tr><tr><td>MCF2L</td><td>41.16</td><td>3.82E-03</td><td>29.04</td><td>NS(0.068)</td></tr><tr><td>ST8SIA3</td><td>16.54</td><td>4.4E-03</td><td>74.5</td><td>*</td></tr><tr><td>TNIP3</td><td>41.79</td><td>1.1E-03</td><td>21.74</td><td>NS(0.072)</td></tr><tr><td>GABRP</td><td>28.89</td><td>1.79E-07</td><td>1369.57</td><td>*</td></tr><tr><td>IFNE</td><td>-999<sup>2</sup></td><td>5.18E-05</td><td>0.002</td><td>*</td></tr><tr><td>TCF23</td><td>0.096</td><td>5.29E-08</td><td>0.051</td><td>***</td></tr><tr><td>TLL2</td><td>0.118</td><td>3.00E-11</td><td>0.148</td><td>*</td></tr><tr><td>EFS</td><td>0.086</td><td>5.39E-06</td><td>0.24</td><td>***</td></tr><tr><td>WINT4</td><td>0.095</td><td>4.63E-03</td><td>0.106</td><td>*</td></tr><tr><td>ITIH5</td><td>0.142</td><td>4.97E-03</td><td>0.101</td><td>***</td></tr><tr><td>HTRIF</td><td>0.007</td><td>2.21E-06</td><td>0</td><td></td></tr><tr><td>TMEM167A</td><td>0.012</td><td>8.16E-04</td><td>0.025</td><td></td></tr><tr><td>DEFB112</td><td>0.087</td><td>1.60E-12</td><td>0.005</td><td>*</td></tr></tbody></table>

enrichment (Additional files 6 and 7). In the oocyte compartment (2 297 genes), we found a high statistical enrichment of genes involved in the cell cycle (with 14 genes contributing to the meiotic arrest), gene regulation, cell morphology and organization, and the reproductive system. The cell cycle function (mainly cell cycle progression) was also enriched in GCs. Last, we observed different involvement of oocytes and GCs in the reproductive system. On one hand, oocytes over-expressed 66 genes that have been shown to play important roles in fertility. The genes over-expressed in GCs (2 833 genes) were found to be mainly related to cellular proliferation and atresia ("cellular growth and proliferation", "cell death" IPA functions), cell morphology (148 genes involved in the "formation of protrusions") and movement, lipid metabolism and the reproductive system such as *BMP15*, *GDF9*, *BOLL*, *FIGLA*, *MOS*, *PIWIL1*, *INSL6*. On the other hand, GCs over-expressed 107 genes involved in ovarian tumors of which eight were related to GC cancer (*AMH*, *FOXL2*, *GATA4*, *GNAS*, *INHA*, *KRAS*, *NR5A1*, *ZFPM2*). Figure 3A is a graphic representation of the DEG biological functions in which the 17 main categories (FDR <0.05) are shown.

### Canonical pathways

Based on the results of the IPA canonical pathway analysis, we identified 36 significant pathways that encompassed metabolism, regulatory and cell signaling. These

pathways are shown in Figures 3B-D. This analysis underlined two different cell cycle checkpoint regulations according to the compartment: enrichment of genes involved in the G1/S checkpoint for GCs and in the G2/M DNA damage checkpoint for oocytes (Figure 3B). We observed significant cellular proliferation and apoptosis mechanisms in GCs, and gene enrichment in nuclear signaling receptors such as PPARG, LRX/RXR and PRX/RXR involved in lipid metabolism. We also identified GC-specific enrichment pathways such as cAMP-mediated, p38 MAPK, wnt/β-catenin and PDGF signaling pathways. Figures 3C-D provide further evidence for the over-expression of growth hormone-encoding genes in GCs compared to oocytes during early ovarian folliculogenesis, and for the involvement of PTEN, calcium, eNos and SHH signaling pathways in GCs. Last, we observed enrichment of genes involved in AMPK signaling related to GC energy regulation.

### Cell-cell interactions

In mammals, bidirectional communications between oocyte-GC and GC-GC are required for meiotic maturation, acquisition of oocyte competence, and follicle development. Analysis of the canonical pathways identified molecular events involved in these communications including IGF1, VEGF, FGF, Notch and tight junction signaling pathways. Figure 4 illustrates oocyte-GC cross-talk with the over expression of *IGF1R* gene in oocytes (FC = 3.33) and the over expression of *IGF1* gene in GCs (FC = 3.57).

Regulation of the IGF1 system was then achieved through the over expression of *IGFBP2*, *IGFBP4*, *IGFBP6*, *IGFBP7* genes in GCs. In the same way, *VEGFA* was highly expressed in GCs (FC = 4.87) while its receptor *FLT1* was over expressed in oocytes (FC = 2.86) (Figure 5A). The involvement of *FGF2* in oocyte-GC analog was also evidenced by the over expression of *FGF2* in oocytes and of *FGFR1* and *FGFR2* in GCs. We also revealed for the first time over expression of FGF family members such as *FGF16* (FC = 7.14), *FGF20* (FC = 4.3) in oocytes and *FGF1* (FC = 8.3), *FGF11* (FC = 4.3), *FGF12* (FC = 6) in GCs. Notch signaling was assumed to be involved in two kinds of interactions. Firstly, the abundance of *NOTCH1* transcripts in GCs and of its ligand *JAG1* in oocyte confirmed compartmental interactions (Figure 5B). Secondly, this pathway could be involved in GC interactions via its functional ligand *CNTN1* (Figure 5C).

Finally, we detected the differential expression of genes involved in cellular connections and adhesions. *GJA1* (*cx43*), which is over expressed in GCs, and *GJA4* (*cx37*), which is over expressed in oocytes (Table 2), are known to participate to gap junctions and connect oocytes with GCs. We also observed over expression of

<footnote>*: pval < 0.05; **: pval < 0.005.</footnote>

<footnote n="1">¹: expression only detected in oocytes.</footnote>

<footnote n="2">²: expression only detected in GCs.</footnote>
<figure><img src="image_5.png" /><figcaption>**Figure 3** Functional enrichment and canonical pathways within oocyte and GCs. **A:** Functional enrichment. Genes differentially expressed between the 2 compartments were evaluated in *silico* using Ingenuity Pathway Analysis (IPA). The figure shows the most significant functional groups (p < 0.05, FDR < 5%). For each compartment, the gene list corresponds to up regulated genes compared to the other compartments. The bars represent the p-value at *logarithmic* scale. GC data are in red and oocyte data are in blue. Red and blue numbers correspond to the number of focus genes that contributed to the function concerned. **B:** Canonical pathways in oocytes and GCs. Significantly enriched canonical pathway categories (p-value < 0.05) identified in the 2 compartments. The Y axis corresponds to percentage gene enrichment in the pathway. Red and blue numbers correspond to the number of focus genes that contributed to the pathway concerned. *: FDR < 0.25 (likely to be valid 3 out of 4 times); **: FDR < 0.10; ***: FDR < 0.05. **C:** Growth factor signaling pathways in oocytes and GCs. Significantly enriched growth factor signaling pathways (p-value < 0.05) identified in the 2 compartments. The Y axis corresponds to percentage gene enrichment in the pathway. Red and blue numbers correspond to the number of focus genes that contributed to the pathway concerned. *: FDR < 0.25; **: FDR < 0.10; ***: FDR < 0.05. **D:** Signaling pathway in oocytes and GCs. Significantly enriched signaling pathways (p-value < 0.05) identified in the 2 compartments. The Y axis corresponds to percentage gene enrichment in the pathway. Red and blue numbers correspond to the number of focus genes that contributed to the pathway concerned. *: FDR < 0.25; **: FDR < 0.10; ***: FDR < 0.05.</figcaption></figure>

genes involved in the adherent junctions (catenin, cadherins and afadin) and in the tight junctions (TJP1, claudins, cingulin) (Additional file 1: Figure S5).

### Transcription regulators

Oocyte and GC transcription factors (*FIGLA*, *NOBOX*, *FOXL2*, etc.) play an important role in folliculogenesis [14]. They control both oocyte development and GC function. To identify new transcription regulators, we performed IPA upstream regulator analysis on the DEG data set. For each compartment, we predicted the

transcriptional regulators that were over expressed and involved in the regulation of this DEG: 47 GC regulators and 9 oocyte regulators (pval < 0.05: Additional file 8). Among them, the oocyte transcription regulator NOBOX and 12 GC transcription regulators (GATA4, FOS, TP53, etc.) were predicted to have their protein activated.

### Discussion

To our knowledge, this is the first comprehensive exploration of transcriptome from oocyte and granulosa cells separately and at several key stages corresponding
<figure><img src="image_6.png" /><figcaption>**Figure 4 IGF signaling pathway.** Gene expression involved in the IGF1 signaling pathway. Genes in red are over expressed in oocytes. Genes in green are over expressed in GCs. Relative expression of IGF1 and IGF1R from RNAseq data in the graph. The Y axis corresponds to RNA-seq normalized read counts.</figcaption></figure>

to *in vivo* major transitions during early folliculogenesis in sheep. The study resulted in reliable transcriptional profiling in time and space. We explored changes in gene expression between the two follicular compartments and identified specific transcriptional programs in either oocyte or granulosa cells. In addition, we identified pathways involved in bidirectional communication.

Thanks to mouse genetic models, a significant number of targeted molecular data (gene expression and proteins) has been accumulated during the past years. Recently, these data were complemented by mouse and human oocyte transcriptome data that give us more information about the molecular mechanisms involved in early folliculogenesis. Our study expanded transcriptome data to GC and enabled us to determine the precise expression of genes involved in these mechanisms in a mono ovulating species.

## Oocyte and GC Transcriptional profiling

This study documented the expression of 15 349 genes (which represents the majority of ovine genes) in ovarian

follicles during early follicular development. Thanks to the sensitivity of RNAseq, we estimated a larger number of genes expressed in oocytes (14 172 genes expressed in ¾ of replicates) than other microarray studies (Additional file 2; discussion section). In our study, 14 087 genes were also estimated to be expressed in GCs. These high numbers of expressed genes were close to the higher number of genes expressed in testis (15 000) compared to the other human and mouse tissues (11 000 to 13 000 genes from RNA-seq data) [15].

The high number of expressed genes in both oocytes and GC evidenced strong transcriptional activity. Interestingly, this analysis revealed a slight decrease in the number of oocyte expressed genes during follicular development from SC follicles that coincides with the gain in methylation from PM to SC follicle transition (oocyte diameter of at least ~50 mm) [16]. The detection of a large number of ubiquitous genes (Figure 2A) was previously reported by Ramköld et al. [15] suggesting that functional specialization is not associated with a large
<figure><img src="image_7.png" /><figcaption>**Figure 5** (See legend on next page.)</figcaption></figure>
(See figure on previous page.)

**Figure 5 VEGF and NOTCH signaling pathways.** Gene expression involved in the VEGF and NOTCH signaling pathways. Genes in red are over expressed in oocytes. Genes in green are over expressed in GCs. Relative expression of VEGF and FLT1 from RNAseq data. The Y axis corresponds to RNA-seq normalized read counts. **A:** Oocyte-GC communications in the NOTCH signaling pathway. **B:** GC-GC communications in the NOTCH signaling pathway. **C:** Oocyte-GC communications in the VEGF signaling pathway.

population of genes but results in more subtle and com-
plex control.

RNA-seq has made it possible to identify novel tran-
scripts. In this study, we uncovered 221 716 genomic
fragments that now need to be annotated. Some of them
present a high expression level and compartment-
specific expression (FC > 100) suggesting they play an
important role in early folliculogenesis. Further explor-
ation using bioinformatics and gene expression valid-
ation will be necessary to refine these results.

Differential expression between compartments reflects
the fate of the cell types

The development of the follicle and the maturation of
oocyte require orchestrated communication between oо-
cytes and GCs at all stages of early folliculogenesis. They
also require the expression of compartment-specific
genes. To investigate these molecular processes, we fo-
cused our analysis on differences in gene expression
between oocytes and GCs and identified a large number
(5 129) of differentially expressed genes. We identified
functional differences between oocytes and GCs and de-
scribed the characteristics of each compartment that
lead to two specific cell fates. As described in our previ-
ous study [12], the top functions of differentially
expressed genes illustrated the three main cellular events
known to be involved in early follicular development:
the change in the shape of the GC cell from flattened to
cuboidal (IPA categories “cell movement”; “shape change
of tumor cell lines”: 30 genes including *FLNB*, *EDNRA*),
the marked increase in the number of GCs (x 40 in
sheep species up to the secondary stage) and oocyte en-
largement [17] (x 300 between the primordial and antral
stages: IPA category “morphology of germ cells”: 31
genes including *ASZ1*, *TDRD5*). DEG comparison of
genes identified in five previous studies validated this list
(Additional file 2; discussion section). In addition to the
information provided by our previous study [12], IPA
RNA-seq analysis increased available information thanks
to the identification of additional DEG.

Meiotic regulation

In most mammals, meiosis is initiated during fetal life
and stops prenatally in the diplotene stage of prophase I
at the G2/M transition until development into a Graa-
fian follicle (before ovulation).

This study highlights oocyte over expression of genes
involved in the regulation of meiotic progression,

acquisition of meiotic maturation, and maintenance of
meiotic arrest.

For instance, SPO11 is involved in initiating meiotic
recombination by catalyzing DNA double strand breaks
[18]. BTRC loss of function in mice results in slow pro-
gression through meiosis and abnormal divisions in
spermatocytes [19]. We also observed a significant
increase in sheep OVOL1 mRNA (FDR = 0.043) and
ASZ1 mRNA (FDR = 0.017) during early folliculogenesis.
Ovol1 is known to be involved in regulating the pachy-
tene progression of male germ cells [20] and shown to
be expressed during sheep early folliculogenesis [21]. In
mice, the ASZ1 protein is a germ cell-specific protein
located in the cytoplasm of oocytes at all stages of de-
velopment that may act as an important signaling pro-
tein and/or transcriptional regulator during germ cell
maturation [22].

The function of these genes is poorly understood in
mammalian ovary development but their expression in
early follicular development suggests they play a role in
the regulation of female meiosis.

Finally, we confirmed the expression of genes already
known to be involved in the regulation of the oocyte
and/or spermatocyte meiotic arrest (AURKA-C, PLK1,
PDE3A, BTRC, CDC25B-C, WEE2 etc.).

The establishment of GC steroid synthesis

This study clearly shows that GCs over-expressed genes
involved in steroid synthesis (such as *CYP11A1*, *HSD3B*
and *HSD17*) and in steroid regulation (such as nuclear
receptor genes (NRs)). NRs include orphan receptors
(NR5A1, NROB2), steroid hormone receptors (*ESR1*, *AR*,
*LXR* (NR1H2), *FXR* (NR1H4)), and retinoid acid recep-
tors (*RXRA*), all of which are over expressed in GCs dur-
ing preantral follicle stages.

NR5A1 (SF-1), which has already been described in
sheep GCs [23], is a master gene that regulates the ex-
pression of numerous genes including reproductive and
steroidogenic enzyme-encoding genes [24]. GC specific
SF-1 KO mice are infertile whose ovaries contain a re-
duced number of follicles and reduced expression of
AMH, which partially correlates with the reserve of ovar-
ian follicles [25].

E2 and androgens hormones exert their influence on
early folliculogenesis via their respective receptors, ESR
and AR. Like in mice, E2 probably plays a role in con-
trolling the primordial follicle pool [26] and androgen
signaling is crucial for normal folliculogenesis [27]. We
observed strong expression of *AR* in GCs, like in humans [28]. Finally, for the first time, we identified over expression of *LXRB*, *FXRA* and *RXRA* in GC preantral stages. RXR is a retinoid X receptor that binds as heterodimers (LXR/RXR, FXR/RXR...) and becomes transcriptionally active only in the presence of a ligand (RXR-selective ligand or ligands of the other dimer). LXRs are activated by oxysterols (generated by intermediates of steroid hormone and cholesterol synthesis (CYP21A1, STAR, P450scc products)) and are key regulators of ovarian steroidogenesis at antral stages [29]. FXRs are activated by bile acid, sterol and different lipids and are possibly important regulators of androgen hormoneostasis in male gonads [30].

## Complex regulation of follicular atresia

Atresia is an important process in folliculogenesis and concerns the majority of the follicles. Apoptosis is found in the oocytes of primordial follicles and progressively extends to GCs of growing follicles.

We characterized the expression of different genes of the BCL2 family either over expressed in oocytes (*BCL2L1*, *BCL2L10*, *BCL2L11* (BIM), *BCL2L14* (BCLG)) or in GCs (*BCL2* *BCL2L2*, *BOK*) and that mainly prevent oocyte loss. For example, the *BCL2L1* gene was shown to play a crucial role in the survival of germ cells [31]. BCL2L10 may antagonize BAX activities and may also play other roles related to cell cycle control and oocyte maturation (for a review see [32]). In sheep, the high expression of *BCL2L10* in oocytes and its up-regulation during preantral stages reinforces its already important role in oocyte development. By decreasing GC atresia, the *BCL2* gene family expressed in GCs, such as *BCL2*, was also believed to increase the number of germ cells and enhance their development. Lastly, *BOK* mRNA, was found to be restricted to reproductive tissue and was detected in GCs from PM to antral follicles in rat species [33].

## GC signaling pathways

The molecular pathways were mainly linked to GC proliferation, steroid production, and morphological change.

Sheep GCs over expressed members of the Sonic hedgehog (Shh) signaling pathway (*SMO*, *PTCH1* and *GLI1*-2-3), which is consistent with a potential role for the HH pathway in increased GC proliferation during preantral stages [34]. Last, genes of the WNT family and their receptors (*WNT4*-3A-5A-6, *FDZ1*, *LRP1*-11) and Rho-GTPase activating proteins (*ARHGAP9*-20-21-31-36-44) were shown to be over expressed in GCs. In the ovary, the β-catenin dependant canonical pathway of WNT4 leads to ovarian determining pathways. WNT proteins can also signal via Rho-GTPases, using non-canonical pathways that involve changes in polarized cell

shape and cell migrations [35]. This pathway may promote changes in cell shape that affect GCs during early folliculogenesis.

## Cell-cell communications

Follicular growth depends on close interactions between oocytes and GCs but also between GCs themselves. Two kinds of cell interactions are described in this study demonstrating the existence of complex dialogs between compartments: a molecular dialog (signaling pathways) and physical communications (gap junctions and transzonal projections).

### Signaling pathways: GC regulation of oocyte development

GCs are known to play a role in supporting the growth of oocytes, regulating the progression of meiosis, and modulating global oocyte transcription activity. Nevertheless, so far, few factors have been identified that support these assumptions. To date, the most widely studied ligand-receptor systems (that were also identified in the present study) are the KIT/KITL and BMP pathways. Below we describe the involvement of other signaling pathways.

The Notch signaling pathway plays many roles in cell communication by influencing cell proliferation, differentiation, and apoptosis [36]. In sheep, we found over expression of *NOTCH1*-2-3 in GCs and over expression of the NOTCH ligands *JAG1* and *DLL3* in oocytes that perfectly illustrate oocyte-GC crosstalk (Figure 5A) and that may influence oocyte apoptosis, like in postnatal mice [37]. In addition, we suggest that NOTCH is involved in GC proliferation [38] through NOTCH1/CNTN1 binding (Figure 5B).

Vascular endothelial growth factor (VEGF) was already known to promote early folliculogenesis [39]. In this study, we identified over expression of members of VEGF pathway including *VEGFA* and *NRP1* (a VEGF receptor) in GCs, and over expression of *FLT1* (another VEGF receptor) in oocytes. These gene expressions suggest a role for the VEGF pathway in GC-oocyte and GC-GC crosstalks (Figure 5C). As already described, the VEGF pathway may protect GC against atresia [40].

The crucial role of the IGF system is largely described during terminal follicular growth [41] but is poorly documented in preantral stages. Ovarian IGF1 production was shown to vary with the species. IGF1 is mainly expressed by GCs in rat and pig species [42] but its protein was also detected in the oocyte of rat preantral follicles [43]. Moreover, the *IGF1* transcript was not detected in adult sheep [44] or cattle ovaries [42] but its receptor *IGF1R* was expressed from primary follicles and its expression increased with follicular growth [44,45]. Our study provides a precise description of members of the IGF system in sheep species at preantral stages (*IGF1*-2, *IGFBPs*, *IGF1R*).
The IGF1 pathway (Figure 4) suggests a role for IGF1 in oocytes and GCs in agreement with studies on mice. Indeed, although insulin/IGF signaling is not essential in oocyte conditional knockout mice [46], IGF-1 significantly increases oocyte maturation [47]. In GCs, the expression of *IGFR1* and the over expression of the downstream signaling molecules *IRS1* and *AKT2* are also in agreement with the involvement of IGF1 in the regulation of GC survival [48]. Additional immunohistochemistry to detect IGF1 and IGF1R proteins and *in vitro* early follicle cultures could be used to evaluate the effect of IGF1 on sheep preantral growth. *In vitro* early follicle cultures are described in Mochida et al [49].

### Signaling pathways: Oocyte regulation of GC function

Numerous studies have shown that oocytes are potent stimulators of GC proliferation. A compilation of oocyte-secreted factors that regulate GC function was published in 2008 [50]. Among these factors, we describe the over expression of *FGF2* in sheep oocytes and *FGFR1-2-L1* receptors in GCs that may influence primordial follicle development in accordance with *in vitro* goat and rat studies [51,52]. We also report the over expression of *FGF16-20* in oocytes for the first time in mammals. A recent study on Nile tilapia species suggests that FGF16-20 are involved in early oocyte development in the female [53]. In other mammal tissues, FGF16 and FGF20 promote cell proliferation and differentiation.

### Physical interactions

Transzonal projections (TZP) are GC specialized extensions corresponding to long cell projections and microvilli that are in close contact with the surface of the oolemma or invaginate the cortical area of the oocyte. These projections establish communicative (gap junctions) and adhesive junctions (tight junctions) both with the oocyte and neighboring transzonal projections. In this study, we observed an enrichment of genes involved in the formation of cellular protrusions in both GCs and oocytes (Additional files 5 and 6). We confirmed over expression of two gap junction genes that are indispensable for the regulation of folliculogenesis and oogenesis: *GJA1* (*Cx43*) in GCs and *GJA4* (*Cx37*) in oocytes. Moreover in our study, *ZP2-4* were seen to significantly increase during preantral stages as shown in human fetal ovaries before the formation of *zona pellicida* [54].

### Regulators

In the present experiment, we identified a set of upstream transcriptional regulators able to regulate follicle gene expression and ovarian development.

In sheep oocytes, *NOBOX*, which is known to promote follicular growth in mice, was detected in this study at

medium expression level in preantral oocytes and was predicted to increase gene expression of *BMP15*, *DNMT1*, *GDF9*, *HIFOO*, *MOS*, *ZAR1* genes [55]. Analysis of IPA upstream regulators also predicted activation of AHR in oocytes (Additional file 9). AHR protein may act on oocyte gene expression through the receptor nuclear translocator *ARNTL2* over expressed in oocytes. These findings are consistent with AHR regulation of follicular growth by affecting germ cell death during female gametogenesis [56].

In GCs, transcriptional regulators referred to sterol regulation (*SREBF1*, *SREBF2*) hormone response (*EGR1*, *CREB*, *SKIL*, *ETS1*), cellular growth (*TP53*, *XBP1*, *MYC*, *NFIA*, *HIF1A*), cellular differentiation (*MITF*), and increased transcription (*MYC*, *NFIA*). Here we describe for the first time, over expression of *NFIA* and *HIF1A* in GCs of preantral follicles. NFIA is a transcription factor involved in the control of cell growth in both humans and model systems [57]. HIF1A is an ARNT interacting protein that activates the transcription of target genes involved in energy metabolism, angiogenesis, and apoptosis [58]. We were able to predict the regulation of 59 genes by *HIF1A* and we observed an increased in its expression during early folliculogenesis, assuming it plays a role in preantral development.

### Genes preferentially expressed according to cell type

Oocyte-specific genes and their role in folliculogenesis have been extensively studied in murine knockout models [59] but there is a relative paucity of information concerning GC-specific genes. In this study, we completed existing data by highlighting the expression of 216 preferentially expressed genes (161 O + 55 GCs), including genes whose function in ovarian folliculogenesis is unknown. Nevertheless, GCs exhibited fewer preferentially expressed genes than oocytes.

In GCs, the *DEFB112* gene is described for the first time as preferentially expressed. This gene is poorly documented but is reported to have a regionalized expression in the epididymis in rat and in human and to be related to sperm maturation [60]. As already described in sheep, the preferential expression of *BMP1*, which promotes BMP signaling via chordinase activities [61], was confirmed. Finally, the molecular chaperone *TCP1*, which assists the folding of proteins upon ATP hydrolysis (such as actin and tubulin), was specifically detected at a slight but constant level of expression during the preantral stage. TCP1 appears to play a key role in the cytodifferentiation of spermatids that relies on the high plasticity of microtubules [62]. Our result may be associated to the remodeling of the microtubule cytoskeleton during the change in cell shape.

In oocytes, we identified the preferential expression of genes such as *BTG4*. Although much about the function of BTG4 remains unknown, it is known to be involved
in cell cycle arrest and to be endowed with antiproliferative properties. In mice, its expression is restricted to the olfactory epithelium, testis, and oocyte [63]. The role of *BTG4* during preantral stages in sheep requires further investigation. Another gene of interest is *PATL2*, a RNA-binding protein (called P100) that was originally characterized as being oocyte specific in *Xenopus laevis*, where it acts as a translational repressor in a CPEB-containing complex [64]. In *Xenopus laevis*, *PATL2* is maternally expressed in immature oocytes, but disappears upon meiotic maturation. Some authors suggest that PATL2 plays a role in regulating the translation of specific maternal mRNAs required for the progression of *Xenopus laevis* oocyte maturation [65]. The expression of *PATL2* has not yet been documented in mammals; but here we describe for the first time increased specific expression of this gene in oocytes during preantral stages in sheep. In comparison with the role of *PATL2* in *Xenopus laevis*, our finding suggests it plays an important role in the regulation of the translation in preantral oocytes in sheep. Finally, we identified over expression of *MUSK* in oocytes at preantral stages. MUSK is a receptor tyrosine kinase (RTK) which is induced during skeletal muscle differentiation, and is dramatically down-regulated in mature muscle [66].

This study highlighted involvement of new gene ex-
pression in early folliculogenesis. Further investigations
are underway to estimate the implication of these genes
in mammalian folliculogenesis.

Conclusions

This work generated a significant amount of expression data that accurately document early folliculogenesis in sheep. The abundance of expressed genes and the differential expression profiles found in oocytes confirm the hypothesis that oocytes are dynamic cells. The decrease in the number of expressed genes following the development of oocytes is in accordance with a general decrease in transcription.

This study identified specific functions and molecular
mechanisms involved in each cell type during early folli-
culogenesis and in the dialog between the different cell
types. We confirm and incorporate a large number of in-
dividual studies in different species. Last, we highlight
new gene expressions that support and extend current
knowledge on early folliculogenesis.

This study is the first thorough, integrated overview of
the molecular mechanisms involved in early folliculogenen-
esis in a mono-ovulating species. It provides a valuable
data base for future functional studies in mono-ovulating
species including humans.

Methods

Laser capture microdissection

Ovaries were sampled during the preliminary experiment
[12] which was conducted in compliance with institutional

guidelines for research studies (animal experiment-
ation authorization no 31-297, prefecture de la Haute
Garonne). As described in that study, the ovaries were re-
moved from lambs at birth and embedded in O.C.T. em-
bedding matrix, frozen in liquid nitrogen and stored at
-80°C until use. LCM capture was performed using the
Arcturus XT apparatus (Arcturus, Applied Biosystems®)
from eight-micrometer frozen sections that were stained
individually using Cresyl Violet.

Each follicle stage was selected as a function of follicu-
lar diameter under the microscope taking the observa-
tion of series of sections into account and according to
the classification of Lundy *et al*. [17]. Briefly, the follow-
ing criteria were applied for the selection of the follicle
stage : i) primordial follicles (PD <35 µM) were defined
as an oocyte surrounded by 2-3 flattened GCs and
no cuboidal GCs, ii) primary follicles (PM: 35-50 µM)
were defined as a monolayer of cuboidal GCs whatever
the size of the oocyte, iii) secondary follicles (SC: 60-
120 µM) were defined as two layers of GCs and iv) small
antral follicles (SA: 250-500 µM) had a diameter of less
than 500 µm. The follicular compartments (i.e. granulosa
cells (GCs)) and oocytes (O) for each follicular stage
were collected in separate CapSure™ HS LCM caps (Arc-
turus, Ref LCM 0214), treated with 10-15 µl of extrac-
tion buffer from Picopure RNA Isolation kit (Arcturus,
ref: KIT0202) and stored at -80°C until use.

RNA extraction

As very few cells were captured per cap, each RNA LCM-
derived sample comprised pooled caps and was extracted
using the PicoPure RNA Isolation Kit according to the
manufacturer's instructions, including on-column DNase
treatment (Qiagen, ref 79254, Courtaboeuf, France) [12].
Finally, RNA of three/four independent replicates was ex-
tracted per biological condition (O and GC compartment
from PD, PM, SC, SA stages).

The 13 tissue samples (fetal ovary, pituitary gland, hypothalamus, muscle, skin, heart, lung, intestine, stomach, kidney, spleen, liver, and theca) were sampled in triplicate at the local slaughterhouse (3 individual animals), frozen in liquid nitrogen and stored at -80°C until RNA extraction. Total RNA was isolated using the Nucleospin RNA II kit (Macherey-Nagel GmbH & Co, ref 740 955 50, HOERDT, France). The 12 total RNA tissue samples (except fetal ovary) were pooled in three independent groups.

T7 linear amplification

A total of 31 RNA LCM-derived (3PDO, 4PMO, 4SCO,
4SAO, 4PDG, 4PMG, 4SCG and 4SAG) and three multi-
tissue RNA samples were subjected to two rounds of T7
linear amplification using the RiboAmp®HS PLUS kit
(Arcturus, ref KIT0525), as previously described [12].
After in vitro transcription, the optical density of anti-sense RNA (aRNA) was measured at 260 and 280 nm. The integrity of the aRNA samples was checked by the Agilent Bioanalyzer 2100 with a RNA 6000 Nano Lab Chip. To complete the check, we added four *Bacillus subtilis* control transcripts with decreasing quantities and transcript length (6.54 kb; 3.36 kb; 2.4 kb and 1.6 kb) (provided by Affymetrix, GeneChip Eukaryotic Poly-A RNA Control Kit, ref 900433) to each LCM-derived RNA and multi-tissue RNA sample before amplification (1 µl of 5*10⁻⁷ dilution). These control transcripts represented internal controls of gene expression for the experimental process.

## Library preparation and massive parallel sequencing

All RNA-seq libraries were constructed from 200 ng of aRNA using Illumina TruSeq RNA Sample Prep (Illumina, San Diego, CA USA) according to the manufacturer's instructions (Additional file 1: Figure S1). Each biological replicate was identified by a different index. Libraries were quantified using the QPCR NGS Library Quant Kit (Agilent). Three libraries were pooled and loaded randomly per lane at a concentration of 8 pM. The PDG4 library was loaded twice as a technical replicate (PDG4, PDG4b). The 11 flow cell lanes were split into five runs (four runs for LCM derived RNA samples and one run for RNA multi-tissue samples). Sequencing was performed on an Illumina HiSeq2000 using the Illumina TruSeq SBS kit v2 (209 cycles including the index) to obtain paired-end reads (2x100 bp). The quality of sequencing results was summarized and plotted using the NG6 storage platform (http://ng6.toulouse.inra.fr) [67].

## Mapping strategy

The sheep genome sequence recently became available in public databases but the 3' UTR non-coding sequences of ovine genes were still poorly documented. Due to the specificity of our data (enrichment of 3' UTR non-coding sequences [68]), we evaluated two different assembly strategies to maximize transcript identification (Additional file 1: Figure S6A).

## Genome strategy

For each sample, the reads were mapped with bwa aln to the sheep genome sequence [69] (CSIRO Oarv2.0 released March 2011: http://www.livestockgenomics.csiro.au/sheep/). Bam files were filtered (keeping only reads in a proper pair with mapq > 30) and converted into BED files. All BED files were individually merged (overlapped genomic fragments or which shared reads from same pair) and then grouped before being merged to produce a collection of genomic fragments. The number of reads per genomic fragment was computed and only genomic fragments with more than 20 reads were kept to perform annotations.

## De novo transcriptome strategy

Reads were assembled per condition using SOAPdenovo-Trans-31kmer. After each assembly, GapCloser was used to resolve N-spacers introduced in the scaffolding process. All contig files were thereafter concatenated and a meta-assembly was performed using TGICL. All contigs and singlets highlighted were used to create a reference to map back all input reads. Bam files were filtered to keep only reads in a proper pair with no more than four mismatches. A counting file was generated for each contig and each condition using all count files.

## Annotation strategy

In comparison to sheep, more gene annotations are available for cattle. In agreement with Jager [70], our preliminary analysis showed that the use of bovine transcripts for the annotation substantially increased the number of genes identified. Consequently, genomic fragments and *de novo* contigs were aligned with the bovine genome sequence [71,72] (Ensembl *Bos taurus* UMD3.1) using Blat. Blat results were filtered to keep only hits with genomic fragment or contig coverage > = 80%, identity > = 90% and best hit score 3 fold greater than other hits. Blat targets (i.e. cattle slices) were crossed with the Ensembl *Bos taurus* GTF file (*Bos taurus*.UMD3.1.65.gtf). Iteratively, for *Bos taurus* slices not linked to a described gene, we searched for overlaps with downstream gene regions extended by 500 bp, 1 kb and 3 kb (Additional file 1: Figure S6B).

All the non-annotated genomic fragments and contigs (not mapped to the bovine transcript) were discarded from the data set. Finally, the annotated contigs (from the *de novo* strategy) were mapped to the sheep genome and contigs whose sequence was unknown in the public sheep genome database were added to the final data set.

## Dataset filtration

As a gene was represented by a mean of eight genomic fragments (Additional file 1: Figure S8) and because most of the RNA-seq reads are located at the end of the genes (amplification bias), the genomic dataset was filtered to conserve only the best genomic fragment per gene (the genomic fragment with the highest count and the most 3' UTR localization) with more than 20–30 reads per condition (20 reads for PDO in the three replicates and 30 reads for the other conditions in the four replicates). This counting method was the most informative for this experiment and did not require gene length normalization. This genomic data set was completed with the annotated contigs from the *de novo* strategy missing in the genome strategy data set.

## Experimental validation

Increasing amounts of four *Bacillus subtilis* control transcripts were added to all LCM-derived RNA samples
before amplification to control the amplification and
RNA-seq processes. For each sample, the reads were
mapped with bwa aln on the 4 *B. subtilis* sequences. The
RPKM (reads per kilobase of exon per million) was com-
puted for each control transcript. To make the samples
comparable for experimental validation, the amount of
each *B. subtilis* control transcript was calculated for each
sample in proportion to the highest one (DAP).

Statistical analysis

The significance of differential gene expression was deter-
mined using DESeq package for multifactor design [13] in
R software (R 2.14.0; DESeq release 1.6.1) (Additional file
1: Figure S6C). Our experimental design included 2 fac-
tors: stage (with 4 levels: PD, PM, SC and SA) and com-
partment (2 levels: O, GCs) with interactions. After
DESeq normalization, the best dispersion estimation of
our data set was obtained using the following arguments:
method ‘pooled’, sharingMode ‘fit-only’, fitType ‘local’. The
GLM method was then used to estimate the coefficient
and deviance for each gene. For inference, we specified
three models: the full model regressing gene expression
on both compartments and follicular stages and taking
interaction effects into account (count ~ stage\*compart-
ment), the additive model (count ~ stage + compartment)
and the reduced model (count ~ 1). We tested the effects
of the compartments, follicle stages, and their interactions
on gene expression using nbinomGLMTest. Differential
gene expression between compartments and differential
expression profiles between compartments (interaction ef-
fect) were selected with FDR < 0.5% (Benjamini-Hochberg
procedure) and a fold change >2. To establish the differen-
tial expression profile during early follicular development,
we selected genes whose intra-compartment variance was
differentially expressed in at least one of the stages (nbi-
nomGLMTest). Thereafter, for each compartment, the
expression value of each gene at each follicle stage
was compared using pairwise comparison (nbinomTest).
Differential gene expression during early follicular devel-
opment was selected with a FDR <5% for global ‘intra-
compartment’ effect (nbinomGLMTest) and pval <1% for
pairwise comparison (nbinomTest) with a fold change >2.

The most expressed genes in a specific compartment (O, GCs) were identified using pairwise comparisons with FDR <0.5% and a fold change >10 for the oocyte compartment (O/GCs + MT) and >5 for GC compartment (GCs/O + MT).

Biological trends

Ingenuity® Pathway Analysis software (IPA; http://www.ingenuity.com) was used to examine functional and molecular pathway enrichment for differentially expressed genes based on Fisher's exact test. This software combines functional annotations of differentially expressed

genes (focus genes) and the corresponding bibliographic
data to generate significant signaling pathways and regu-
lation networks. The biological analysis focused on the
list of differentially expressed genes. The Ingenuity Path-
ways Knowledge Base was used as the reference set.

Functional category analysis identified the functions in the
IPA library that were most significant to the input data set.
The significance of the association between the data set and
the functional category was determined by a P value calcu-
lated using Fischer's exact test corrected for multiple testing
(p-value < 0.05; FDR < 0.05 (Benjamini-Hochberg test)).

Canonical pathway analysis identified the pathways in the IPA library of canonical pathways that were most significant to the input data set. The significance of the association between the data set and the canonical pathway was determined based on two parameters: (1) the ratio of the number of genes from the data set that mapped to the pathway divided by the total number of genes that mapped to the canonical pathway and (2) a P value calculated using Fischer's exact test corrected for multiple testing to determine the probability that the association between the genes in the data set and the canonical pathway is due to chance alone (p-value < 0.05; FDR < 0.05-0.25 (FDR threshold <0.25 is proposed as appropriate for enrichment gene sets in BROAD GSEA analyses (http://www.broadinstitute.org/gsea/doc/GSEAUserGuideFrame.html). In these conditions, the result is likely to be valid 3 out of 4 times).

Downstream effects analysis was used to identify the ex-
pected effect of the observed change in gene expression on
the functions (expected to increase or decrease). This ana-
lysis examined genes in our data set that are known to
affect functions, compared the gene direction of change to
expectations derived from the literature, and provided a
prediction for each function based on the direction of
change (in the list of differential gene expression). IPA uses
the regulation z-score algorithm to make predictions.

In the same way, IPA upstream regulator analysis was
used to understand the cause of the observed change in
gene expression. This analysis identifies the cascade of
upstream transcriptional regulators (any molecules that
can affect the expression of other molecules) that can
explain the observed change in gene expression in the
dataset. The overlap pvalue measures statistically signifi-
cant differences between genes in the data set and the
genes that are regulated by a transcriptional regulator.
The activator z-score infers the activation state to the
transcriptional regulator (activator z-score < -2 and >2).
Without identified bias, this z-score can always be used
as an independent test to call downstream regulators. A
bias means that the activation z-score and pvalue must
be used to call downstream regulators. The expression
and/or differential expression of the transcriptional regu-
lator in the data set provides more evidence for the bio-
logical mechanism.
Combining oocyte and GC upstream effect/ regulator
analyses provided clues to the putative regulation in-
volved in the oocyte-GC dialog.

Quantitative real time RT-PCR analysis for gene
expression

Gene primers were designed from RNA-seq sequences
(Additional file 10). The intron–exon organization of
ovine genes was deduced by comparison with the Hu-
man genome using the NCBI database. Gene primers
were designed preferentially in the 3′ UTR (last 1000
base pairs) or the last exon using LightCycler Probe De-
sign2 software (Roche Diagnostics). The primer pairs
were confirmed by Primer3 (http://frodo.wi.mit.edu/pri-
mer3/). Sequences are available in Additional file 8.

The LCM-derived aRNA samples were reverse tran-
scribed as previously described [12]. The resulting LCM-
derived cDNA samples were completed to 50 µl and
diluted 1/20 before PCR. The assay for each gene consisted
of four replicates per condition (except for PDO = 3) and
negative controls.

Total tissue RNA samples (2 µg) were reverse tran-
scribed using Primer p(dT)15 for cDNA synthesis
(Roche Diagnostic, ref 10814270001, Meylan, France)
and 200 U of SuperScript II (Invitrogen, ref: 18064–014,
Cergy Pontoise, France) according to the manufacturer's
instructions. The resulting tissue cDNA samples were
completed to 50 µl and diluted 1/50 before PCR. The
assay for each gene consisted of a pool of three inde-
pendent replicates per tissue.

QPCR for compartment specific genes (17 genes) was
performed using 3 µl of the final dilution (~1.5 ng LCM-
derived cDNA samples (31 LCM-derived samples) and
2.4 ng tissue cDNA samples (12 tissue + fetal ovary)
using SYBR green fluorescence detection during amplifi-
cation on an ABI Prism 7900 Sequence Detection Sys-
tem 2.1 (Applied Biosystems®) as previously described
[12]. The efficiency of real-time PCR amplification was
calculated for each primer pair using five serial dilution
points from an sheep fetal ovary cDNA sample (1:6
(~16 ng cDNA/point); 1:3; 1:3; 1:2; 1:2).

The expression of genes involved in canonical path-
ways (28 genes) was analyzed using 96.96 Dynamic
Array™ IFCs and the BioMark™ HD System from Flui-
digm. Two specific target amplifications (STA) were per-
formed on the 36 cDNA samples (31 LCM-derived
aRNA samples, 3 multi-tissue aRNA samples, a calibra-
tor sample (pool of all the LCM-derived aRNA samples)
and a pool of cDNA from total fetal ovary RNA (for the
calculation of PCR amplification efficiency) in 96-well
PCR plates. A mix was prepared containing 132 µL 2X
TaqMan® PreAmp Master Mix plus (ref 4391128, Ap-
plied Biosystem) 26.4 µL 10X Preamplification Primer
Mix (500 nM each primer (for 60 genes + 2 reference

genes (β-actin and *TMED4*)), and 3.75 µL of this mix
was added to each cDNA sample from aRNA (2 µl: 1 ng
cDNA) and fetal ovary cDNA pool (2 µl =11 ng). Follow-
lowing a brief vortex and centrifugation, the plate was
transferred to a thermal cycler and subjected to the fol-
lowing thermal protocol: 95°C for 10 min; 14/17 cycles
of (95°C for 15 s; 60°C for 4 min); 4°C hold. *RLP19* gene
primers were not included in the STA Primer Mix. The
resulting STA cDNA samples were then treated with
Exonuclease I to digest the primers. A mix was prepared
containing 24 µL 20 units/µL Exonuclease I (M02935,
NEB), 12 µL 10 x Exonuclease I Reaction Buffer and
84 µL H₂O. Two µL of this mix was added to 5 µl of
each STA cDNA sample. Following a brief vortex and
centrifugation, the plate was transferred to a thermal cy-
cler and subjected to the following thermal protocol: 37°C
for 30 min; 80°C for 15 min; 4°C hold. 3.5 µl of each STA
cDNA sample were diluted by adding 9 µL buffer consist-
ing of 10 mM Tris–HCl, pH 8.0; 0.1 mM EDTA (1/5 dilu-
tion) and stored at –20°C.

For the determination of the PCR amplification effi-
ciency specific to each gene, the STA cDNA sample of
the fetal ovary cDNA pool (11 ng cDNA) was serial di-
luted (1, 1:3; 1:3; 1:2; 1:2).

In order to load the chip, a sample pre-mix and primer
pre-mix were prepared. The sample pre-mix consisted of
440 µL 2X Taqman Gene Expression Master Mix (ref
4369510, Applied Biosystem), 44 µL 20X DNA Binding
Dye Sample Loading Reagent (ref 100–0388, Applied Bio-
system), 44 µl 20X EvaGreen (ref BI1790, Interchim) and
132 µL 1X TE. A 6 µL aliquot of this mix was dispensed to
each well of a 96-well assay plate. A 2 µL aliquot of diluted
STA cDNA samples (14 cycle STA cDNA samples and
17 cycle STA cDNA samples) was added to each well and
the assay plate was briefly vortexed and centrifuged. The
primer pre-mix consisted of 440 µl of 2X assay loading re-
agent and 110 µl of 1X TE and 6 µl of this mix was dis-
pensed to each well of a 96-well primer plate. Aliquots
(2 µl, 20X) of individual primer assay were added to each
well and the primer plate was briefly vortexed and centri-
fuged. Following priming of the IFC in the IFC Controller
HX, 5 µL from each well of the 96-well assay plate (cDNA
sample + reagent mix) were dispensed to each sample inlet
of 96.96 IFC and 5 µL of each well of the 96-well primer
plate were dispensed to each detector inlet of the 96.96
IFC. After loading, the chip was transferred to the BioMark
HD and PCR was performed using the thermal protocol:
UNG and Hot Start at 50°C for 120 s; 95°C for 600 s, PCR
Cycles: 35 cycles of (95°C for 15 s; 60°C for 60 s) followed
by a melting phase. Data was analyzed using Fluidigm
Digital PCR Analysis software using the Linear (Derivative)
Baseline Correction Method, the User (Global) Ct Threshold
Method with the threshold set at 0.01, and a Ct Range
of 12 to 25 cycles.
After determination of the threshold cycle (Ct) for each LCM-derived aRNA sample, the PFAFFL method was used to calculate the relative expression of each gene [73]. The relative expression was normalized by the corresponding geometric average of three reference genes using geNorm v3.4 [74]: *β-actin*, *TMED4* and *RPL19* genes, which were respectively slightly, moderately and highly expressed and not regulated during follicle development or between follicular compartments.

The significance of compartment-specific and bio-
markers differential gene expression was evaluated using
Student's test after fourth square transformation.

The significance of the differential expression of genes involved in canonical pathways was tested using the one-way ANOVA model in the R statistical software system (the Comprehensive R Archive National, http://www.cran.r-project.org) after fourth square transformation of the data. For each gene, an ANOVA model was fitted, with the 2 factor stage (4 levels) and compartment (2 levels) with interactions. A backward variable selection procedure was applied as previously described [12].

Data access

All the raw RNA-seq data have been deposited in [EMBL-
EBI ArrayExpress http://www.ebi.ac.uk/arrayexpress/ under
accession number E-MTAB-1587].

**Additional files**

**Additional file 1: Supplemental Figures.** This section includes 9 Supplemental Figures that illustrate the RNAseq design, bioinformatics processes, validation, and signaling pathways.

**Additional file 2: Supplemental Results and Discussion.** This section provides a detailed description of LCM, RNA-seq, expression level preservation and reproducibility results followed by a comparative study of gene expression in 5 previous studies.

**Additional file 3: Differentially expressed genes (DEG) between compartments.**

**Additional file 4: Gene expression of known genes.** \*: pval < 0.05; \*: pval < 0.01; \*:*: pval < 0.005.

**Additional file 5: Compartment specific gene analysis.**

**Additional file 6: Oocyte functional enrichment.**

**Additional file 7: CC functional enrichment.**

**Additional file 8: Transcriptional regulators.** Ingenuity identified the cascade of upstream transcriptional regulators that explain the observed changes in gene expression in our data set. Pval is based on a significant overlap between genes in our data set and known target regulated and transcriptional regulators. The activator z-score infers the activation state to the transcriptional regulator. Without identified bias, this z-score can always be used as an independent test to select downstream regulators. Bias means activation z-score and pval must be used to select downstream regulators. The fold change column lists FC in the data set: fold change =1 corresponds to a differential expression between compartments (>2); fold change ≠ 1 refers to the differential fold change during follicular development.

**Additional file 9: Upstream regulator analysis.**

**Additional file 10: Primer sequences for real-time PCR.**

Competing interests

The authors declare that they have no competing interests.

Authors' contributions

AB conceived and designed the study, performed the LCM experiments, analyzed the RNA-seq data, carried out the RT-PCR and statistical analyses, and was responsible for much of the writing. OB, NM and JS performed RNA sequencing. CC and SF carried out the RNA-seq bioinformatics analysis. FB and FW helped obtain embedded ovary and staining sections. PM and BMP contributed to the conception and design of the study, and were involved in the drafting and revising the manuscript. All authors have read and approved the final version of the manuscript.

Acknowledgements

We thank the Langlade INRA experimental farm that provided the lambs. LCM was processed at FR3450 set, TRI platform Toulouse, Genopole (http://tri.genotoul.fr/index.php?id=36). RNA sequencing was processed at the genomic platform Toulouse, Genopole (http://genomique.genotoul.fr/). We thank Katia Feve for her help with Fluidigm experiment. We thank our colleagues at the Biochips platform Toulouse Genopole (http://biopuce.insa-toulouse.fr/) for their help on this work.

This work was supported by an INRA "Bioressources 2010" grant.

Author details

¹INRA, UMR444 Génétique Cellulaire, Auzeville, BP52627, F-31326 Castanet-Tolosan, France. ²ENVT, UMR444, Laboratoire de Génétique Cellulaire, F-31326 Castanet-Tolosan, France. ³INRA, SIGENAE, UR83 Recherches Avicoles, F-37380 Nouzilly, France. ⁴GeT-PlaGe, Genotoul, INRA Auzeville, F31326 Castanet-tolosan, France. ⁵Université de Toulouse; INSA, UPS, INP; LISBP, 135 Avenue de Rangueil, F-31077 Toulouse, France. ⁶INRA, UMR792, Ingénierie des Systèmes Biologiques et des Procédés, F-31400 Toulouse, France. ⁷CNRS, UMR5504, F-31400 Toulouse, France. ⁸INRA, UMR1198 Biologie du Développement et de la Reproduction, F-78350 Jouy-en-Josas, France.

Received: 6 June 2013 Accepted: 10 December 2013

Published: 19 December 2013

**References**

1. McNatty KP, Smith P, Hudson NL, Heath DA, Tisdall DJ, WS O, Braw-Tal R: Development of the sheep ovary during fetal and early neonatal life and the effect of fecundity genes. *J Reprod Fertil Suppl* 1995, **49**:123-135.

2. Choi Y, Rajkovic A: Genetics of early mammalian folliculogenesis. *Cell Mol Life Sci* 2006, **63**(5):579-590.

3. Fabre S, Pierre A, Mulsant P, Bodin L, Di Pasquale E, Persani L, Monget P, Monniaux D: Regulation of ovulation rate in mammals: contribution of sheep genetic models. *Reprod Biol Endocrinol* 2006, **4**:20.

4. Monget P, Bobe J, Gougeon A, Fabre S, Monniaux D, Dalbies-Tran R: The ovarian reserve in mammals: a functional and evolutionary perspective. *Mol Cell Endocrinol* 2012, **356**(1-2):2-12.

5. Su Y-Q, Sugiura K, Eppig JJ: Mouse oocyte control of granulosa cell development and function: paracrine regulation of cumulus cell metabolism. *Semin Reprod Med* 2009, **27**(1):32-42.

6. Binelli M, Murphy BD: Coordinated regulation of follicle development by germ and somatic cells. *Reprod Fertil Dev* 2010, **22**(1):1-12.

7. Gittens JE, Barr KJ, Vanderhyden BC, Kidder GM: Interplay between paracrine signaling and gap junctional communication in ovarian follicles. *J Cell Sci* 2005, **118**(Pt 1):113-122.

8. Pangas SA, Choi Y, Ballow DJ, Zhao Y, Westphal H, Matzuk MM, Rajkovic A: Oogenesis requires germ cell-specific transcriptional regulators Sohlh1 and Lhx8. *Proc Natl Acad Sci USA* 2006, **103**(21):8090-8095.

9. Yoon SJ, Kim KH, Chung HM, Choi DH, Lee WS, Cha KY, Lee KA: Gene expression profiling of early follicular development in primordial, primary, and secondary follicles. *Fertil Steril* 2006, **85**(1):193-203.

10. Markholt S, Grøndahl ML, Ernst EH, Andersen CY, Ernst E, Lykke-Hartmann K: Global gene analysis of oocytes from early stages in human folliculogenesis shows high expression of novel genes in reproduction. *Mol Hum Reprod* 2012, **18**(2):96-110.

11. Grøndahl ML, Borup R, Vikeså J, Ernst E, Andersen CY, Lykke-Hartmann K: The dormant and the fully competent oocyte: comparing the transcriptome of human oocytes from primordial follicles and in metaphase II. *Mol Hum Reprod* 2013, **19**(9):600-617.
12. Bonnet A, Bevilacqua C, Benne F, Bodin L, Cotinot C, Liaubet L, Sancristobal M, Sarry J, Terenina E, Martin P, et al: Transcriptome profiling of sheep granulosa cells and oocytes during early follicular development obtained by laser capture microdissection. *BMC Genomics* 2011, **12**:417.

13. Anders S, Huber W: Differential expression analysis for sequence count data. *Genome Biol* 2010, 11(10):R106.

14. Pangas SA, Rajkovic A: Transcriptional regulation of early oogenesis: in search of masters. *Hum Reprod Update* 2006, 12(1):65–76.

15. Ramskold D, Wang ET, Burge CB, Sandberg R: An abundance of ubiquitously expressed genes revealed by tissue transcriptome sequence data. *Plos Computational Biology* 2009, 5(12):e1000598.

16. Denomme MM, White CR, Gillio-Meina C, Macdonald WA, Deroo BJ, Kidder GM, Mann MR: Compromised fertility disrupts Peg1 but not Snrpn and Peg3 imprinted methylation acquisition in mouse oocytes. *Front Genet* 2012, 3:129.

17. Lundy T, Smith P, O'Connell A, Hudson NL, McNatty KP: Populations of granulosa cells in small follicles of the sheep ovary. *J Reprod Fertil* 1999, 115(2):251–262.

18. Keeney S: Spo11 and the formation of DNA double-strand breaks in meiosis. *Genome Dyn Stab* 2008, 2:81–123.

19. Guardavaccaro D, Kudo Y, Boulaire J, Barchi M, Busino L, Donzelli M, Margottin-Goguet F, Jackson PK, Yamasaki L, Pagano M: Control of meiotic and mitotic progression by the F box protein beta-Trcp1 in vivo. *Dev Cell* 2003, 4(6):799–812.

20. Li BA, Nair M, Mackay DR, Bilanchone V, Hu M, Fallahi M, Song HQ, Dai Q, Cohen PE, Dai X: Ovo1v1 regulates melanic pachytene progression during spermatogenesis by repressing Id2 expression. *Development* 2005, 132(6):1463–1473.

21. Mandon-Pepin B, Oustry-Vaiman A, Vigier B, Piumi F, Cribiu E, Cotinot C: Expression profiles and chromosomal localization of genes controlling meiosis and follicular development in the sheep ovary. *Biol Reprod* 2003, 68(3):985–995.

22. Yan W, Rajkovic A, Viveiros MM, Burns KH, Eppig JJ, Matzuk MM: Identification of Gasz, an evolutionarily conserved gene expressed exclusively in germ cells and encoding a protein with four ankyrin repeats, a sterile-alpha motif, and a basic leucine zipper. *Mol Endocrinol* 2002, 16(6):1168–1184.

23. Logan KA, Juengel JL, McNatty KP: Onset of steroidogenic enzyme gene expression during ovarian follicular development in sheep. *Biol Reprod* 2002, 66(4):906–916.

24. Parker KL, Rice DA, Lala DS, Ikeda Y, Luo XR, Wong M, Bakke M, Zhao LP, Frigeri C, Hanley NA, et al: Steroidogenic factor 1: An essential mediator of endocrine development. *Recent Progr Horm Res Reprod Horm Human Health* 2002, 57:19–36.

25. Pelusi C, Ikeda Y, Zubair M, Parker KL: Impaired follicle development and infertility in female mice lacking steroidogenic factor 1 in ovarian granulosa cells. *Biol Reprod* 2008, 79(6):1074–1083.

26. Chen Y, Breen K, Pepling ME: Estrogen can signal through multiple pathways to regulate oocyte cyst breakdown and primordial follicle assembly in the neonatal mouse ovary. *J Endocrinol* 2009, 202(3):407–417.

27. Sen A, Hammes SR: Granulosa cell-specific androgen receptors Are critical regulators of ovarian development and function. *Mol Endocrinol* 2010, 24(7):1393–1403.

28. Saunders PTK, Millar MR, Williams K, Macpherson S, Harkiss D, Anderson RA, Orr B, Groome NP, Scobie G, Fraser HM: Differential expression of estrogen receptor-alpha and -beta and androgen receptor in the ovaries of marmosets and humans. *Biol Reprod* 2000, 63(4):1098–1105.

29. Beltowski J, Semczuk A: Liver X receptor (LXR) and the reproductive system - a potential novel target for therapeutic intervention. *Pharmacol Rep* 2010, 62(1):15–27.

30. Koutsounas I, Giaginis C, Theocharis S: Farnesoid X Receptor (FXR) from normal to malignant state. *Histol Histopathol* 2012, 27(7):835–853.

31. Rucker EB, Dierisseau P, Wagner KU, Garrett L, Wynshaw-Boris A, Flaws JA, Hennighausen L: Bcl-x and bax regulate mouse primordial germ cell survival and apoptosis during embryogenesis. *Mol Endocrinol* 2000, 14(7):1038–1052.

32. Boumela I, Assou S, Aouacheria A, Haozuï D, Dechaud H, De Vos J, Handyside A, Hamamah S: Involvement of BCL2 family members in the regulation of human oocyte and early embryo survival and death: gene expression and beyond. *Reproduction* 2011, 141(5):549–561.

33. Hsu SY, Kaipia A, McGee E, Lomeli M, Hsueh AJW: Bok is a pro-apoptotic Bcl-2 protein with restricted expression in reproductive tissues and heterodimerizes with selective anti-apoptotic Bcl-2 family members. *Proc Natl Acad Sci USA* 1997, **94**(23):12401–12406.

34. Russell MC, Cowan RG, Harman RM, Walker AL, Quirk SM: The hedgehog signaling pathway in the mouse ovary. *Biol Reprod* 2007, 77(2):226–236.

35. Schlessinger K, Hall A, Tolwinski N: Wnt signaling pathways meet Rho GTPases. *Genes Dev* 2009, 23(3):265–277.

36. Artavanis-Tsakonas S, Rand MD, Lake RJ: Notch signaling: cell fate control and signal integration in development. *Science* 1999, 284(5415):770–776.

37. Xu J, Gridley T: Notch2 is required in somatic cells for breakdown of ovarian germ-cell nests and formation of primordial follicles. *BMC Biol* 2013, 11(1):13.

38. Zhang CP, Yang JL, Zhang J, Li L, Huang L, Ji SY, Hu ZY, Gao F, Liu YX: Notch signaling is involved in ovarian follicle development by regulating granulosa cell proliferation. *Endocrinology* 2011, 152(6):2437–2447.

39. Yang MY, Fortune JE: Vascular endothelial growth factor stimulates the primary to secondary follicle transition in bovine follicles in vitro. *Mol Reprod Dev* 2007, **74**(9):1095–1104.

40. Greenaway J, Connor K, Pedersen HG, Coomber BL, Lamarre J, Petrik J: Vascular endothelial growth factor and its receptor, Flk-1/KDR, are cytoprotective in the extravascular compartment of the ovarian follicle. *Endocrinology* 2004, **145**(6):2896–2905.

41. Mazerbourg S, Overgaard MT, Otvig C, Christiansen M, Conover CA, Laurendeau I, Vidaud M, Tosser-Klopp G, Zapf J, Monget P: Pregnancy-associated plasma protein-A (PAPP-A) in ovine, bovine, porcine, and equine ovarian follicles: involvement in IGF binding protein-4 proteolytic degradation and mRNA expression during follicular development. *Endocrinology* 2001, **142**(12):5243–5253.

42. Mazerbourg S, Bondy CA, Zhou J, Monget P: The insulin-like growth factor system: a key determinant role in the growth and selection of ovarian follicles? A comparative species study. *Reprod Domest Anim* 2003, **38**(4):247–258.

43. Zhao J, Taverne MA, Van Der Weijden GC, Bevers MM, Van Den Hurk R: Insulin-like growth factor-I (IGF-I) stimulates the development of cultured rat pre-antral follicles. *Mol Reprod Dev* 2001, **58**(3):287–296.

44. Perks CM, Denning-Kendall PA, Gilmour RS, Wathes DC: Localization of messenger ribonucleic acids for insulin-like growth factor I (IGF-I), IGF-II, and the type 1 IGF receptor in the ovine ovary throughout the estrous cycle. *Endocrinology* 1995, **136**(12):5266–5273.

45. Wandji SA, Pelletier G, Sirard MA: Ontogeny and cellular localization of 125I-labeled insulin-like growth factor-I, 125I-labeled follicle-stimulating hormone, and 125I-labeled human chorionic gonadotropin binding sites in ovaries from bovine fetuses and neonatal calves. *Biol Reprod* 1992, **47**(5):814–822.

46. Pitetti JL, Torre D, Conne B, Papaioannou MD, Cederroth CR, Xuan S, Kahn R, Parada LF, Vassalli JD, Efstratiadis A, et al: Insulin receptor and IGF1R are not required for oocyte growth, differentiation, and maturation in mice. *Sex Dev* 2009, **3**(5):264–272.

47. Kiapekou E, Loutradis D, Drakakis P, Zapanti E, Mastorakos G, Antsaklis A: Effects of GH and IGF-I on the in vitro maturation of mouse oocytes. *Hormones (Athens)* 2005, **4**(3):155–160.

48. Mani AM, Fenwick MA, Cheng Z, Sharma MK, Singh D, Wathes DC: IGF1 induces up-regulation of steroidogenic and apoptotic regulatory genes via activation of phosphatidylinositol-dependent kinase/AKT in bovine granulosa cells. *Reproduction* 2010, **139**(1):139–151.

49. Mochida N, Akatani-Hasegawa A, Saka K, Ogino M, Hosoda Y, Wada R, Sawai H, Shibahara H: Live births from isolated primary/early secondary follicles following a multistep culture without organ culture in mice. *Reproduction* 2013, **146**(1):37–47.

50. Gilchrist RB, Lane M, Thompson JG: Oocyte-secreted factors: regulators of cumulus cell function and oocyte quality. *Hum Reprod Update* 2008, **14**(2):159–177.

51. Matos MH, van den Hurk R, Lima-Verde IB, Luque MC, Santos KD, Martins FS, Bao SN, Lucci CM, Figueiredo JR: Effects of fibroblast growth factor-2 on the in vitro culture of caprine preantral follicles. *Cells Tissues Organs* 2007, **186**(2):112–120.

52. Nilsson E, Parrott JA, Skinner MK: Basic fibroblast growth factor induces primordial follicle development and initiates folliculogenesis. *Mol Cell Endocrinol* 2001, **175**(1-2):123–130.

53. Sun YL, Zeng S, Ye K, Yang C, Li MH, Huang BF, Sun LN, Zhou LY, Wang DS: Involvement of FGF9/16/20 subfamily in female germ cell development of the Nile tilapia, Oreochromis niloticus. *Fish Physiol Biochem* 2012, **38**(5):1427–1439.

54. Törmälä RM, Jääskeläinen M, Lakkakorpi J, Liakka A, Tapanainen JS, Vaskivuo TE: Zona pellucida components are present in human fetal ovary before follicle formation. *Mol Cell Endocrinol* 2008, **289**(1-2):10–15.
55. Rajkovic A, Pangas SA, Ballow D, Suzumori N, Matzuk MM: NOBOX deficiency disrupts early folliculogenesis and oocyte-specific gene expression. *Science* 2004, **305**(5687):1157–1159.

56. Robles R, Morita Y, Mann KK, Perez GI, Yang S, Matikainen T, Sherr DH, Tilly JL: The aryl hydrocarbon receptor, a basic helix-loop-helix transcription factor of the PAS gene family, is required for normal ovarian germ cell dynamics in the mouse. *Endocrinology* 2000, **141**(1):450–453.

57. Gronostajski RM: Roles of the NF1/CTF gene family in transcription and development. *Gene* 2000, **249**(1–2):31–45.

58. Zagórska A, Dulak J: HIF-1: the knowns and unknowns of hypoxia sensing. *Acta Biochim Pol* 2004, **51**(3):563–585.

59. Sun QY, Liu K, Kikuchi O: Oocyte-specific knockout: a novel in vivo approach for studying gene functions during folliculogenesis, oocyte maturation, fertilization, and embryogenesis. *Biol Reprod* 2008, **79**(6):1014–1020.

60. Patil AA, Cai Y, Sang Y, Blecha F, Zhang G: Cross-species analysis of the mammalian beta-defensin gene family: presence of syntenic gene clusters and preferential expression in the male reproductive tract. *Physiol Genomics* 2005, **23**(1):5–17.

61. Canty-Laird E, Carré GA, Mandon-Pépin B, Kadler KE, Fabre S: First evidence of bone morphogenetic protein 1 expression and activity in sheep ovarian follicles. *Biol Reprod* 2010, **83**(1):138–146.

62. Soues S, Kann ML, Fouquet JP, Melki R: The cytosolic chaperonin CCT associates to cytoplasmic microtubular structures during mammalian spermiogenesis and to heterochromatin in germline and somatic cells. *Exp Cell Res* 2003, **288**(2):363–373.

63. Buanne P, Corrente G, Micheli L, Palena A, Lavia P, Spadafora C, Lakshmana MK, Rinaldi A, Banfi S, Quarto M, et al: Cloning of PC3B, a novel member of the PC3/BTG/TOB family of growth inhibitory genes, highly expressed in the olfactory epithelium. *Genomics* 2000, **68**(3):253–263.

64. Marnef A, Maldonado M, Bugaut A, Balasubramanian S, Kress M, Weil D, Standart N: Distinct functions of maternal and somatic *Pat1* protein paralogs. *RNA* 2010, **16**(11):2094–2107.

65. Nakamura Y, Tanaka KJ, Miyauchi M, Huang L, Tsujimoto M, Matsumoto K: Translational repression by the oocyte-specific protein P100 in Xenopus. *Dev Biol* 2010, **344**(1):272–283.

66. Kim N, Burden SJ: MuSK controls where motor axons grow and form synapses. *Nat Neurosci* 2008, **11**(1):19–27.

67. Mariette J, Escudié F, Allias N, Salin G, Noiret C, Thomas S, Klopp C: NG6: Integrated next generation sequencing storage and processing environment. *BMC Genomics* 2012, **13**:462.

68. Li P, Ponnala L, Gandotra N, Wang L, Si Y, Tausta SL, Kebrom TH, Provart N, Patel R, Myers CR, et al: The developmental dynamics of the maize leaf transcriptome. *Nat Genet* 2010, **42**(12):1060–1067.

69. Archibald AL, Cockett NE, Dalrymple BP, Faraut T, Kijas JW, Maddox JF, McEwan JC, Hutton Oddy V, Raadsma HW, Wade C, et al: The sheep genome reference sequence: a work in progress. *Anim Genet* 2010, **41**(5):449–453.

70. Jager M, Ott CE, Grunhagen J, Hecht J, Schell H, Mundlos S, Duda GN, Robinson PN, Lienau J: Composite transcriptome assembly of RNA-seq data in a sheep model for delayed bone healing. *BMC Genomics* 2011, **12**:158.

71. Zimin AV, Delcher AL, Florea L, Kelley DR, Schatz MC, Puiu D, Hanrahan F, Pertea G, Van Tassell CP, Sonstegard TS, et al: A whole-genome assembly of the domestic cow. *Bos taurus*. *Genome Biol* 2009, **10**(4):R42.

72. Elsik CG, Tellam RL, Worley KC, Gibbs RA, Muzny DM, Weinstock GM, Adelson DL, Eichler EE, Elnitski L, Guigó R, et al: The genome sequence of taurine cattle: a window to ruminant biology and evolution. *Science* 2009, **324**(5926):522–528.

73. Pfaffl MW: A new mathematical model for relative quantification in real-time RT-PCR. *Nucleic Acids Res* 2001, **29**(9):e45.

74. Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, Speleman F: Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. *Genome Biol* 2002, **3**(7). RESEARCH0034. Epub 2002 Jun 18.

doi:10.1186/1471-2164-14-904

**Cite this article as:** Bonnet *et al*. An overview of gene expression dynamics during early ovarian folliculogenesis: specificity of follicular compartments and bi-directional dialog. *BMC Genomics* 2013 **14**:904.

Submit your next manuscript to BioMed Central and take full advantage of:

• Convenient online submission

• Thorough peer review

• No space constraints or color figure charges

• Immediate publication on acceptance

• Inclusion in PubMed, CAS, Scopus and Google Scholar

• Research which is freely available for redistribution

Submit your manuscript at
www.biomedcentral.com/submit

BioMed Central