Expression profiling and function analysis identified new microRNAs regulating cumulus expansion and apoptosis in cumulus cells

Abstract

Although microRNAs (miRNAs) expressed in cumulus cells (CCs) may be used to select competent oocytes/embryos, only a limited number of such miRNAs has been reported. To identify more miRNAs that regulate cumulus expansion (CE) and CC apoptosis, we first established that mouse cumulus–oocyte complexes (COCs) cultured in expansion-supporting medium supported full CE while undergoing mild apoptosis, whereas mouse oocytectomized COCs (OOXs) cultured in apoptosis-triggering medium underwent severe apoptosis while supporting no CE. RNA- and miRNA-sequencing and bioinformatics using CCs from these cultured COCs/OOXs identified candidate apoptosis- and/or CE-regulating miRNAs. Transfection of COCs/OOXs with miRNA mimic or inhibitor validated that miR-212-5p and 149-5p promoted CE by facilitating Has2 expression; miR-31-5p and 27a-3p promoted CE by increasing both Has2 and Ptx3 expression; and miR-351-5p and 503-5p inhibited CE by suppressing Ptx3 expression. Furthermore, miR-212-5p, 149-5p and Nov798 inhibited CC apoptosis, involving both Bcl2/Bax and Fas signaling. Analysis using in vivo matured COCs further verified the above apoptosis- and/or CE-regulating miRNAs, except for miR-149-5p. In conclusion, this study identified and validated new CE- and apoptosis-regulating miRNAs in CCs, which could be used as biomarkers to select competent oocytes/embryos and for elucidating how the oocyte-derived factors regulate CE and CC apoptosis.

Introduction

In the antral follicles of mammalian ovaries, there are two types of anatomically and functionally different granulosa cells. Whereas the mural granulosa cells line the wall of the follicle, the cumulus cells (CCs) enclose the oocyte forming a functional unit, the cumulus–oocyte complex (COC). In vivo, CCs undergo apoptosis with the advancement of atresia of the antral follicles, and in vitro, a negative correlation has been observed between apoptosis of CCs and the developmental potential of oocytes (Bosco et al., 2017; Almeida et al., 2018). Furthermore, in the preovulatory follicles, the cumulus oophorus expands dramatically in response to the preovulatory gonadotrophin surge, and it has been demonstrated that an optimal cumulus expansion (CE) is essential for normal oocyte maturation, ovulation, fertilization and embryo development (Yokoo and Sato, 2004).

Although IVM can provide large numbers of mature oocytes for human and animal ART, the developmental competence of IVM oocytes is significantly lower than that of the in vivo matured oocytes (Trounson et al., 1994; Barnes et al., 1996; Cha et al., 2000). Thus, the optimization of IVM protocols is of great importance for human and animal ART. Selection of competent oocytes for IVM or competent embryos for transfer is crucial for optimization of IVM protocols. For example, because the most frequent source of ovaries for animal IVM is the slaughterhouse, and many factors, such as age, the stage of the estrous cycle, nutritional status, genetic potential and the presence of a reproductive disorder in donors, can influence oocyte quality, it is almost impossible to avoid the retrieval of a heterogeneous population of oocytes with distinct developmental competences (Aguila et al., 2020). Furthermore, selecting embryos with the highest implantation potential is crucial for successful human IVF, as single embryo transfer is becoming increasingly common to avoid multiple pregnancies (Assou et al., 2010).

Although the current oocyte quality assessment systems, such as CC apoptotic indexes (Bosco et al., 2017) and CE morphology (Moulavi and Hosseini, 2018), have played an important role in the efficiency of IVF, their precision is still insufficient. Thus, it is very important to develop objective, accurate, fast and affordable tests for pre- and/or post-IVM assessment of oocyte competence. Because the oocyte plays a dominant role in regulating CC functions during follicle development, the functions of CCs are believed to indirectly reflect oocyte competence (Gilchrist et al., 2004; Shepel et al., 2016). Therefore, gene expression analysis in CCs should be able to provide a non-invasive method for identification of competent oocytes or embryos. For example, RT-PCR and DNA microarrays have shown that genes expressed in human CCs can be used to predict embryo quality and pregnancy outcomes following IVF (Assou et al., 2010). Furthermore, although it is accepted that oocytes can promote CE by secreting the CE enabling factors (Buccione et al., 1990; Zhang et al., 2008; Li et al., 2012) and protect CCs from apoptosis by releasing growth differentiation factor 9 and/or bone morphogenetic protein 5 (Hussein et al., 2005; Han et al., 2017), the identities of these oocyte-secreted factors and their mechanisms of action are largely undefined.

Being detectable in circulatory bio-fluids and even in culture media, microRNAs (miRNAs) have been proposed as biomarkers for disease diagnosis (Kumar et al., 2017; Matsuzaki and Ochiya, 2017) and embryo selection (Rosenbluth et al., 2014; Capalbo et al., 2016). Expression of miRNAs in CCs was reported in human (Assou et al., 2013), bovine (Abd El Naby et al., 2013), mice (Yao et al., 2014) and pig (Pan et al., 2015). Functional analyses of some of the miRNAs expressed in CCs indicated that they were involved in regulation of CE (Yao et al., 2014; Pan et al., 2015; Li et al., 2017), CC apoptosis (Assou et al., 2013; Han et al., 2017) or both (Cui et al., 2018; Liu et al., 2018). Thus, these miRNAs may be used to select competent oocytes or embryos and to reveal the mechanisms by which the oocyte-secreted factors regulate CE and CC apoptosis. However, the number of miRNAs reported so far is far from being enough to serve these purposes. Furthermore, most of the above studies on miRNAs in CCs, including our previous study (Han et al., 2017), were focusing on a single miRNA in their functional analysis. This can be very challenging and may not result in significant measurable biological changes, as most miRNAs act in clusters and are fine tuners of cellular functions (Maalouf et al., 2016). Thus, profiling studies should be conducted to assess the function of a cluster of miRNAs that are co-expressed and more likely to affect the functions of each other.

The objective of this study was to identify and validate new miRNAs that regulate CE and CC apoptosis. miRNA/mRNA sequencing and functional analysis were conducted using CCs after maturation of mouse COCs and oocytectomized COCs (OOXs) in vitro. The new miRNAs identified can be used as biomarkers for selection of not only competent embryos for transfer but also of competent oocytes for IVM. They are also important for elucidating how the oocyte-derived factors regulate CE and CC apoptosis, which will potentially contribute to the progress of IVM/IVF practice.

Materials and methods

All chemicals and reagents used in this study were purchased from Sigma-Aldrich Corp. (St. Louis, MO, USA) unless otherwise specified.

Experimental design

In this study, miRNA/mRNA sequencing and functional analysis were conducted to identify more miRNAs that regulate CC apoptosis and CE. Mouse oocytectomized COCs (OOXs) and COCs were cultured for 18 h in apoptosis-triggering medium (ATM) and CE-supporting medium (ESM), respectively, before isolation of mRNAs and miRNAs (Fig. 1). The mRNAs and miRNAs isolated were subjected to RNA- and miRNA-seq to identify differentially expressed (DE) genes and miRNAs, respectively. Target genes were predicted from the DE miRNAs and were overlapped with the DE genes identified by RNA-seq to determine the candidate genes. The candidate genes were subjected to Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment and protein–protein interaction (PPI) analysis to obtain the top enriched and top degree genes (genes with the highest degree of connectivity), respectively. The top enriched and top degree genes obtained were then overlapped to identify the candidate miRNAs, from which the would-be apoptosis- and/or expansion-related (A/E) miRNAs were selected according to fold changes (FCs). For functional analysis, COCs and OOXs were transfected with mimics or inhibitors of the A/E miRNAs, and the transfected COCs and OOXs were cultured in ESM and ATM before assessment of CE and apoptosis, to determine the apoptosis-regulating (AR) and expansion-regulating (ER) miRNAs, respectively.

Animals and oocyte collection

Mice (Kunming strain) were raised in a room with 14-h light and 10-h dark cycles, with the dark starting from 20:00. The experimental procedures for mouse care and handling were approved by the Shandong Agricultural University Animal Care and Use Committee, PR China, with permit number SDAUA-2010-004. Female mice were injected i.p. with 10 IU eCG at the age of 8–10 weeks. At 48 h after eCG injection, the mice were killed to recover ovaries for oocyte collection. Only COCs with more than three layers of CCs, containing oocytes larger than 70 µm in diameter and with a homogenous cytoplasm were used for experiments.

To recover in vivo matured COCs from zearalenone (ZEN)-treated mice, some female mice were injected with ZEN (1 mg/kg body weight) daily for 5 days, as described previously (Pan et al., 2021), and age-matched mice that were sham-injected were used for controls. At 24 h after the last ZEN injection, the mice were injected with 10 IU eCG and, 48 h later, they were injected with 10 IU hCG. At 6 h after hCG injection, mice were killed to recover in vivo matured COCs.

Preparation of OOXs

To produce OOXs, the cytoplasm of each oocyte was micro-surgically removed from the COCs using micromanipulators. The OOXs produced consisted of a hollow zona pellucida surrounded by more than three layers of intact CCs.

Culture of COCs and OOXs

COCs and OOXs were cultured for 18 h in microdrops (20–30 COCs or OOXs per 100 µl drop) at 37.5°C under a humidified atmosphere with 5% CO₂ in air. The OOXs were cultured in an ATM medium, and the COCs were cultured in an ESM medium. The ATM was a simplified MEM (s-MEM) alone, and the ESM consisted of s-MEM supplemented with FSH and fetal bovine serum. The s-MEM contained inorganic salts (1.8 mM CaCl₂, 0.81 mM MgSO₄, 5.3 mM KCl, 26.2 mM NaHCO₃, 117.2 mM NaCl, 1.0 mM NaH₂PO₄), 2 mM glutamine, 5.56 mM glucose, 1 mM sodium pyruvate, 4 mg/ml bovine serum albumin, 0.03 mM phenol red, 50 IU/ml penicillin and 50 µg/ml streptomycin.

ELISA for hyaluronic acid in CCs from OOXs or COCs

Following culture of COCs and OOXs for 18 h in ESM and ATM, respectively, 50 COCs or OOXs were placed in a freezing tube containing 60 µl M2 medium. Meanwhile, another group of 50 COCs were pipetted in M2 medium containing 0.1% hyaluronidase to prepare cumulus-denuded oocytes, and the denuded oocytes obtained were placed in the freezing tube containing OOXs. We thus obtained COCs and OOXs samples that contained the same number of oocytes and CCs without hyaluronidase treatment. Then the cells were lysed by the freeze–thaw method to release proteins. The ELISA assays were conducted using a HA ELISA kit (BLUE GENE) as follows: 50 µl standards or samples were added in duplicate to wells of a microtiter plate coated with mouse monoclonal antibodies; 100 µl conjugate were added to each well and incubated at 37°C for 60 min; the microtiter plate was washed using the wash solution and dried using paper towels; and 50 µl substrate A and 50 µl substrate B were added to each well and incubated for 10–15 min at 37°C. Optical density at 450 nm was measured using a plate reader (BioTek-ELx808; BioTek Instruments Inc., Winooski, VT, USA) within 15 min of the reaction termination by application of 50 µl stop solution. The concentrations of hyaluronic acid in the samples were calculated against the standard curves.

Assessment of CC apoptosis by Hoechst staining

The CCs released from 30 OOXs or COCs were dispersed by pipetting in M2 medium with a fine pipette. The dispersed CCs were centrifuged at 200×g for 5 min. Then, the pellets were resuspended in 50 µl M2 containing 0.01 mg/ml of Hoechst 33342 and stained in the dark for 5 min. The stained cells were resuspended in M2 and centrifuged at 200×g twice (5 min each). Finally, a 5 µl drop of suspension was smeared on a slide and examined under a Leica DMLB fluorescence microscope at a magnification of 400×. Six to eight fields were randomly observed on each smear. The apoptotic cells showed pyknotic nuclei full of heterochromatin that was heavily stained and gave bright fluorescence, whereas healthy cells showed normal nuclei with sparse heterochromatin spots. To reduce the subjectivity, percentages of the apoptotic cells were always calculated, double blind, by two investigators.

Sequencing of mRNAs and miRNAs in CCs by illumina HiSeq

Extraction of total RNA

Total RNA was extracted from CC samples using an RNeasy Mini Kit (Qiagen, Valenda, CA, USA) or miRNeasy Mini Kit (Qiagen) for mRNA and miRNA sequencing, respectively. Each CC sample for mRNA sequencing contained CCs from about 300 OOXs/COCs, whereas each sample for miRNA sequencing contained CCs from about 900 OOXs/COCs. On each experimental day, about 120 OOXs/COCs were collected, of which 30 were used for mRNA and 90 for miRNA extraction. The OOXs/COCs recovered were centrifuged (2000×g for 5 min) and stored in liquid nitrogen before mRNA/miRNA extraction.

RNA integrity evaluation

Total RNA extracted from each sample was qualified and quantified by Agilent 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA, USA) and NanoDrop (Thermo Fisher Scientific Inc., Waltham, MA, USA). For library preparation, 1 μg of total RNA with an RNA integrity number (RIN) value above 7 was used for mRNA library preparation, and 2 μg of total RNA with a RIN value above 7.5 was used for miRNA library preparation.

Library preparation and sequencing

Library preparation and sequencing were carried out by Genewiz (Suzhou, China). The next-generation sequencing library preparations were constructed according to the manufacturer’s protocol (NEBNext Ultra RNA Library Prep Kit from Illumina for mRNAs, or Multiplex Small RNA library Prep Set from Illumina for miRNAs). Then, libraries with different indexes were multiplexed and loaded on an Illumina HiSeq instrument according to manufacturer’s instructions (Illumina, San Diego, CA, USA). Sequencing was carried out using 2 × 150 paired-end (PE) configuration. Image analysis and base calling were conducted by the HiSeq Control Software (HCS) + Off-Line Basecaller software (OLB) + GAPipeline-1.6 (Illumina) on the HiSeq instrument. The HCS + OLB + GAPipeline-1.6 (Illumina) can control the flow tank platform, fluid system and flow tank temperature and capture cluster images. It can automatically generate image analysis, base reading and base reading quality reports.

Bioinformatic analysis

Quality control

The files with the nucleotide sequences and the quality description of the bases generated by the Illumina Hiseq X ten system (Illumina, San Diego, CA, USA) were analyzed with the program FastQC v0.10.1 (www.bioinformatics.babraham.ac.uk/projects/fastqc). Approximately 90% of the reads were mapped, with phred-score >30 (Q30 > 85% is considered qualified). Cutadapt (version 1.9.1) (Martin, 2011) and Trimmomatic (version 0.30) (Parra et al., 2007) were used to remove the adapter sequences for the RNA-Seq and miRNA-Seq, respectively.

Mapping

The mapping of generated fragments (reads) was performed with Hisat2 (v2.0.1) (Powell et al., 2016), using Ensembl (mm10) as the reference genome and NCBI Release 81 as gene annotation. The index build command is Code: Hisat2-build Genome Genome. For this analysis the parameters used were: -x <bt2-idx> -1 <m1> -2 <m2> [-S <sam>]. After mapping, the files with the final alignment information of the reads in the ‘bam’ binary format were converted to the ‘sam’ format using SAMtools v0.1.19-44 428cd.

Differential expression analysis

Differential expression analysis for mRNAs was conducted using the Cuffdiff (v 2.2.1) (Trapnell et al., 2012) and that for miRNAs was conducted using the DESeq2 (v 1.6.3) (Love et al., 2014), a model based on the negative binomial distribution. The normalization and differential expression between the groups from the raw counts were performed in the statistical environment R (v 3.1.2) using the DESeq2 (v 1.6.3) package, with Bonferroni P-value adjustment method. In the same statistical environment, principal component analysis was performed among the groups to check the clustering of the samples according to their overall expression. The DE genes had values of FC >2 and adjusted q-value <0.05, and the DE miRNAs had values of FC >1.5 and adjusted P-value <0.05.

Cluster analysis of DE mRNAs and miRNAs

The FPKM (fragments per kilo bases per million reads) and TPM (transcripts per million) values were used to represent the level of mRNA and miRNA expression, respectively, in CCs from different treatments, to carry out the hierarchical clustering analysis. The normalization of miRNA expression levels between OOX and COC was carried out using the normalization formula: microRNA read count ×10⁶/total reads; and the normalization of mRNA expression was conducted using the normalization formula: Total exon fragments/mapped reads (millions) × exon length (KB). Different cluster grouping information was expressed as different colors. Heat maps of DE mRNAs and miRNAs were drawn using the MeV 4.9.0 software (Saeed et al., 2003).

RT-qPCR validation of selected DE mRNAs and miRNAs

For mRNA assay, CCs from 100 OOXs or COCs were isolated for total RNA isolation. Total RNA isolation was performed using an RNAqueous-Micro Total RNA Isolation Kit (AM1931, Ambion, Austin, TX, USA). The RNA isolated was resuspended in dH₂O treated with diethylpyrocarbonate (DEPC-dH₂O). Reverse transcription was performed in a total volume of 20 µl using transcriptor reverse transcriptase (03531287001 or 03335399001, Roche, Basel, Switzerland). Briefly, 2 µl of each RNA sample was mixed in a 0.2 ml reaction tube with 1 µl Oligo dT18 (SO132 or R0192, Fermentas, St. Leon-Rot, Germany), and 10 µl of DEPC-dH₂O, and the mixture was incubated in a PCR instrument at 65°C for 10 min. At the end of incubation, the reaction tube was cooled on ice for 2 min and centrifuged (200×g, 4°C) for 10 s. Then, 4 µl of 5× RT buffer, 0.5 µl RNase inhibitor (03335402001, Roche), 2 µl dNTP (R0192, Fermentas) and 0.5 µl transcriptor reverse transcriptase were added to the reaction tube. The mixture was then incubated at 55°C for 30 min, followed by incubation at 85°C for 5 min before storage at −20°C until use.

Gene-specific primers for real-time qPCR are listed in Table I. Quantification of mRNA was conducted using the Mx3005P real-time PCR instrument (Stratagene, Valencia, CA, USA). Amplification reactions were performed in a 10 µl reaction volume containing 1 µl of cDNA, 5 µl of 2× SYBR Green Master Mix (600882, Agilent Technologies, Palo Alto, CA, USA), 0.15 µl of ROX (reference dye), 3.25 µl of DEPC-dH₂O and 0.3 µl each of forward and reverse gene-specific primers (10 µM). Cycle amplification conditions comprised an initial denaturation step at 95°C for 3 min followed by 40 cycles at 95°C for 20 s and 60°C for 20 s. Immediately after amplification, PCR products were analyzed by sequencing, dissociation curve analysis and gel electrophoresis to determine specificity of the reaction. Gene expression was normalized to the Gapdh internal control. All values were then expressed relative to the calibrator samples using the 2^–(ΔΔCT) method.

For miRNA assay, CCs from 240 OOXs or COCs were used for total RNA isolation, which was performed using an Ambion mirVana miRNA Isolation Kit (AM1561, Ambion, Austin, TX, USA). The RNA isolated was resuspended in DEPC-dH₂O. Reverse transcription was performed using a Takara Bio Reverse-transcription kit (638313, Takara Bio, Dalian, China). Briefly, 5 µl 2× mRQ buffer, 3.75 µl total RNA (0.25–8 µg/ml), and 1.25 µl mRQ enzyme were mixed and incubated at 65°C for 1 h followed by heating at 85°C for 5 min. The cDNA obtained were then diluted to 50 µl for use in the following real-time qPCR experiment.

Reactions for real-time qPCR of miRNAs were performed using a Takara Bio RT-qPCR kit (RR820A, Takara Bio, Dalian, China). A 10-µl reaction mixture was prepared according to the manufacturer’s instructions, which contained 5 µl of 2× TB Green Premix Ex TaqII, 0.4 µl each of forward and reverse primers, 0.2 µl of 50× ROX Reference Dye, 1 µl of cDNA and 3 µl H₂O. The sequences of primers used are listed in Table I. Reactions were run in Mx3005P (Stratagene) as follows: initial denaturation at 95°C for 30 s followed by 40 cycles of 95°C for 5 s and 60°C for 34 s. The ROX fluorescence was collected after each cycle. Relative expression of miRNAs was calculated by MxPro software with 2^–(ΔΔCT) and U6 as normalizer.

MicroRNA target mRNA analysis

The software programs miRanda (version 3.3a, https://anaconda.org/bioconda/miranda), TargetScan (http://www.targetscan.org/mamm_31/) and RNAhybride (version 2.1.2, https://bibiserv.cebitec.uni-bielefeld.de/rnahybrid) were used to predict the targets of miRNAs. Only genes that were predicted positive by two or three software programs were selected as the target genes.

KEGG pathway enrichment and PPI analyses

Candidate genes generated from overlap between DE genes and target genes (Fig. 1) were subjected to KEGG (http://en.wikipedia.org/wiki/KEGG) or PPI (String Database, http://string-db.org/) analyses. For KEGG analysis, the software of DAVID Bioinformatics Resources 6.8 was used and P < 0.05 was considered as significantly enriched. To generate the PPI network on the STRING database, we selected the ‘none/query proteins only’ option in the 1st shell, and selected the ‘none’ option in the 2nd shell. Only PPIs with a medium confidence score of 0.4 were considered, and the software Cytoscape v3.6.1 (https://www.cytoscape.org/) was used to visualize the protein interaction network. The hub proteins were identified based on the node degree that represents the number of interactions (edges) the node has, and only proteins with a node degree of ≥20 were selected as top degree genes.

Transfection of COCs or OOXs with miRNA inhibitor or mimic

Knockdown or overexpression of miRNAs was realized by transfection of COCs or OOXs with miRNA inhibitors or mimics, respectively. The miRNA mimics, inhibitors and their negative controls (Table II) were synthesized by Ribobio Co., Ltd, Guangzhou, Guangdong, China. The oligonucleotides of miRNA mimics and inhibitors were transfected at a final concentration of 50 and 100 nM, respectively. To prepare a 100 mM stock solution, the oligonucleotides were dissolved in DEPC-dH₂O, and the stock solution was stored at −20°C before use. miRNA mimics or inhibitors or their controls were transfected by adding to each culture well 0.75 µl of Lipofectamine RNAiMAX Reagent (13778500, Invitrogen, Carlsbad, CA, USA), 25 µl Opti-MEM I Reduced Serum Medium (51985042, Life technologies, Gaithersburg, MD, USA), 10 µl of miRNA-lipid complex and 90 µl ESM or ATM, along with 25–30 COCs or OOXs. After 18 h of transfection culture, the CCs cultured in ATM were used for assessment of apoptosis, and the CCs cultured in ESM were used to examine CE indices.

Data analysis

Unless otherwise specified, each treatment was repeated three times. Percentage data were arc sine transformed before being analyzed using independent-sample Student’s t-test. The Statistics Package for Social Science (SPSS 20, SPSS, Inc., Chicago, IL, USA) was used to conduct the analysis. Data were expressed as mean ± SEM and were considered significant when the P-value was <0.05.

Results

The establishment of CC models with significantly different levels of CE and apoptosis

The OOXs were cultured in ATM, and COCs were cultured in ESM for 18 h before examination for CE and apoptosis. CE in individual COCs/OOXs was classified into 0 to 4 grades. 0 indicates no expansion with all CCs forming a compact mass; 1 indicates expansion only in the peripheral two layers of CCs; 2 indicates expansion involving several CC layers; 3 indicates expansion involving all CC layers but corona radiata cells; and 4 indicates expansion of the whole cumulus including corona radiata cells (Fig. 2A and B). None of the OOXs expanded, but 96% of the COCs showed 3- or 4 grades of CE (Fig. 2C). While only 15% of the CCs from COCs underwent apoptosis, significant apoptosis was observed in 34% of the CCs from OOXs (Fig. 2D–F). The CCs from COCs contained significantly more hyaluronan than CCs from OOXs did (Fig. 2G). The ratio of Bcl2/Bax mRNAs and the mRNA level of the CE-supporting genes, Has2 and Ptx3, were significantly higher in COCs than in OOXs (Fig. 2H and I). Furthermore, our flow cytometry showed that although the percentage of healthy cells was higher, the percentage of late apoptotic and necrotic CCs was lower significantly in COCs than in OOXs (Fig. 2J–L). Thus, while CCs from OOXs cultured in ATM undergo severe apoptosis but with no potential to support any CE, those from COCs cultured in ESM suffer mild apoptosis but with full potential to support CE. Therefore, comparing DE genes between the two CC models will reveal new genes and miRNAs regulating CE and/or CC apoptosis.

Identification of DE genes and miRNAs in CCs between OOXs and COCs

RNA- and miRNA-sequencing was conducted to identify DE genes and miRNAs in CCs between OOXs and COCs. The RNA-seq identified a total of 2773 DE genes (q < 0.05, FC > 2), and among them, 1333 were upregulated and 1440 were downregulated in COCs relative to OOXs (Fig. 3A). The miRNA-seq identified a total of 90 DE miRNAs (P < 0.05, FC > 1.5), and among them, 48 were upregulated and 42 were downregulated in COCs relative to OOXs (Fig. 3B). The 90 DE miRNAs identified by miRNA-seq targeted 11 242 genes. Four distinct zones are apparent in both the mRNA and miRNA heat maps. Interestingly, the zones with downregulated mRNAs in Fig. 3A happen to correspond to the zones with upregulated miRNAs in Fig. 3B, and vice versa. This suggested that, in general, our DE miRNAs targeted our DE genes, validating the effectiveness of our OOXs and COCs as models to find new ER/AR genes and miRNAs.

To validate the sequencing results, RT-qPCR was performed on selected mRNAs and miRNAs. The results revealed a similar expression pattern to that detected by the RNA-Seq or miRNA-seq in 10 of the 12 genes (Fig. 3C) and in 12 of the 16 miRNAs (Fig. 3D). Thus, our RT-qPCR results were generally in agreement with the sequencing data, and they further confirmed significantly different patterns of gene and miRNA expression between OOXs and COCs.

Pathway enrichment analysis

An association analysis between the 11 242 genes targeted by the DE miRNAs and the 2773 DE genes identified by RNA-seq gave rise to 1597 candidate genes (Fig. 4A). The 1597 candidate genes were significantly enriched in 57 KEGG pathways, from which the top 20 significant pathways were selected in the P-value order (Fig. 4B). According to the literature, most of the top 20 significantly enriched pathways have relations with either CE or apoptosis or both. The top 20 significant KEGG pathways included 215 top enriched genes (Fig. 5A).

PPI analysis

When the 1597 candidate genes were subjected to PPI analysis, 159 proteins with a degree ≥20 were obtained (Fig. 5B). The top 20 genes, which showed a PPI degree of >40, are shown in Table III. According to the literature, all the top 10 genes are involved in apoptosis and CE, and most of them showed a correct tendency of expression in our OOX/COC models.

Identification of the would-be A/E miRNAs

An overlap between the 215 top enriched genes from KEGG and the 159 top degree genes from PPI produced 90 genes, which were targeted by 39 miRNAs (Fig. 5A). Twenty miRNAs were selected from the 39 candidate miRNAs according to their FC. After excluding four miRNAs that showed a same direction of upregulation or downregulation as their targeted genes, we obtained 16 would-be A/E miRNAs (Fig. 5C).

Functional analysis of A/E miRNAs to determine the ER miRNAs

From the 16 would-be A/E miRNAs, nine were selected for functional analysis based on the following criteria: no involvement in CE had been reported; conflicting results were reported on their role in cell apoptosis; new miRNAs discovered in this study; and FC order. To validate the ER miRNAs, COCs were transfected with inhibitors or mimics of the upregulated or downregulated A/E miRNAs, respectively, and the transfected COCs were cultured for 18 h in ESM to observe CE or mRNA expression level of expansion-related genes. When COCs were transfected with inhibitors of miR-212-5p, 149-5p or 31-5p, or with mimics of miR-351-5p or 503-5p, the percentages of oocytes showing 3- or 4 grades of CE decreased significantly compared to those in control COCs transfected with inhibitor control or mimic control (Fig. 6A). Transfection with miR-27a-3p inhibitor also decreased percentages of oocytes with 3- or 4 grades of CE although statistically the difference was not significant (P = 0.152). The results suggested that miR-212-5p, 149-5p, 31-5p and 27a-3p were CE-promoting, whereas miR-351-5p and 503-5p were CE-inhibiting.

As expected, the mRNA level of Has2 decreased after transfection with inhibitors of the CE-promoting miR-212-5p, 149-5p, 31-5p or 27a-3p but, unexpectedly, it increased following transfection with mimics of the CE-inhibiting miR-351-5p or 503-5p (Fig. 6B). As expected, the mRNA level of Ptx3 decreased significantly after transfection with inhibitors of CE-promoting 31-5p or 27a-3p, or with mimics of CE-inhibiting 351-5p or 503-5p but, unexpectedly, it increased after transfection with inhibitors of CE-promoting 212-5p or 149-5p (Fig. 6C). Thus, the results verified miR-212-5p, 149-5p, 31-5p, 27a-3p, 351-5p and 503-5p as ER miRNAs. While miR-212-5p and 149-5p promoted CE mainly by facilitating Has2 expression, miR-31-5p and 27a-3p facilitated expression of both Has2 and Ptx3. miR-351-5p and 503-5p inhibited CE mainly by suppressing Ptx3 expression.

Functional analysis of A/E miRNAs to determine the AR miRNAs

To validate the AR miRNAs, OOXs were transfected with mimics or inhibitors of the upregulated or downregulated A/E miRNAs, respectively, and the transfected OOXs were cultured for 18 h in ATM to observe CC apoptosis or mRNA expression level of apoptosis-related genes. When OOXs were transfected with mimics of miR-212-5p, 149-5p or Nov798, the percentages of apoptotic CCs decreased significantly compared to those in control OOXs transfected with mimic control (Fig. 6D). However, transfection of OOXs with mimics or inhibitors of other A/E miRNAs showed no significant effects on percentages of apoptotic CCs. Our RT-qPCR showed that while the ratio of Bcl2/Bax mRNAs was higher (Fig. 6E), the level of Fas mRNA was significantly lower (Fig. 6F) after transfection of OOXs with mimics of miR-212-5p, 149-5p or Nov798 than with respective mimic controls. The results validated miR-212-5p, 149-5p and Nov798 as AR miRNAs, which inhibited apoptosis through Bcl2/Bax and Fas signaling.

Validation of ER and AR miRNAs using in vivo matured oocytes

Our previous study had demonstrated that in vivo ZEN exposure significantly compromised mouse oocyte competence (Pan et al., 2021), therefore in vivo matured COCs from ZEN-treated mice were used to further verify the ER and AR miRNAs we had identified. The results indicated that although 57% COCs from control mice showed 2- or 3 grades of CE, the percentage decreased significantly to 30% in ZEN-treated mice (Fig. 7A). Levels of Has2 and Ptx3 mRNAs were significantly lower in ZEN-treated than in control CCs (Fig. 7B). While the Bcl2/Bax ratio was lower (Fig. 7G), percentages of apoptotic cells were higher significantly in ZEN-treated than in control CCs as revealed by both Hoechst staining (Fig. 7C) and flow cytometry (Fig. 7D–F). Furthermore, our RT-qPCR demonstrated that while the CE-promoting miR212-5p, 31-5p and 27a-3p, and the apoptosis-inhibiting NovelmiR798 were downregulated, the CE-inhibiting miR351-5p and 503-5p were upregulated significantly in ZEN CCs compared to those in control CCs (Fig. 7H). However, the expression of miR149-5p did not differ significantly between control and ZEN CCs. Together, the results first confirmed that CCs from the in vivo matured COCs from ZEN-treated mice can be used to validate the ER and AR miRNAs and, second, validated six of the seven ER/AR miRNAs that had been identified using our CC models from cultured COCs/OOXs.

Discussion

Our Hoechst labeling for apoptotic CCs, ELISA for hyaluronan content, and RT-qPCR for the mRNA levels of Bcl2/Bax and CE-supporting genes all validated that OOXs cultured in ATM could be used as CC models that undergo severe apoptosis but with no potential to support CE, whereas COCs cultured in ESM could be used as CC models that suffer mild apoptosis but possess full potential to support CE. Both our RNA-seq and RT-qPCR demonstrated that the mRNA expression level of the four important CE genes, including Has2, Ptx3, Tnfaip6 and Ptgs2, was significantly upregulated in COCs relative to OOXs. Our miRNA sequencing indicated that among the seven miRNAs that had been reported as regulating CE (Li et al., 2016a, 2017; Cui et al., 2018; Pan et al., 2018; Vahdat-Lasemi et al., 2019), five showed correct trends of expression. When our RNA-seq data on the well-known pro-apoptotic genes were compared between OOX and COC samples, mRNA levels of Caspase 2 (FC = 0.59), Caspase 3 (FC = 0.61), Caspase 6 (FC = 0.55), Caspase 8 (FC = 0.89) and Fas (FC = 0.38) were found to be significantly (P < 0.05) downregulated in COCs relative to OOXs. All the top 10 genes from our PPI analysis were reported as regulating apoptosis and CE, with most of them showing a correct tendency of expression in the present study (see below for detailed discussion). Furthermore, our analysis using in vivo matured COCs from ZEN-treated mice, which showed a significantly lower CE potential but a higher level of CC apoptosis than controls, validated six of the seven ER and/or AR miRNAs we had identified using our CC models from cultured COCs/OOXs.

According to the literature, most of the top enriched pathways that our KEGG pathway enrichment analysis revealed have relations with either CE or apoptosis or both. For example, the essence of CE is the formation of extracellular matrix (ECM) (Russell and Salustri, 2006), which is involved in regulation of apoptosis (Gilmore et al., 2009). The major ECM components, including ‘proteoglycans’ (Russell and Salustri, 2006), play a critical role in both ‘axon guidance’ (Masu, 2016; Saied-Santiago and Bülow, 2018) and CE (Nagyova, 2018). Several molecules involved in ‘axon guidance’ are regulated by p53 (Arakawa, 2005), which induces apoptosis (Polyak et al., 1997). ‘RAS’ and ‘RAP1’ are required for cAMP signaling to ‘MAPK signaling’ (Li et al., 2016b), which mediated the gonadotrophin-induced CE (Su et al., 2002; Di Giacomo et al., 2016). The ‘ECM receptor interaction’, ‘focal adhesion’ and ‘PI3K/AKT’ pathways were involved in granulosa cell luteinization (He et al., 2016). Inhibition of ‘focal adhesion’ kinase suppressed CE (Kitasaka et al., 2018) and enhanced apoptosis of colon cancer cells (Golubovskaya et al., 2003). miR-21 inhibited CC apoptosis by activation of the ‘PI3K/AKT pathway’ (Han et al., 2017). ‘Pathways in cancer’ is related to apoptosis (Peng et al., 2014; Pandurangan et al., 2018), and it involves the transforming growth factor-β signaling that has a relation with CE (Dragovic et al., 2005). Lactoferrin promoted ‘lipolysis in adipocytes’ by controlling the activity of ‘cAMP/ERK signaling’ (Ikoma-Seki et al., 2015), suggesting that the ‘lipolysis in adipocytes’ pathway might regulate CC apoptosis and CE. The ‘protein digestion and absorption’ pathway is often reported in association with the ‘ECM-receptor interaction’ and ECM organization pathways (Chen et al., 2019; He et al., 2019).

According to the literature, all the top 10 genes from our PPI analysis are involved in apoptosis and CE, and most of them showed a correct tendency of expression in our OOX/COC models. Among the genes upregulated in our COCs, IL-6 inhibits cancer cell apoptosis (Guo et al., 2012), sustains ECM turnover (Yang et al., 2017) and promotes COC expansion (Wang et al., 2014). Vascular endothelial growth factor A showed antiapoptotic effects in cultured granulosa cells (Gao et al., 2020), and promoted CE of porcine oocytes (Liu et al., 2020). Prolyl 4-hydroxylase, β subunit (P4HB) inhibited apoptosis of human HT29 cells (Zhou et al., 2019) and was essential for assembly of functional ECM (Balasubramanian et al., 2018). Hsp70 member 5 (Hspa5) inhibited apoptosis in hepatoma cells (Wang et al., 2012) and was secreted into ECM (Delpino and Castelli, 2002). Furthermore, neural cell adhesion molecule 1 (NCAM1) prevents cell death (Sasca et al., 2019) and is involved in cell-ECM interactions (Berezin et al., 2014). Among genes downregulated in our COCs, hepatocyte growth factor elicited apoptosis and ECM degradation (Mizuno et al., 2005). IL-1β is involved in apoptosis and matrix destruction of intervertebral disc degeneration (Wang et al., 2020). Thrombospondin-1 promoted apoptosis in granulosa cells (Zhu et al., 2019) and it acted as an extracellular mediator of matrix mechano-transduction during vascular remodeling (Yamashiro et al., 2020). Insulin receptor substrate 1 (IRS1) decreased DNA fragmentation in podocytes of kidney (Mima et al., 2020) and activation of the IR/IRS-1/Akt pathway reduced accumulation of ECM in cardiac fibroblasts (Zhao et al., 2020).

Using CC models from cultured COCs/OOXs, we identified and validated six new ER miRNAs, and showed that while miR-212-5p, 149-5p, 31-5p and 27a-3p promoted, miR-351-5p and 503-5p inhibited CE. Furthermore, we validated that miR-212-5p, 149-5p and Nov798 inhibited CC apoptosis involving both Bcl2/Bax and the Fas signaling. Most of the miRNAs reported so far inhibited CE and few were reported promoting CE. For example, among the six miRNAs that had been reported regulating CE, five (Li et al., 2016a, 2017; Cui et al., 2018; Pan et al., 2018; Vahdat-Lasemi et al., 2019) were found inhibiting CE, and only one (Sinha et al., 2017) was found promoting CE. Thus, one significant contribution of this study was the identification and validation of four new CE-promoting miRNAs.

Our study using CC models of cultured OOXs/COCs indicated that miR-212-5p and 149-5p promoted CE while inhibiting apoptosis of CCs. Several miRNAs have been reported to promote CE while preventing apoptosis of CCs or to promote CC apoptosis while inhibiting CE. For instance, miR-21 is antiapoptotic in CCs (Han et al., 2017), and it facilitates CE by down-regulating tissue inhibitor of metalloproteinase 3 (TIMP3) (Pan and Li, 2018). High levels of miR-375 decreased expression of PTX3, HAS2 and PTGS2 while promoting apoptosis by inhibiting p-SMAD2/3 and p-SMAD1/5/8 in bovine CCs (Liu et al., 2018). The expression level of miR-145-5p in human CCs was correlated negatively with oocyte maturation but positively with CC apoptosis (Cui et al., 2018). miR-224 impaired mouse ovulation by decreasing the expression of Ptx3 and Smad4 while increasing apoptosis (Vahdat-Lasemi et al., 2019). Overexpression of miR-378 in porcine CCs impaired cumulus expression by decreasing Has2 and Ptgs2 expression (Pan et al., 2015), and CMV-miR-378 lentivirus transduction in mouse ovaries increased apoptosis rates (Sun et al., 2018).

The current results indicated that different miRNAs regulate CE by modifying expression of different molecules. Thus, while miR-212-5p and 149-5p promoted CE mainly by facilitating Has2 expression, miR-31-5p and 27a-3p promoted CE by facilitating expression of both Has2 and Ptx3. Furthermore, miR-351-5p and 503-5p inhibited CE mainly by suppressing Ptx3 expression. As Yokoo and Sato (2011) reviewed, CE is the result of synthesis and accumulation of hyaluronan around the CC, and hyaluronan is synthesized by hyaluronan synthase (HAS). The Has2 mRNA is strongly expressed in CCs during FSH- or eCG-stimulated CE in mice, pigs and cattle. Hyaluronan needs extracellular hyaluronan-associated proteins, such as inter-a-trypsin inhibitor (ITI), tumor necrosis factor alpha-induced protein 6 (TNFAIP6 or TSG6) and pentraxin 3 (PTX3), to form and maintain the stability of hyaluronan matrix of COCs. Mice deficient in Tnfaip6 (Fulop et al., 2003) or Ptx3 (Salustri et al., 2004) are infertile owing to their inability to organize hyaluronan into a stable matrix, although they synthesized a normal amount of HAS2.

In this study, both RNA-seq and RT-qPCR indicated that the levels of Bcl2 and Bax mRNAs in COCs were either lower than, or not significantly different from, those in OOXs. This is inconsistent with our expectation that the anti-apoptotic BCL2 should be higher, but the pro-apoptotic BAX should be lower in COCs than in OOXs. However, the Bcl2/Bax ratio revealed by our RT-qPCR was always significantly higher in COCs than in OOXs. According to Hengartner (2000), BCL2 and BAX meet at the surface of mitochondria to compete for regulation of cytochrome c release. If the pro-apoptotic molecules prevail, a series of molecules, including cytochrome c, are released from mitochondria to induce apoptosis. Thus, it is the ratio of BCL2/BAX but not their contents that determines the susceptibility of a cell to apoptogenic stimuli.

In summary, to identify new ER and AR miRNAs, we first demonstrated that COCs cultured in ESM supported full CE while suffering mild apoptosis, whereas OOXs cultured in ATM underwent severe apoptosis and supported no CE. Then, we conducted RNA- and miRNA-sequencing and bioinformatics using CCs from the cultured COCs/OOXs, and identified candidate ER and/or AR miRNAs. We then performed functional analysis on the candidate miRNAs by transfecting COCs and OOXs with mimic or inhibitor, and validated six new ER miRNAs, and three new AR miRNAs. Our further analysis indicated that different ER miRNAs regulate CE by modifying CC expression of Has2, Ptx3 or both, and that the three new AR miRNAs inhibited apoptosis through the Bcl2/Bax and Fas signaling. Furthermore, our experiment using in vivo matured COCs from ZEN-treated mice further verified the above ER/AR miRNAs, except for miR-149-5p. The new AR and ER miRNAs we identified can be used not only as biomarkers for the selection of competent oocytes and embryos, but also for studying the mechanisms by which the oocyte-derived factors regulate CE and CC apoptosis.

Data availability

The datasets presented in this study can be found in the online repository (https://www.ncbi.nlm.nih.gov/) with accession number of GSE189959.

Authors’ roles

X.H., M.Z., H.-J.Y., G.-L.W., X.-Y.Z., Z.-B.L. and S.G. conducted the experiments; X.H., M.Z. and S.G. analyzed the data; and J.-H.T. designed the experiments and wrote the manuscript. All authors reviewed the manuscript.

Funding

This study was supported by grants from the National Natural Science Foundation of China (31902160, 31772599 and 32072738), the National Key R&D Program of China (2017YFC1001601 and 2017YFC1001602), the Natural Science Foundation of Shandong Province (ZR2017BC025) and the Funds of Shandong Double Tops Program (SYL2017YSTD12).

Conflict of interest

The authors declare that they have no conflicts of interest.

References

Author notes