Circular RNA expression profiling of human granulosa cells during maternal aging reveals novel transcripts associated with assisted reproductive technology outcomes

Circular RNAs (circRNAs) are a unique class of endogenous RNAs which could be used as potential diagnostic and prognostic biomarkers of many diseases. Our study aimed to investigate circRNA profiles in human granulosa cells (GCs) during maternal aging and to uncover age-related circRNA variations that potentially reflect decreased oocyte competence. CircRNAs in GCs from in vitro fertilization (IVF) patients with young age (YA, ≤ 30 years) and advanced age (AA, ≥ 38 years) were profiled by microarray, and validated in 20 paired samples. The correlation between circRNAs expression and clinical characteristics was analyzed in additional 80 samples. Chip-based analysis revealed 46 up-regulated and 11 down-regulated circRNAs in AA samples (fold change > 2.0). Specifically, circRNA_103829, circRNA_103827 and circRNA_104816 were validated to be up-regulated, while circRNA_101889 was down-regulated in AA samples. After adjustment for gonadotropin treatment, only circRNA_103827 and circRNA_104816 levels were positively associated with maternal age (partial r = 0.332, P = 0.045; partial r = 0.473, P = 0.003; respectively). Moreover, circRNA_103827 and circRNA_104816 expressions in GCs were negatively correlated with the number of top quality embryos (r = -0.235, P = 0.036; r = -0.221, P = 0.049; respectively). Receiver operating characteristic (ROC) curve analysis indicated that the performance of circRNA_103827 for live birth prediction reached 0.698 [0.570–0.825], with 77.2% sensitivity and 60.9% specificity (P = 0.006), and that of circRNA_104816 was 0.645 [0.507–0.783] (P = 0.043). Bioinformatics analysis revealed that both circRNAs were potentially involved in glucose metabolism, mitotic cell cycle, and ovarian steroidogenesis. Therefore, age-related up-regulation of circRNA_103827 and circRNA_104816 might be potential indicators of compromised follicular micro-environment which could be used to predict IVF prognosis, and improve female infertility management.

Circular RNAs (circRNAs) are a unique class of endogenous RNAs which could be used as potential diagnostic and prognostic biomarkers of many diseases. Our study aimed to investigate circRNA profiles in human granulosa cells (GCs) during maternal aging and to uncover age-related circRNA variations that potentially reflect decreased oocyte competence. CircRNAs in GCs from in vitro fertilization (IVF) patients with young age (YA, 30 years) and advanced age (AA, ! 38 years) were profiled by microarray, and validated in 20 paired samples. The correlation between circRNAs expression and clinical characteristics was analyzed in additional 80 samples. Chip-based analysis revealed 46 up-regulated and 11 down-regulated circRNAs in AA samples (fold change > 2.0). Specifically, circRNA_103829, circRNA_103827 and circRNA_104816 were validated to be upregulated, while circRNA_101889 was down-regulated in AA samples. After adjustment for gonadotropin treatment, only circRNA_103827 and circRNA_104816 levels were positively associated with maternal age (partial r = 0.332, P = 0.045; partial r = 0.473, P = 0.003; respectively). Moreover, circRNA_103827 and circRNA_104816 expressions in GCs were negatively correlated with the number of top quality embryos (r = -0.235, P = 0.036; r = -0.221, P = 0.049; respectively). Receiver operating characteristic (ROC) curve analysis indicated that the performance of circRNA_103827 for live birth prediction reached 0.698 [0.570-0.825], with 77.2% sensitivity and 60.9% specificity (P = 0.006), and that of circRNA_104816 was 0.645 [0.507-0.783] (P = 0.043). Bioinformatics analysis revealed that both circRNAs were potentially involved in glucose metabolism, mitotic cell cycle, and ovarian steroidogenesis. Therefore, age-related up-regulation of circRNA_103827 and circRNA_104816 might be potential indicators of compromised follicular micro-environment which could be used to predict IVF prognosis, and improve female infertility management. a1111111111 a1111111111 a1111111111 a1111111111 a1111111111

Introduction
Advanced maternal age was closely associated with decreased pregnancy rates and increased spontaneous miscarriages [1]. Even seeking assisted reproductive technology (ART), only 12.2% of women aged 41-42, and 4.2% of women older than 42 years achieved live births with their own oocytes, which were remarkably lower than that of women younger than 35 years (40.1%) [2]. Due to the postponement of childbearing in modern society [3] and low success rate of ART in older women [2], age-related decline in female fertility has become a health issue urgent to be solved. Elucidating underlying molecular signatures that link maternal age and fecundity decline is particularly important for developing predictive markers and therapeutic targets for improving reproductive outcomes in aging women.
Ovarian aging or diminished ovarian function, which is generally attributed to progressive depletion of follicular reserve, decreased oocyte competence and poor granulosa cells (GCs) [4][5][6][7][8], has been recognized as the main driving force underlying age-related female subfertility and adverse ART outcomes [9]. Since GCs provide an essential follicular micro-environment that influences oocyte competence and early embryo development [10][11][12], many studies have tried to use the transcriptome of GCs as a non-invasive approach to gain insight into microenvironment surrounding oocytes and predict oocyte developmental potential [13][14][15][16]. However, there is a lack of overlap among different studies possibly owing to population heterogeneity and the instability of linear RNAs [15]. Therefore, exploring other contributors with better stability appears to be necessary. Circular RNAs (circRNAs) are a unique class of non-coding RNAs with a stable structure formed by special loop splicing. To date, thousands of highly expressed circRNAs have been successfully identified across different eukaryotic organisms [17]. These endogenous RNAs often exhibit tissue-and developmental-phase-specific expression patterns [18,19] and can be sorted into exosomes [20]. Compared with traditional linear RNAs, circRNAs do not have 5' to 3' polarity or a polyadenylated tail, but are structured as a covalently closed ring instead, thus resisting RNA exonucleases mediated degradation and maintaining a high stability [21,22]. In addition, circRNAs contribute to post-transcriptional regulation of gene expression by binding specifically to microRNAs (miRNAs) and sequestering miRNAs to terminate suppression of their targets [23]. Increasing evidences indicate that circRNAs play important roles in aging [24,25] and many age-related diseases, such as circ-Foxo3 in cardiac senescence [26], ciR-7 in Alzheimer's disease [27], and cir-ITCH in colorectal cancer [28]. Moreover, some circRNAs have been suggested as novel diagnostic biomarkers and potential therapeutic targets, such as hsa_circ_0005075 for prediction of hepatocellular carcinoma [29], circRNA_101222 for prediction of pre-eclampsia [30], and circRNA-CER for the treatment of osteoarthritis [31]. However, the role of circRNAs in reproduction and their contribution to ovarian aging remains to be discovered.
The aim of this study was to characterize circRNA expression profiles of GCs according to maternal age, and identify age-related changes that could potentially reflect oocyte competence, thus providing new insight into the molecular signatures of ovarian aging.

Patients recruitment and study design
This prospective study included 126 women enrolled in ART treatment at Reproductive Medicine Center of Tongji Hospital (Wuhan, Hubei, China) between February 2015 and August 2015. Only one cycle from each patient was included in this study. Exclusion criteria were patients with polycystic ovary syndrome (PCOS), endometriosis, premature ovarian failure, preimplantation genetic screening/preimplantation genetic diagnosis cycles, and young women ( 30 years old) with diminished ovarian reserve (DOR, indicated by serum anti-Müllerian hormone (AMH) < 1.77ng/ml [32], or antral follicle count (AFC) < 5 on day 3 of menstrual cycle).
This study was designed as three stages. In stage 1, GCs samples from women of advanced age (AA, ! 38 years, n = 3) and young age (YA, 30 years, n = 3) were analyzed by a human circRNA microarray (8×15K, Arraystar Inc.). In stage 2, significantly differentially expressed circRNAs discovered in the microarray were validated in additional 20 AA and 20 YA samples by real-time quantitative reverse transcription-polymerase chain reaction (qRT-PCR). In stage 3, GC samples from another independent cohort (80 women range from 22 to 48 years old) were analyzed by qRT-PCR to further verify the association of diminished ovarian function with candidate circRNAs which were validated in stage 2. Clinical and demographic characteristics of patients were collected from electronic patient charts and summarized in Table 1, S1 and S2 Tables.
This study was approved by Institutional Review Board of Tongji Hospital (Wuhan, Hubei, China) and carried out in accordance with the principles expressed in Declaration of Helsinki. All patients who took part in the prospective observational study (NCT02294500) were informed about GCs sample collection/analysis, and written informed consents were given.

ART procedures
Routine ovarian stimulation, oocyte retrieval, in vitro fertilization (IVF) / intracytoplasmic sperm injection (ICSI), embryo culture and transfer were performed at Reproductive Medicine Center of Tongji Hospital. Oocytes were considered as normally fertilized if two pronuclei and two polar bodies were observed 18-20 hours after microinjection or insemination. Embryo quality was assessed by microscopic morphological observation. Embryos that contained 7-8 regular blastomeres and less than 20% fragments on day 3 were considered as top quality. Top quality embryos were selected for transfer or freezing, whereas others were cultured up to day 5. Blastocysts on day 5 were classified according to the scoring system developed by Gardner [33]. Clinical pregnancy was confirmed by the observation of at least one gestational sac and embryonic heart activity on ultrasound examination. These patients were followed up until delivery. The detailed IVF outcomes of the participants were reported in Table 2 and S1 Table.

GCs isolation
On oocyte retrieval day, all follicular fluid (FF) samples from the same patient were pooled after cumulus-oocyte complexes (COCs) were isolated for conventional IVF or ICSI procedures. GCs were individually purified from FF samples with similar methods as previously described [10,34]. Briefly, FF sample of each person was centrifuged and resuspended in phosphate-buffered saline (PBS). Then, it was slowly placed on 50% Percoll gradient and centrifuged at 400 × g for 30 min at 4˚C. The cells in the middle layer were carefully collected, washed and resuspended in PBS. Next, the modified cell strainer methodology was applied to purify GCs. In order to confirm the purity of GCs, cells were stained with antibodies against CD45 (leukocyte specific antigen, not expressed by GCs) [34]. Flow cytometry revealed that the levels of CD45+ cells contamination in the samples were less than 5%. The harvested GCs were stored in TRIzol reagent (Invitrogen, Karlsruhe, Germany) at -80˚C until RNA extraction.

Total RNA extraction
All GCs samples were processed for total RNA extraction using TRIzol reagent following manufacturer's instructions. The quantity (ng/mL) and purity (260/280 and 260/230 ratios) of RNAs were measured by a NanoDrop 2000 spectrophotometer and denaturing agarose gel electrophoresis. For RNase R digestion, 2μg RNA were incubated with 2IU RNase R (RNR-07250, Epicentre) at 37˚C for 10 min followed by heat inactivation at 95˚C for 3 min.

CircRNA microarray hybridization and analysis
Sample preparation and microarray hybridization were performed according to the Arraystar's standard protocols. Briefly, total RNAs were digested with RNase R to remove linear RNAs  and enrich circRNAs. The enriched circRNAs were then amplified and transcribed into fluorescent cRNAs utilizing random primers according to Arraystar Super RNA Labeling protocol (Arraystar, Inc.). The labeled cRNAs were purified by RNeasy Mini Kit (Qiagen). As with quantification method for other RNAs [35], the concentration and specific activity of labeled cRNAs (pmol Cy3/μg cRNA) were measured by NanoDrop ND-1000. 1μg of the labeled cRNAs were hybridized onto the Arraystar Human circRNA Arrays (8x15K, Arraystar), and incubated for 17 hours at 65˚C in an Agilent Hybridization Oven. The hybridized arrays were washed, fixed and scanned using the Agilent Scanner G2505C. Scanned images were imported into Agilent Feature Extraction software for raw data extraction. A series of data processing including quantile normalization were performed using R software packages. The statistical significance of the difference was estimated by t-test. CircRNAs having fold changes ! 2 and P-values < 0.05 were recognized as significantly differentially expressed. The Arraystar's home-made computer program based on TargetScan and miRanda was applied to predict their miRNA targets and the circRNA/ miRNA interaction. The miRNA support vector regression (mirSVR) algorithm was utilized to score the efficiency of predicted miRNA targets.

Characterization and validation of candidate circRNAs by RT-PCR
Based on fold change > 3, we selected the top 10 among 17 up-regulated circRNAs, and all 4 down-regulated circRNAs for subsequent validation. qRT-PCR was performed using Fast SYBR Green Master Mix (4385612, Thermo Fisher Scientific) according to the manufacturer's instructions. Three paired divergent primers encompassing circRNA-specific back-splice sites were designed for each candidate circRNA. Sequences of the primers achieving a single peak in melting curve were provided in S3 Table. Relative expression levels of candidate circRNAs were normalized to glyceraldehyde-3-phosphate dehydrogenase (GAPDH) and were calculated according to 2 -ΔΔCt method [29,30,36]. In addition, the circular characteristics of differentially expressed circRNAs were assessed by two methods: testing their enrichment in poly(A)-fractions and detecting their resistance to RNase R digestion. RNA with and without poly(A) tail were fractionated using oligo-d (T) 25 -Magnetic Beads (E7490L, New England Biolabs Inc.). Then, poly(A)+ RNA, poly(A)-RNA, digested RNA and control total RNA were reversely transcribed into cDNA using the PrimeScript RT reagent kit (RR037A, TaKaRa) with random primers. Next, candidate cir-cRNAs were amplified from cDNA by KOD Xtreme Hot Start Polymerase Kit (71975-3, Novagen). Eventually, PCR products were detected by agarose gel electrophoresis and Sanger sequencing.
Mann-Whitney U test and Kruskal-Wallis test were used to compare circRNA expression levels between two groups and multiple groups in stage 3, respectively. Spearman correlation analysis was used to investigate the correlation of circRNA expressions and clinical characteristics. The ability of circRNA levels in GCs to predict pregnancy outcomes was determined by constructing the Receiving Operator Curve (ROC) curve and calculating the area under the curve (AUC) with 95% confidence intervals (CI). The sensitivity and specificity for the optimal cut-off were calculated. Statistical tests were performed using SPSS 16.0 software. Results were considered significant when P < 0.05.

Results
Stage 1: CircRNA expression profiles in GCs according to female age GCs RNA from three YA and three AA samples were analyzed by circRNA microarray. Unsupervised hierarchical clustering of circRNA expression patterns discriminated YA from AA women obviously (Fig 1A), revealing different circRNA expression profiles in GCs during  (Fig 1B). Overall, 57 circRNAs were discovered to be significantly differentially expressed (fold change > 2.0, P-value < 0.05). Compared to YA samples, 46 circRNAs were up-regulated and 11 circRNAs were down-regulated in AA samples (Table 3 and S4 Table). The complete microarray data in this study was available at Gene Expression Omnibus (GEO) database under the accession number GSE97193.

Stage 2: Identification and characterization of differentially expressed circRNAs in GCs during maternal aging
The top 10 up-regulated and top 4 down-regulated circRNAs ranked by fold change were selected for validation by qRT-PCR. As this was the first time to use these divergent primers for circRNAs quantification, it was necessary to test their specificity. Our melting curve analysis revealed that eight candidate transcripts showed a single peak (S1 Fig), including six upregulated circRNAs (circRNA_103829, circRNA_103828, circRNA_103827, circRNA_ 100833, circRNA_104816 and circRNA_104852) and two down-regulated circRNAs (circRNA_101889 and circRNA_103611) in AA samples. Single-peak curves were not obtained from other six circRNAs, although three pairs of divergent primers were designed for each of them.
Then we compared the expression of these eight candidate circRNAs in GCs from additional 20 YA and 20 AA women by qRT-PCR. Among them, circRNA_103829, circRNA_103827 and circRNA_104816 were significantly up-regulated, while circRNA_101889 were remarkably down-regulated in AA women by Mann-Whitney U Test (P = 0.004, P = 0.009, P < 0.001, P = 0.010, respectively) (Fig 2). Likewise, Spearman's correlation analysis showed that these four circRNAs levels were significantly associated with maternal age (r = 0.312, P = 0.050; r = 0.446, P = 0.004; r = 0.597, P < 0.001; r = -0.369, P = 0.019; respectively). In addition to maternal age, GCs transcriptome profiles may be affected by gonadotropin treatment [37,38]. Interestingly, we found that the expression of the four circRNAs was comparable in GCs from women who received agonist or antagonist protocols. Conversely, circRNA_101889 expression in GCs varied significantly according to the type of gonadotropin treatment (recombinant follicle-stimulating hormone (r-FSH) and/or highly purified human menopausal gonadotropin (HP-hMG) (P = 0.024) (S2 Fig). Besides, circRNA_103829 was significantly up-regulated in GCs from women who received high doses of gonadotropins (! 2500 IU/l) compared to those treated with doses < 2500 IU/l (P = 0.022) (S2 Fig). After adjustment for gonadotropin treatment, only circRNA_103827 and circRNA_104816 levels were significantly and positively associated with maternal age in partial correlation analysis (partial r = 0.332, P = 0.045; partial r = 0.473, P = 0.003; respectively).
The unique circular characteristics of circRNA_103827 and circRNA_104816 were assessed by poly(A)-enrichment and RNase R digestion. As shown in S3 Fig, both circRNAs were abundant in poly(A)-fractions, barely detectable in poly(A)+ fractions, and resistant to 3'-5' exoribonuclease RNase R digestion. In addition, their back-splice sequences were successfully detected in RNase R treated RNA (S3 Fig). These results strongly support that they are indeed circRNAs. Considering the relationship between circRNAs and maternal age, we continued to investigate whether the variations of cirRNA_103827 and circRNA_104816 in GCs could reflect reduced ovarian reserve and oocyte quality in another cohort (80 women, aged from 22-48 years). Currently, ovarian reserve status is mainly evaluated by AMH and AFC, while oocyte quality could be preliminary assessed based on embryo development outcomes. As shown in Fig 3,  circRNA_103827 and circRNA_104816 were significantly up-regulated in GCs from women with serum AMH concentration < 2ng/ml than in those with higher AMH levels (! 2ng/ml) (P = 0.028; P = 0.013; respectively). Spearman's correlation analysis indicated a significant negative correlation of circRNA_103827 / circRNA_104816 levels and serum AMH concentration (r = -0.283, P = 0.011; r = -0.311, P = 0.005; respectively) (S5 Table). Moreover, circRNA_ 103827 and circRNA_104816 tended to be up-regulated in GCs from women with low AFC (< 10) than those from women with AFC ! 10 (P = 0.006; P = 0.012; respectively) (Fig 3 and Table 1). Indeed, both circRNA_103827 and circRNA_104816 expression levels were found to be significantly and negatively correlated with AFC (r = -0.305, P = 0.006; r = -0.314, P = 0.005; respectively) (S5 Table).
In addition to ovarian reserve parameters, GCs from women with a low number of retrieved oocytes ( 5) contained significant higher circRNA_103827 levels than those with higher number of oocytes (> 5) (P = 0.018) (Fig 3 and Table 2). Spearman's correlation analysis also showed that circRNA_103827 levels in GCs were negatively correlated with the number of retrieved oocytes (r = -0.345, P = 0.002) (S6 Table), indicating the relationship between up-regulation of circRNA_103827 and poor ovarian response. On the other hand, circRNA_103827 and circRNA_104816 expressions did not vary significantly between different ovarian stimulation protocols (S7 Table).
For embryo outcomes at day 3 post-fertilization, we found that oocytes cohorts which gave rise to a small number of top quality embryos (< 2) were related to GCs with higher cir-cRNA_103827 and circRNA_104816 expression levels than those achieving at least two top quality embryos (P = 0.009, P = 0.037, respectively) (Fig 3 and Table 2). The expression of cir-cRNA_103827 and circRNA_104816 was also shown to be significantly and negatively correlated with the number of top quality embryos (r = -0.235, P = 0.036; r = -0.221, P = 0.049; respectively) (S6 Table). In addition, the ratio between number of top quality embryos and total number of mature oocytes was calculated to estimate global developmental competence of each oocyte cohort. CircRNA_103827 tended to be up-regulated in GCs of embryo cohorts with low proportion of top quality embryos ( 25%; decreased developmental competence) than those with top quality embryos > 25% (P = 0.047) (Fig 3 and Table 2). These data suggest that up-regulated circRNA_103827 and circRNA_104816 may be associated with adverse embryo outcomes.

CircRNA_103827 and circRNA_104816 predictive values for pregnancy outcomes
Given that age-related decline in ovarian function is the most important factor for determining pregnancy success in women of advanced age [9], we further investigate whether the expression of cirRNA_103827 and circRNA_104816 in GCs could predict pregnancy outcomes in assisted reproduction. ROC curve analysis revealed the performance of circRNA_103827 for clinical pregnancy prediction was 0.725 [0.608-0.842] (P = 0.001), with 79.6% sensitivity and 61.3% specificity, and that of circRNA_104816 was 0.672 [0.550-0.795] (P = 0.01). Moreover, the predictive power of GCs circRNA_103827 expression for clinical pregnancy was higher than that of top quality embryo percentage (AUC = 0.674 [0.551-0.797]; P = 0.01) in our population (Fig 4A-4C).
Apart from clinical pregnancy, circRNA_103827 have predictive values for live births as well. ROC curve analysis indicated that the performance of circRNA_103827 for live birth Circular RNA expression profiles of human granulosa cells associated with ovarian aging prediction reached 0.698 [0.570-0.825] (P = 0.006), with 77.2% sensitivity and 60.9% specificity, and that of circRNA_104816 was 0.645 [0.507-0.783] (P = 0.043). In our population, the power of circRNA_103827 expression for live birth prediction was higher than that of top quality embryo percentage (AUC = 0.661 [0.533-0.789], P = 0.027) (Fig 4D-4F).

Prediction and annotation of circRNA_103827 and circRNA_104816 targeted miRNA-gene network
To explore biological function of circRNA_103827 and circRNA_104816, we used bioinformatics to predict their potential targets and interactive networks. Given that circRNAs were reported to function as miRNA sponges [23], we assumed that circRNA_103827 and circRNA_104816 may regulate gene expression by a similar mechanism. The molecular interaction of circRNA_103827 and circRNA_104816 with their five highest ranking candidate miRNA targets (Top 5) based on mirSVR scores were depicted in S4 Fig. Since the complete competing endogenous RNA (ceRNA) relationship could not be presented clearly in a map, we constructed an interactive network with two circRNAs and their respective relevant top 5 miRNAs and top 6 genes. Interestingly, this network revealed that circRNA_103827 and Circular RNA expression profiles of human granulosa cells associated with ovarian aging circRNA_104816 may regulate the same target genes, including ANKRD20A9P, KCNQ1OT1 and XIST (Fig 5A and S5 Fig).
To gain further insights into functions of circRNA_103827 and circRNA_104816, GO and KEGG pathway analysis were utilized based on their predicted targets. CircRNA_103827 and circRNA_104816 showed a strong relationship with biological processes of glucose import, glucose metabolic process, glucose transport, G-protein coupled receptor signaling pathway, and G2/M transition of mitotic cell cycle ( Fig 5B). Pathway analysis revealed that circRNA_103827 and circRNA_104816 were predicted to be involved in ovarian steroidogenesis ( Fig 5C).

Discussion
Ovarian aging, termed as decline of ovarian function with increasing age, is considered as the trigger of female aging and subfertility, and it is associated with many age-related diseases such as osteoporosis and cardiac vascular diseases [5,39]. Prolonged human lifespan and postponed Circular RNA expression profiles of human granulosa cells associated with ovarian aging childbearing prompt researchers to investigate mechanisms underlying ovarian aging. Extensive studies have focused on characterizing mRNA and miRNAs profiles in human GCs according to maternal age [10,15], and expanded our knowledge of molecular processes occurring in follicular cells. However, the application of new sequencing technologies provides opportunities to study other novel molecules with promising prospects.
In this study, we introduced circRNAs to shed new light on ovarian aging mechanisms. There are three important properties for circRNAs: first, they are highly conserved and tissuespecific sequences; second, they can tolerate exonuclease degradation, and exhibit high degree of stability in mammalian cells; third, they have unique ceRNA characteristics [17,23,40]. Compared with other linear RNAs, these properties provide circRNAs with the potential to be used as ideal biomarkers and therapeutic targets. Therefore, it is worth exploring circRNAs in human reproduction.
Our study is the first to profile circRNA expression patterns of human GCs according to female age, and uncovered a general trend of circRNA up-regulation in GCs during aging. This age-accumulation trend of circRNA was previously reported in mouse brain tissues [24]. Specifically, circRNA_103827 and circRNA_104816 in GCs were found to be up-regulated in older women. Moreover, the expression of the two circRNAs were negatively correlated with serum AMH levels and AFC, which have been demonstrated to be excellent biomarkers of ovarian reserve [41,42]. This relationship indicates the important role of circRNAs in ovarian aging. Considering the high dispersion of circRNA_103827 and circRNA_104816 expressions in populations, they may not be suitable for use as predictors of ovarian function alone, but serve as complementary indicators.
Oocyte quality, which can be influenced and reflected by follicular micro-environment, has been well recognized as a key limiting factor in female fertility [9]. Many investigations have tried to use proteins, mRNA and miRNAs as surrogate biomarkers for predicting oocyte developmental competence and pregnancy outcomes after ART [10,14,[43][44][45]. Despite some progress, the universally applicable and non-invasive test for the determinant of pregnancy outcomes has not yet been established. For example, the performance of serum AMH for pregnancy prediction was 0.634 (95% CI, 0.618-0.650), with 44.0% sensitivity and 66.5% specificity [46]. In addition, the predictive power of FF miR-29a levels for clinical pregnancy was 0.68 (95% CI, 0.55-0.79) with a sensitivity of 83.3%, but a low specificity (53.5%) [47]. In our study, circRNA_103827 and circRNA_104816 expressions in GCs were found to be negatively associated with embryonic development. The predictive performance of circRNA_103827 expression in GCs for clinical pregnancy was 0.725, with 79.6% sensitivity and 61.3% specificity, and that of circRNA_104816 was 0.672, with 74.1% sensitivity and 57.2% specificity. Thus, these indicators have their own advantages in terms of sensitivity and specificity. Combination of proteins, mRNAs or microRNAs in previous studies and circRNAs in this study may be beneficiary for clinical evaluation of ovarian function and pregnancy outcomes. Our study adds new information and promising biomarkers to the field. Given the effect of population heterogeneity, the results need to be validated in larger number size and other populations.
For potential biological functions, circRNA_103827/circRNA_104816-miRNA-mRNA networks suggest that both circRNAs may regulate common target genes, ANKRD20A9P, KCNQ1OT1 and XIST. ANKRD20A9P has been identified as a novel loci associated with longevity in Han Chinese GWAS [48]. XIST (X-inactive specific transcript) plays critical roles in the initiation of imprinted X-chromosome inactivation. In females, lack of XIST leads to genome-wide transcriptional misregulation in early blastocysts and post-implantation lethality [49]. KCNQ1OT1 is a paternally expressed allele involved in the transcriptional silencing by regulating histone methylation, and its misregulation could lead to disastrous physical and genetic abnormalities [50]. Epigenetic dysregulation has been identified as an important mechanism underlying ovarian aging and abnormal embryonic development [51,52]. Accordingly, the critical roles of XIST and KCNQ1OT1 in epigenetics support our hypothesis that agerelated up-regulation of circRNA_103827 and circRNA_104816 may be involved in ovarian aging.
Investigations of GCs could contribute to understanding compromised follicular microenvironment in the aging ovary, and provide insight into the mechanisms underlying reduced oocyte quality. Based on previous studies, GCs from women aged > 38 years showed decreased proliferation, reduced steroidogenic potential, increased apoptotic changes, abnormal glucose metabolism and more damaged mitochondria compared to younger women [7,8]. Proteomic analysis revealed that genes involved in steroid hormone synthesis, oxidative phosphorylation and post-transcriptional splicing processes were down-regulated with advancing female age [16]. Our in silico prediction suggests that circRNA_103827 and circRNA_104816 were potentially involved in glucose metabolism, mitotic cell cycle, and ovarian steroidogenesis. Therefore, age-related up-regulation of circRNA_103827 and circRNA_104816 may be responsible for follicular micro-environment deterioration with consequences for poor oocyte and embryo quality. Next, we will explore corresponding circRNA functions using overexpression and knockdown experiments.
Despite these novel and promising findings in this study, there were still some weaknesses we recognized. Firstly, one limitation is the small sample size in microarray analysis which may decrease statistical power. Considering this limitation, we expanded the sample size to 20 patients in each arm for subsequent validation to obtain credible results. Another limitation is heterogeneity in the study population. However, we controlled other confounding factors such as PCOS and endometriosis to minimize the differences, and adjusted ovarian stimulation protocols by statistical analysis. Thirdly, the physiological significance of the circRNAs in reproduction is only predicted. The novel identified circRNAs associated with ovarian aging and the prediction of their potential functions laid a solid foundation for further investigation.
In conclusion, our study was the first to characterize global circRNA expression patterns in human GCs based on female age, and demonstrated the significant difference between young and older women. Elevated expression of circRNA_103827 and circRNA_104816 were closely related to declining ovarian reserve and adverse reproductive outcomes. Intriguingly, circRNA_103827 has predictive performance for pregnancy outcomes after ART cycles. This study sheds new light on altered follicular micro-environment in aging ovary and adds novel transcripts associated with reproductive outcomes. Deciphering molecular mechanisms and biological functions of circRNAs would improve personalized ART strategies and identify new biomarkers and therapeutic targets in infertility female with advanced age.  Table. Clinical characteristics and assisted reproductive technology outcomes of the patients used for microarray analysis in stage one. YA, young age; AA, advanced age; BMI, body mass index; FSH, follicle-stimulating hormone; LH, luteinizing hormone; E2,17βestradiol; T, testosterone; PRL, prolactin; AMH, anti-Müllerian hormone; AFC, antral follicle count; IVF, in vitro fertilization; ICSI, intracytoplasmic sperm injection; Gn, gonadotropin; OPU, oocyte pick-up; COCs, cumulus oocyte complexes; MII, oocyte blocked in meiotic metaphase II; 2PN, 2-pronuclear fertilization; Good quality embryos, embryos with 7-8 regular blastomeres and less than 10% fragments on day 3. (DOCX) S2 Table. Clinical characteristics of women used for candidate circRNAs validation in stage two. YA, young age; AA, advanced age; E2,17β-estradiol; PRL, prolactin; T, testosterone; AMH, anti-Müllerian hormone; AFC, antral follicle count; IVF, in vitro fertilization; ICSI, Intracytoplasmic sperm injection. Gn, gonadotropin; P4, progesterone; NS, no significance (P ! 0.05). Numerical values are in the form of Mean ± SD. a Data were analyzed by two-tailed t test; b Chi-square test. Ã P < 0.05; ÃÃ P < 0.01. (DOCX) S3  Table. CircRNA_103827 and circRNA_104816 expression levels in granulosa cells relative to ovarian stimulation protocols. r-FSH, recombinant follicle-stimulating hormone; HP-hMG, highly purified human menopausal gonadotropin; E2,17β-estradiol. Ã P < 0.05; ÃÃ P < 0.01. (DOCX) S1 File. STROBE statement-Checklist of items that should be included in reports of casecontrol studies. (DOC)