Skip to main content
Download PDF

Maternal age affects equine day 8 embryo gene expression both in trophoblast and inner cell mass

BMC Genomics volume 23, Article number: 443 (2022) Cite this article

Abstract

Background

Breeding a mare until she is not fertile or even until her death is common in equine industry but the fertility decreases as the mare age increases. Embryo loss due to reduced embryo quality is partly accountable for this observation. Here, the effect of mare’s age on blastocysts’ gene expression was explored. Day 8 post-ovulation embryos were collected from multiparous young (YM, 6-year-old, N = 5) and older (OM, > 10-year-old, N = 6) non-nursing Saddlebred mares, inseminated with the semen of one stallion. Pure or inner cell mass (ICM) enriched trophoblast, obtained by embryo bisection, were RNA sequenced. Deconvolution algorithm was used to discriminate gene expression in the ICM from that in the trophoblast. Differential expression was analyzed with embryo sex and diameter as cofactors. Functional annotation and classification of differentially expressed genes and gene set enrichment analysis were also performed.

Results

Maternal aging did not affect embryo recovery rate, embryo diameter nor total RNA quantity. In both compartments, the expression of genes involved in mitochondria and protein metabolism were disturbed by maternal age, although more genes were affected in the ICM. Mitosis, signaling and adhesion pathways and embryo development were decreased in the ICM of embryos from old mares. In trophoblast, ion movement pathways were affected.

Conclusions

This is the first study showing that maternal age affects gene expression in the equine blastocyst, demonstrating significant effects as early as 10 years of age. These perturbations may affect further embryo development and contribute to decreased fertility due to aging.

Peer Review reports

Background

Unlike other farm animals, the main purpose for the use of horses is athletic performance and leisure. In many cases, mares produce foals regularly to remain profitable after their sport career so that their reproductive career may last until advanced age, even though mares’ fertility declines with age (for review [1]). Indeed, in Thoroughbreds, pregnancy rates per estrous at Day 42 post-ovulation were shown to be reduced from 63.0% in young mares to 58.5% in mares older than 9 years of age, and to reach 43.8% in mares older than 18 years old [2]. It is difficult to determine an age limit defining when a mare is considered too old for breeding as aging is a progressive process during the life. One study on 535,746 fertility records in mixed breeds in France, however, suggested that if a 10% of fertility decline is chosen as cutoff value for retiring breeding mares, then mares older than 8 should not be bred [3].

Embryo loss seems to be the principal cause for the fertility decline in old mares. Using pregnancy data obtained from 898 cycles in Thoroughbred mares over period of 11 years, the cutoff for embryo mortality to reach values above average was calculated to be 10 years of age [4]. In practice, embryo loss frequency is generally observed between the first (14–22 days) and the second (42–45 days) pregnancy check in old than in young mares [2, 4, 5]. Embryo loss may, however, be unobserved in earlier stages. Indeed, although experimental oviductal flushes on Day 2 post-ovulation showed that the fertilization rate was similar between young and fertile mares vs old and subfertile mares (respectively, 10/14 vs 11/14 embryos recovered per group on Day 2), when the same mares were inseminated during another estrus period, the pregnancy rate was significantly higher at 14 days post ovulation in young and fertile mares compared to older ones (respectively, 12/15 vs 3/15 pregnant mares on Day 14 by ultrasonography) [6].

Although uterine degeneration, including large inflammatory cell infiltration, fibrosis and reduced gland density [7,8,9] is observed in old mares, oocyte quality seems to play a major role in the observed age-related fertility decrease: the transfer of Day 4 embryos collected from old and subfertile mares resulted in lower pregnancy rates at 14 days post-ovulation than when embryos collected from young and fertile mares were transferred. This remained true when only the subset of good quality embryos was considered [10].

Controversial effects of maternal age on oocyte morphology have been observed in several studies [11,12,13,14,15] but there is a clear effect of maternal age on oocyte quality and function. Old mares’ oocytes are less likely to reach metaphase II after in vitro maturation than oocytes from younger mares [11]. Alterations of the expression of several genes have been also observed during in vivo and in vitro oocyte maturation [16, 17]. Moreover, after in vitro maturation, higher spindle instability is observed in old versus young mares’ oocytes, resulting in increased risks of chromosome misalignment leading to embryo aneuploidy [18,19,20]. Controversial effects of mare age on mitochondria degeneration have been also reported [16, 21,22,23].

The effect of maternal age on embryo size is not so clear, as studies demonstrated a reduction of embryo diameter in old mares at several developmental stages [24,25,26,27,28] while others did not detect any effect [9, 13, 25, 29] of maternal aging on embryo size. In these studies, however, the developmental stage was not uniform and the method to measure embryo size could be different (direct measurement or using ultrasonography). Furthermore, the size of the equine embryo is highly variable for the same gestational age [30]. The only study that directly measured more than 200 Day 7 to 10 embryos as part of an embryo transfer program observed that maternal aging induced a reduction in embryo size (in average, more than 740 μm in diameter for embryos produced by mares younger than 15 vs less than 620 μm in diameter for embryos produced by mares older than 16) [31].

Maternal age, however, seems to affect embryo morphology from the earliest stages of development, before the embryonic genome activation started (around the 8–16 cells stage) [32, 33]. Zygotes collected in vivo on Day 1.5 post ovulation, indeed, had the same number of cells, but a poorer morphology grade (uniformity of cell size and texture, imperfections and abnormalities) when they had been collected from mares older than 15 years of age compared to young ones [34]. In addition, at Day 3 post ovulation, zygotes collected from mares older than 20 years of age combined poorer morphology grade and reduced cell number [34]. In contrast, in other studies, no difference according to maternal age was observed in terms of number, size or texture of blastomeres of embryos collected in oviducts on Days 2 and 4 post ovulation [6, 10]. Finally, at blastocyst stage, a higher incidence of morphological abnormalities was observed in embryos that were collected from subfertile and/or 18 years old mares compared to 4 years old nulliparous mares [24].

Only a few studies have examined the effect of maternal age on equine embryo quality. As observed in oocytes, metabolic function and capacity were reduced in Day 2 embryos collected from old mares (≥20 years old) compared to those collected from young mares (≤14 years old) [22]. Moreover, a reduced quantity of mitochondrial DNA, but not genomic DNA, has been reported in Day 7 blastocysts from mares older than 16 years of age compared to younger counterparts, suggesting a reduction in mitochondria/cell ratio in old mares’ blastocysts [23]. Taken together, these alterations could lead to the alterations in embryo mobility and/or day of fixation that have been observed in some studies [9, 28].

Altogether, these data indicate that the quality of embryos is reduced in older mares. Considering the Developmental Origins of Health and Diseases (DOHaD) [35], dam age could thus also affect feto-placental and post-natal growth as well as long-term offspring health. To the authors’ knowledge, the effect of maternal age on gene expression in the equine embryo, however, has not yet been explored.

The aim of this study was to determine the effect of maternal age in the mare on embryo gene expression at the blastocyst stage. Young (< 10 years) and older (> 10 years) multiparous mares were inseminated with the same stallion. Blastocysts were collected on Day 8 post-ovulation and bisected to separate the pure trophoblast (TE_part) from the inner cell mass enriched hemi-embryo (ICMandTE). Gene expression was analyzed by RNA-seq in each compartment. Deconvolution was used to isolate gene expression of inner cell mass from trophoblast cells in ICMandTE part.

Results

Embryo recovery rates, diameter, total RNA content and quality

Altogether, 36 uterine flushings were performed (12 in YM < 10 years old and 24 in OM > 10 years old) to obtain one embryo per mare. For 5 mares, after a negative embryo collection, a second attempt was performed. Characteristics of mares are detailed in Table 1. One YM and 3 OM had a double ovulation but no twin embryo was obtained. Post-breeding uterine fluid accumulation was checked at the time of ovulation confirmation and none of the mares presented sign of pathological fluid accumulation.

Table 1 Mares’ characteristics at embryo collection time
Full size table

Embryo recovery rates (ERR) per mare and per ovulation were, respectively, 42 and 39% in YM and 25 and 22% in OM, with no significant difference between groups (p = 0.45 for ERR/mare and p = 0.46 for ERR/ovulation). All embryos were expanded blastocysts, grade I or II according to the embryo classification of McKinnon and Squires [36].

Embryo diameter, sex and RNA yield are presented in Additional Table 1. Embryo diameter ranged from 409 μm to 2643 μm, with no effect of group on embryo diameter (p = 0.97). RNA yield per embryo ranged from 16 ng to 2915.5 ng and was not related to mares’ age (p = 0.92). There is a clear exponential relation between RNA yield/embryo and embryo size (Fig. 1a).

Fig. 1
figure 1

Effect of embryo size on RNA yield and gene expression in equine embryos. a Green circles represent embryos from young mares (YM) and yellow squares represent embryos from old mares (OM). The used relation is exponential. b Biplot graphics representing principal component analyses of differentially expressed genes in inner cell mass and trophoblast of embryos collected in YM or OM. Embryo were grouped according to maternal age (YM or OM) and their diameter (< 1100 μm = small; > 1100 μm = large). The top 20 variables contributing to the axis formation were represented

Full size image

The median RNA Integrity Number (RIN) was 9.8 (8.7–10 range). Between 29.4 and 69.5 million reads per sample were obtained after trimming. On average, 83.61% of reads were mapped on the modified EquCab 3.0 using STAR and 66.66% of reads were assigned to genes by featureCounts.

Deconvolution of gene expression to discriminate ICM and TE gene expression on ICMandTE hemi-embryos

After keeping genes with more than 3 non-zero count values at least in one group (OM or YM) per hemi-embryo (ICMandTE or TE_part), 16,741 genes were conserved for deconvolution. In addition, twenty-seven genes were removed because their variance was null in the TE_part. However, for these genes, the mean counting in ICMandTE samples was above 10 counts, excepted for one which had a mean expression of 34 counts. On the other hand, this gene was only expressed in 5 ICMandTE over the 10 counts threshold.

During the analysis, 3 additional genes were excluded because the deconvolution quality for these genes was not sufficient. Therefore, at the end of deconvolution algorithm, 16,711 genes were available for differential analysis.

Before deconvolution, 341 genes were differentially expressed (FDR < 0.05) between the ICMandTE and the TE_part (Fig. 2a). After deconvolution, comparison between DeMixT_ICM_cells and DeMixT_TE_cells yielded 7569 differentially expressed genes while the comparison DeMixT_ICM_cells vs TE_part yielded 6085 differentially expressed genes, with 5925 in common (77% of the total number of unique DEG in all comparisons). Moreover, except one of the initially 341 differentially expressed genes before deconvolution, were also identified as differentially expressed in both post-deconvolution analyses. On the PCA graphic of individuals, ICMandTE and TE_part were not really separated before deconvolution (Fig. 2b). DeMixT_TE_cells and TE_part were partly superposed, suggesting that datasets before and after deconvolution have a similar global gene expression; whereas the DeMixT_ICM_cells group is clearly separated from both on Axis 2 (14.3% of variance), indicating that the deconvolution effectively enabled the separation of gene expression in the two cell types.

Fig. 2
figure 2

Gene expression in ICM and TE before and after deconvolution using DeMixT. A Venn diagrams of the differential gene expression in ICMandTE vs TE part (before deconvolution), DeMixT_ICM_cells vs DeMixT_TE_cells (after deconvolution) and DeMixT_ICM_cells vs TE_part (gene expression of ICM after deconvolution vs gene expression in TE_part without deconvolution); B Principal Component Analysis of gene expression of DeMixT_ICM_cells, DeMixT_TE_cells, ICMandTE and TE part datasets. Deconvolution was used to isolate gene expression of ICM and TE cells in ICMandTE hemi-embryos. ICMandTE: inner cell mass + trophoblast; TE part: pure trophoblast. Here trophoblast represents trophectoderm + endoderm

Full size image

Six of the 12 genes previously identified as more expressed in ICM [37] were also more expressed in ICMandTE vs TE_part comparison (Table 2). After deconvolution (comparison DeMixT_ICM_cells vs TE_part), 3 more genes previously identified were also differentially expressed. However, for two of them (POU Class 5 Homeobox 1, POU5F1 and Undifferentiated Embryonic Cell Transcription Factor 1, UTF1), their expression was increased in TE_part after deconvolution. In TE, any genes previously identified were differentially expressed in the comparison ICMandTE vs TE_part, i.e., before deconvolution. After deconvolution, the expression of 2 genes was increased in TE_part compared to DeMixT_ICM_cells.

Table 2 Comparison of selected genes expression before and after deconvolution
Full size table

For further analysis, DeMixT_ICM_cells and TE_part gene expression datasets are presented. Nevertheless, the comparison of maternal age on ICMandTE dataset (graphics, differential gene expression analysis and GSEA) was also performed and this analysis is available in Additional Figure 1 and Additional Tables 2 and 3.

Differential gene expression in deconvoluted ICM cells

After filtering out of genes with an average expression < 10 counts/maternal age group/hemi-embryo, 13,735 genes were considered as expressed in the ICM cells from OM or YM embryos including 255 differentially expressed genes (227 downregulated and 28 upregulated in OM, with YM as the reference group) (Fig. 3a and Additional Table 4). Only 218 and 12 genes out of the down- and upregulated genes in OM compared to YM, respectively, were associated to a known protein in human. These downregulated genes in ICM of embryos from OM compared to YM mainly belonged to cellular process, metabolic process and biological regulation in GO biological process gene sets (Fig. 3a). From downregulated genes in OM ICM compared to YM ICM, 84 pathways were statistically overrepresented in GO Biological process leading to 16 most specific subclass pathways according to hierarchy sorting by PANTHER (Fig. 3b). Overrepresented pathways in downregulated genes in OM ICM compared to YM ICM were linked to cell division, embryo morphogenesis, histone methylation, protein modification processes and cytoskeleton organization.

Fig. 3
figure 3

Analysis of differentially expressed genes (DEG) in embryos according to maternal age. A representation of down (blue) and upregulated (red) DEG in ICM (from DeMixT_ICM_cells data obtained after deconvolution of ICMandTE using DeMixT R package [38, 39]) and TE (from TE_part dataset) of embryos from OM. Nota bene: YM is the reference group. Functional classifications (for down regulated genes in ICM and TE and for up-regulated genes in ICM) using GO Biological Process obtained with PANTHER online software are presented as pie charts. For upregulated DEG in TE, no functional classification was made. B Bar chart presenting p-adjusted of most specific subclass pathways provided by PANTHER online software after statistical overrepresentation test (Fisher’s Exact test and Benjamini-Hochberg correction) with the Human GO Biological Process annotation on down-regulated DEGs in ICM of OM. Nota bene: YM is the reference group. DEG: Differentially Expressed Genes (FDR < 0.05); TE: Trophoblast; ICM: Inner Cell Mass; OM: Old mares; YM: young mares

Full size image

Upregulated differentially expressed genes (DEG) in OM ICM compared to YM ICM were mainly part of biological regulation, response to stimulus, cellular process, signaling and metabolic process in GO biological process gene sets (Fig. 3a). No pathway was statistically overrepresented in upregulated genes in ICM of embryos from old mares compared to ICM of embryos from young mares using PANTHER overrepresentation testing.

Differential gene expression in the TE part

In the TE, 13,178 genes were considered as expressed in OM or YM. Fourteen were differentially expressed (Additional Table 5) with 11 genes being downregulated and 3 being up regulated in OM, YM being the reference group (Fig. 3a). Only 12 genes were associated to a known protein in human. Downregulated genes in OM compared to YM mainly belonged to both cellular and metabolic process in GO biological process gene sets.

Analysis of the variability of differentially expressed genes

As the embryo size was highly variable within each group, a PCA analysis was made on DEG in ICM and TE to verify that observed differences were effectively due to maternal age more than to embryo size (Fig. 1b). Respectively 3 YM and 3 OM were classified in small (< 1100 μm in diameter) while the remaining 2 YM and 3 OM were classified in large group (> 1100 μm in diameter).

In the ICM, the first axis of the PCA explained 58.3% of the variance of the data and clearly separated OM from YM, indicating than more than half of the variance in DEG is due to maternal aging. The second axis slightly separated genes according to embryo size but explained only 11.9% of the variance, indicating that in the ICM, embryo size explains a minor part of the variance. These observations were confirmed by the Partial redundancy analysis. Indeed, in ICM, 43% of the variance was explained by maternal age while only 14% was explained by embryo size. Moreover, the 20 variables that contributed the most to the separation are mostly explaining differences related to maternal age.

In the TE part, both the first and the second axis (explaining respectively, 67.3 and 16.5% of data variability) separated both maternal age and embryo size groups. The partial redundancy analysis also showed a less clear impact of maternal age in DEG from TE. Indeed, mares’ age explained 31% of the variance of DEGs whereas 22% was explained by embryo size. Moreover, although most of genes showed a clear contribution in the construction of the first axis, some such as cystathionine beta-synthase (CBS) and CD69 molecule (CD69) were clearly involved in the construction of the second axis.

Gene set enrichment analysis in deconvoluted ICM cells

After Entrez Gene ID conversion, 12,344 genes were considered expressed in ICM. In the GO Biological Process database, 106 identified gene sets differed according to maternal age in the ICM (11 enriched in OM and 95 enriched in YM) (Additional Table 6). In the KEGG database, 14 gene sets were perturbed (Additional Table 6).

Using SUMER analysis (Fig. 4), gene sets enriched in OM ICM compared to YM ICM were mainly grouped under terms “Respiratory electron transport chain” and “Translational elongation”. The first term represents mitochondrial function while the second represents translation function. Moreover, terms relative to “Systemic lupus erythematosus” are conserved by SUMER analysis and enriched in OM ICM compared to YM ICM. Genes involved in these enrichments were also present in dysregulated metabolic pathways.

Fig. 4
figure 4

GSEA clustering of the most perturbed terms in ICM and TE according to mares’ age. Nodes represent altered gene sets in ICM and TE (FDR < 0.05). Node size represents the gene set size. Node shape represents the gene set database: GO BP (circle) or KEGG (diamond). YM was the reference group: blue nodes represent enriched gene sets in YM (NES < 0) while red nodes represent enriched gene sets in OM (NES > 0). Edges represent the level of connection between representative gene sets. This graph was performed using SUMER R package [40] and modified using cytoscape 3.8.2 [41]. ICM: Inner Cell Mass; TE: trophoblast; FDR: False Discovery Rate; GO BP: Gene Ontology Biological Process; Kyoto Encyclopedia of Genes and Genomes; NES: Normalized Enrichment Score; YM: Young mares; OM: Old mares

Full size image

Gene sets enriched in YM ICM compared to OM ICM were involved in different biological processes such as “Mitotic nuclear division” which is relative to mitosis. Under this term, “Spindle organization” and “Negative regulation of supramolecular fiber organization” terms represent the assembly of spindle during mitosis. “Chromosome segregation” is also found in this grouping. Furthermore, terms under “Golgi vesicle transport” seem to be involved in cytokinesis and vesicle trafficking. The catabolism of protein and mRNA is represented by “Ubiquitin mediated proteolysis”. Translation function is represented by several terms grouped under “covalent chromatin modification” in the mitosis part.

“Leukocyte transendothelial migration” term is also enriched in YM ICM compared to OM ICM and seems to represent cell signaling and adhesion. Indeed, this term regroups several pathways such as “Focal adhesion”, “Adherens junction” and “WNT signaling pathway”.

Finally, further terms related to embryo development could be found in different groups such as “Dendrite development”, “Cell differentiation in spinal cord” or “endothelial cell development”.

Gene set enrichment analysis in the TE part

After Entrez Gene ID conversion, 11,874 genes were considered expressed in TE from OM or YM embryos. In GO BP database, 38 gene sets were perturbed by mares’ age in the TE of embryos (11 enriched in YM and 27 enriched in OM) (Additional Table 7). In the KEGG database, 5 pathways were enriched in the TE of OM compared to TE from YM embryos (Additional Table 7).

As in ICM part, using SUMER analysis, translation (represented by both “Translational termination” and “Ribosome” terms) and mitochondrial function (represented by “Oxidative phosphorylation” associated terms) were enriched in OM embryos compared to YM embryos. “Proteasome” term is enriched in OM embryos, compared to YM representing protein catabolism.

Moreover, transcription (represented by “Peptidyl lysine trimethylation” term) and embryo development (represented by “spinal cord development”) were, as in ICM, also, enriched in YM embryos compared to OM. Finally, in TE, “Regulation of action potential” term and associated terms were enriched in YM compared to OM. These terms could be related to ion movement.

Discussion

Maternal aging did not affect embryo recovery rate (ERR), embryo size nor total RNA quantity. It is important to notice that, as the experimental protocol required only 5 embryos per group, a limited number of embryo collection attempts were performed, resulting in a limited statistical power for these analyses. ERR, especially in young mares, was not as good as in commercial embryo transfer programs [31, 42]. Moreover, total RNA content in equine embryo was exponentially correlated to embryo diameter and was not affected by maternal age.

As ovulation check was only performed 48 hours after induction of ovulation, there is a bias in the estimation of the embryo age in the present study. This certainly partly explains the large variation in embryo size observed here. Different studies, however, have shown that embryo size is highly variable [43,44,45], even when ovulation determination is accurate to the nearest 6 hours [46]. Accurately determining fertilization time in vivo is not possible and the size of 142 equine embryos recovered at 8 days post-ovulation was reported to range from 120 to 3980 μm with an average of 1132 μm [30]. In the present project, all embryos remain within these limits.

Although embryo size could be considered as normal for Day 8 embryos, embryo sizes were highly variable and not equally distributed between groups. Embryo size was considered as a cofactor in the differential analysis. In addition, partial redundancy analysis showed that embryo size explained less than 25% of the variability of DEGs in both ICM and TE. As effects of maternal age on equine embryo growth are not so clear [24,25,26,27,28], an interaction between both effects on embryo gene expression can not be excluded. Maternal age, however, explained 43 and 31% of the total variability of DEGs in both ICM and TE, respectively, suggesting that identified DEGs were mainly due to maternal age effects. It is important to note that in ICM, the proportion of variability explained by maternal age was 3 times greater than that of embryo size, whereas in the TE, this difference was smaller. This could be explained by the fewer number of DEGs identified in TE compared to ICM but also because trophoblast expansion was possibly more correlated to embryo size than the ICM.

As in tumors, today, the only way to analyze separate expression of ICM and TE cells in mammalian embryos is to use micromanipulation but micromanipulation is expensive in time and requires adapted skills and material. In this study, micromanipulation to properly separate ICM and TE was not possible. Therefore, deconvolution appeared to be a good option to analyze the effects of maternal ICM gene expression in mixed ICM and TE samples. Deconvolution seemed to have been effective on embryo to separate ICM and TE gene expression. Indeed, the 341 genes differentially expressed in ICMandTE vs TE_part were also identified by comparing TE with post-deconvolution ICM cells gene expression dataset. Moreover, after deconvolution, more genes were differentially expressed between TE and ICM datasets and among them, more genes identified in literature were differentially expressed. Finally, PCA showed that deconvolution separate DeMixT_ICM_cells from other datasets while before deconvolution ICMandTE and TE_part were partly superposed. Only POU5F1 and UTF1 differed in terms of allocation between previously published data [37] and the present study. In data before deconvolution, expression of POU5F1 and UTF1 were extremely variable according to embryo diameter in both compartments. Deconvolution reduced the variability of both gene expression in ICM but not in TE. Studies showed that POU5F1 is expressed in both compartments until Day 10 post ovulation in equine in vivo produced embryos [47] but its expression decreases progressively from Day 7 in TE cells. Thus, the statistically increased expression of POU5F1 in TE on Day 8 is compatible with previously published data.

Therefore, it was chosen to only present the analysis of the effects of age in DeMixT_ICM_cells dataset. In all cases, deconvolution did not change the principal results of maternal aging on embryo gene expression because the comparison of maternal age in ICMandTE and DeMixT_ICM_cells gave comparable results.

Genes with higher expression in both the ICM and the TE of OM embryos were overrepresented for functional terms related to mitochondrial function and translation with ribosome biogenesis. Those with higher expression in YM embryos were involved in transcription and chromatin modification (Fig. 5). Mitosis, signaling and adhesion pathways were additionally altered and the expression of genes involved in embryo development was particularly reduced in the ICM hemi-embryos from old mares. Finally, in the TE, ion movement related genes were affected.

Fig. 5
figure 5

Schematic representation of the observed effects of maternal aging on ICM and TE in equine blastocysts. Normal blastocyst development from day 7 to day 9 is represented inside black square using literature [48, 49]. Light blue boxes represent biological processes that are enriched in YM embryos and light red boxes represent biological processes that are enriched in OM embryos

Full size image

Several studies on mature oocytes in women and mice, using microarrays and RNAseq, had shown a differential expression of genes involved in mitochondrial, cell cycle, cytoskeleton functions and related pathways according to maternal age [50,51,52], as observed in the present study in equine blastocysts, especially in the ICM. In equine follicles and oocytes, the expression of genes essential for oocyte maturation, as studied by RT-qPCR, quantitatively and temporally differed with maternal age, suggesting a desynchronization [16]. Moreover, a study on the expression of 48 genes that have been previously identified as related to aneuploidy and maternal aging in humans showed differential expression in cumulus-oocytes complexes from young (< 12 years of age) and old (> 18 years old) mares [17]. Genes related to mitochondrial replication and function were also down-regulated in oocytes of old mares [23]. Even effects of maternal aging on gene expression in equine oocytes, using high throughput techniques, are not yet available, results on equine Day 8 blastocysts suggest that mitochondrial perturbations present in the oocyte persist until the blastocyst stage. In humans, single-embryo RNA-seq analysis showed that maternal aging alters the transcriptome of blastocysts derived from ICSI. Genes altered by maternal age were related to segregation of chromosomes, embryonic growth and cell signaling [53, 54], as observed in our study on equine blastocysts.

Unlike the human embryo, the equine embryo remains spherical, grows very fast and is surrounded by its glycoprotein capsule up to 22 days post ovulation [43]. At 8 days post ovulation, equine embryos are therefore growing and the endodermal layer entirely delineates the future yolk sac (Fig. 5) [48]. Moreover, whereas in Day 7 embryos, the epiblast is lined by a continuous layer of trophoblast named the polar trophoblast. As the blastocyst grows, epiblast cells phagocyte polar trophoblast cells that become discontinuous. During this time, epiblast and polar trophoblast cells divided and are bound by tight junctions [49]. The turnover in focal adhesion is essential for cell movement and phagocytosis [55].

Cell signaling and adhesion function were altered in OM embryos compared to YM. Indeed, KEGG pathways “focal adhesion” and “adherens junction” regrouped in “Leukocyte transendothelial migration” on Fig. 3 are disturbed in old mares’ embryos. Focal adhesions are subcellular structures observed in adherent cells, that are clusters of different structural and signaling molecules and integrins that form a physical link between the actin cytoskeleton and the extracellular matrix [56, 57]. Integrins are the major receptors for cell adherence to extracellular matrix proteins in focal adhesions by playing an important role in transmembrane connections to the cytoskeleton and activation of several signaling pathways, as growth factor receptors [58]. Integrin Subunit Alpha V (ITGAV) expression was reduced in the ICM of embryos collected from OM compared to the ones from YM. In different human breast cancer cell lines, cell proliferation, invasion and renewal are restricted when ITGAV is inhibited [59].

Moreover, the activity and localization of small GTPases within the Rho family must be precisely regulated as the communication between them is indispensable to achieve the turnover of actin filaments and cell adhesions necessary for cell migration, that are essential processes for the proper development of the embryo [60]. The Rho GTPase Activating Protein 9 and 35 (ARHGAP9 and ARHGAP35) were downregulated in OM ICM compared to YM ICM. ARHGAP35 codes for P190 Rho GTPase activating protein that mediates signalization through integrin-dependent adhesion in cells. Aberrant neural morphogenesis, related to an accumulation of polymerized actin in embryoneural tube cells, is reported in mice lacking P190 Rho GTPase activating proteins [61]. Furthermore, P21 activated kinases (PAK) are effector kinases for other small Rho GTPases such as Rac that regulate the polymerization of actin and can affect the microtubule organization [62]. Once activated, PAKs promote the turnover of focal adhesions by regulating transcription [63]. PAK3 was downregulated in the OM compared to YM ICM. However, in TE, SLIT-ROBO Rho GTPase Activating Protein 3, SRGAP3 was one of the upregulated genes in OM compared to YM. This gene regulates cytoskeletal reorganization and plays an important role in the formation of cell-cell adhesion [64]. Therefore, maternal aging seems to alter focal adhesion in both ICM and TE of embryos but, while in ICM, it seems to reduce cohesion and signalization, in TE it seems to reinforce them.

Several signaling pathways required for development were altered and numerous genes related to focal adhesions and signaling pathways were down-regulated in OM ICM compared to YM (ARHGAP9 and ARHGAP35; B-Raf Proto-Oncogene, Serine/Threonine Kinase, BRAF; Filamin C, FLNC; FYN Proto-Oncogene, Src Family Tyrosine Kinase, FYN; PAK3; Protein Phosphatase 1 Regulatory Subunit 12A, PPP1R12A; Talin 2, TLN2; X-Linked Inhibitor Of Apoptosis, XIAP). Among all signaling pathways, transforming growth factor beta (TGF-β), Wg et Int (Wnt), mitogen-activated protein kinases (MAPK) and apoptosis signaling were altered by maternal age in equine embryos. Indeed, transforming growth factor beta receptor 1 (TGFBR1) was downregulated in the ICM of embryos from old mares compared to the ones from young mares. TGFBR1 is one of the two dimers for the receptor for transforming growth factors beta that is essential in the generation of axes and cell fate commitment during mammal embryogenesis [65]. Moreover, Casein Kinase 1 Alpha 1, (CSNK1A1) and Protein Phosphatase 2 Regulatory Subunit B’Epsilon (PPP2R5E), involved in Wnt signaling pathway, were downregulated in OM embryos compared to YM. Wnt signaling pathway is also involved in nervous development system and may act on cell determination in mice embryos [66]. Moreover, Doublesex And Mab-3 Related Transcription Factor 3 (DMRT3) was downregulated in the OM ICM compared to YM. DMRT3 is required for the proper development of the brain in prenatal chicken and mice [67].

Cell cycle related pathways and genes were also disturbed by maternal aging. Indeed, the expression of Pleckstrin Homology domain-Interacting Protein (PHIP), also known as replication initiation determinant, was reduced in the OM ICM compared to YM. This gene is required for the initiation of DNA replication. Fibroblasts originated from mice embryos with depleted PHIP expression presented less replication-initiation as well as abnormal replication fork progression events [68]. Moreover, several cyclins (CCN) were less expressed in the ICM of embryos from old mares compared to the ones from young mares (CCNB3, CCNI and CCNT2). Cyclins are key cell cycle regulators and, for instance cyclin B3 degradation allows the passage from metaphase (chromatids aligned on metaphase plate) to anaphase (segregation of sister chromatids) [69]. During metaphase and anaphase, spindle assembly and chromosome segregation could also be compromised. Indeed, genes involved in the formation of microtubules (Tubulin Gamma 2, TUBG2) [70], the generation of the mitotic spindle (Never in Mitosis A related kinase 9, NEK9) [71], the production of cohesion and condensing complexes (Structural Maintenance Of Chromosomes 4, SMC4) [72], the control of centriole length, important for centrosomes and microtubules organization (Centriolar Coiled-Coil Protein 110, CCP110) [73], the amplification of intra-spindle microtubules allowing the formation of a sufficient quantity of kinetochore (HAUS Augmin Like Complex Subunit 6, HAUS6) [74], the stabilization and protection of the cohesin complex association with chromatin (Cell Division Cycle Associated 5, CDCA5; Shugosin 2, SGO2) [75, 76] and the production of microtubule-based motor proteins (Kinesin Family Member 14 and 23, KIF14 and KIF23) [77] were downregulated in OM ICM compared to YM. Furthermore, Centromere protein F (CENPF), an important regulator of chromosome segregation during mitosis, was reduced in OM ICM compared to YM. CENPF is indispensable for murine embryo development from the first stages of development until at least the blastocyst stage [78]. Impaired CENPF expression could induce chromosome mis-segregation in the ICM of old mares’ embryos compared to the ones from young mares.

Cytokinesis seems to also be altered by maternal age. CD2-associated protein (CD2AP) and Vacuolar Protein Sorting 4 Homolog A (VPS4A) were downregulated in OM ICM compared to YM. Both are involved in the late phase of cytokinesis which correspond to the abscission of the cell, producing two daughter cells [79, 80]. Globally, results suggest that cell division is impaired or indicate that there is less division in the ICM of old mares’ embryos compared to the ones from young mares.

Transcription is regulated by epigenetic marks. Several pathways linked to methylation and acetylation were enriched in YM ICM and TE in relation to numerous upregulated genes in YM ICM compared to OM (i.e., Histone Deacetylase 1 and 2, HDAC1 and HDAC2; Helicase, Lymphoid Specific, HELLS; DNA Methyltransferase 3 Alpha, DNMT3A; Lysine Acetyltransferase 6B, KAT6B; Lysine Demethylase 4A, KDM4A; Lysine Methyltransferase 2A, KMT2A; Nuclear Receptor Binding SET Domain Protein 3, NSD3; SET Domain Containing 5, SETD5; Tet Methylcytosine Dioxygenase 1, TET1). At blastocyst stage, de novo methylation is more important in the inner cell mass than in trophoblast cells and is essential for the proper development of mammalian embryos [81]. DNMT3A is essential for de novo methylation and required for embryo viability [82]. As it has been shown that maternal aging increases the acetylation of lysine 12 of histone H4 in mice oocytes, affecting fertilization and embryo development [83], it is likely that the methylation and acetylation dysregulation observed in the embryos of old mares could be inherited from the oocyte rather than directly due to embryo environment.

Both in OM ICM and TE compared to YM, several pathways related to mitochondrial function were enriched. Mitochondria are long known to be the “powerhouse” of cells, producing most of the necessary ATP for cellular function [84]. ATP production is possible through oxidative phosphorylation enabled by the electron transport chain in the inner membrane of mitochondria. This process begins with the oxidation of NADH or FADH that generates a proton gradient across the mitochondrial membrane, driving ATP synthesis through the phosphorylation of ADP via ATPase [85]. Enriched pathways in OM compared to YM embryos demonstrate that the NADH dehydrogenase (complex I), the oxidoreductase acting on a heme group of donors (complex III) with the global electron carrier, and the proton-transporting ATP synthase (1rst part of the complex V) activities were particularly enriched while the second part of complex V is depleted. Since all complexes are required to produce ATP, the impaired complex generation may lead to mitochondrial stress and uncoupled respiratory chain, resulting in poor ATP production. In human IVF embryos, reduced mitochondrial respiration efficiency, with a negative effect on embryo development, is observed when oocytes are collected from old women [86]. In mammalian expanded blastocysts, energy requirement is high, probably because of the need for ion transport enabling quick blastocoel cavity expansion. Defects in metabolic requirements during blastocyst development are associated with impaired implantation [87]. Equine embryos are not an exception [88].

Mitochondria being maternally inherited, uncoupling of the mitochondrial respiratory chain may originate from mitochondrial defects already present in oocytes [89]. Studies have shown that, in human, poor oocyte quality is responsible for the decline in fertility associated with maternal aging. The reduced oocyte developmental competence in older women is related to mitochondrial dysfunction [90]. In horses’ oocytes, before in vitro maturation, no difference between young and old mares is observed while after, reduced mitochondrial DNA copy numbers, increased swelling and reduced mitochondrial cristae are reported in oocytes collected from mares > 12 years. This suggests that oocyte mitochondria from old mares present deficiencies that become problematic when energy needs are high, as during in vitro maturation or subsequent embryo development [21]. Because mitochondrial replication in equine embryos does not begin until the blastocyst stage [91], oocyte mitochondrial numerical or functional deficiencies may compromise early embryogenesis.

In TE, the GO biological process “Regulation of action potential” and “Transmission nerve impulse”, related to ion movement, were enriched in YM compared to OM embryos. These processes are associated with genes involved in the action potential dependent ion channels. These ion channels such as NA+, K+, ATPase are essential for the accumulation of blastocoel and subsequently yolk sac fluid in equine embryos [92]. Their inhibition is detrimental for embryo development [93]. Thus, both reduced ATP production and disturbed ion channels’ function may impact embryo growth, although effects of mares’ age on embryo growth are controversial [24,25,26,27,28].

Conclusions

This is the first study showing an effect of maternal age on gene expression in the equine blastocyst at Day 8 post ovulation. Maternal age, even for mares as young as 10 years old, disturbs several processes related to chromatin segregation, cytokinesis, apposition of epigenetic marks and mitochondrial function and that the ICM appears to be more affected than the TE. These perturbations on Day 8 post ovulation may affect the further development of the embryo and as a result, may contribute to decreased fertility due to aging. So far, maternal age does not appear to have an impact on foal morphology at birth but, as for placental function, more subtle effects may exist. Here, as embryos were collected long before implantation, it is impossible to know if these embryos would have implanted and if they would have produced a live foal or not. Studies on long-term growth and health of foals born to the same mares are currently being performed.

Methods

Ethics

The experiment was performed at the experimental farm of IFCE (research agreement C1903602 valid until March 22, 2023). The protocol was approved by the local animal care and use committee (“Comité des Utilisateurs de la Station Expérimentale de Chamberet”) and by the regional ethical committee (“Comité Régional d’Ethique pour l’Expérimentation Animale du Limousin”, approved under N° C2EA - 33 in the National Registry of French Ethical Committees for animal experimentation) under protocol number APAFIS#14963–2018050316037888 v2. All experiments were performed in accordance with the European Union Directive 2010/63EU. The authors complied with the ARRIVE guidelines.

Embryo collection

Thirty-one non-nursing multiparous Saddlebred mares (French Anglo-Arabian and Selle Français breeds) were included in this study. Multiparous mares were defined as dams that had already foaled at least once. Mares were bred and raised in the “Institut Français du Cheval et de l’Equitation” experimental farm (Chamberet, France, 45° 34′55.17″N, 1°43′16.29″E, 442 m). During the experimental protocol, mares were managed in herds in natural pastures 24 h/day with free access to water. The experiments took place from April 1st to May 23rd, 2019. All mares remained healthy during this period.

Mares were allocated to 2 groups according to their age: young mares (YM, < 10 years, n = 11) and older mares (OM, ≥ 10 years, n = 20). As several studies had shown an increase embryo loss and/or reduced pregnancy rates as early as 10 years old, the age threshold was defined in accordance with the literature [3, 4]. Before experimentation, mare’s withers’ height was measured. During the experimentation, mares were weighted and mares’ body condition score (BCS, scale 1 to 5) was determined by one trained operator [94].

The mares’ estrous period was monitored routinely by ultrasound with a 5 MHz trans-rectal transducer. During estrus, ovulation was induced with a single injection of human chorionic gonadotropin (i.v.; 750 - 1500 IU; Chorulon® 5000; MSD Santé animale, France) as soon as one ovarian follicle > 35 mm in diameter was observed, together with marked uterine edema. Ovulation usually takes place within 48 h, with > 80% occurring 25 to 48 h after injection [95, 96]. At the same time, mares were inseminated once with fresh or fresh overnight cooled semen containing at least 1 billion motile spermatozoa from a single fertile stallion. Ovulation was confirmed within the next 48 hours by ultrasonography.

Embryos were collected by non-surgical uterine lavage using prewarmed (37 °C) lactated Ringer’s solution (B.Braun, France) and EZ-Way Filter (IMV Technologies, France) 10 days after insemination, i.e. approximately 8 days post ovulation. All embryo collections were performed by the same operator, properly trained but not performing embryo collection routinely. Just after embryo collection, mares were treated with luprotiol an analogue of prostaglandin F2α (i.m; 7.5 mg; Prosolvin, Virbac, France).

The aim of the embryo collection was to obtain 5 embryos/group with each embryo coming from a different mare. Therefore, some mares that failed to produce an embryo at their first attempt had been bred again.

Embryo bisection and RNA extraction

Collected embryos were immediately observed using a binocular magnifier and photographed with a size standard. Embryo diameter was subsequently determined on photographs using ImageJ® software (version 1.52a; National Institutes of Health, Bethesda, MD, USA). Embryos were washed 4 times in commercially available Embryo holding medium (IMV Technologies, France) at 34 °C and then bisected with a microscalpel under binocular magnifying glass to separate trophoblast part, at this stage composed of trophectoderm and endoderm (TE_part) from inner cell mass (composed of epiblast cells) enriched hemi-embryos (ICMandTE) (Fig. 6). Each hemi-embryo was transferred to an Eppendorf DNA LoBind tube (Eppendorf, Montesson, France) with extraction buffer (PicoPure RNA isolation kit, Applied Biosystems, France) and incubated for 30 min at 42 °C prior to storage at − 80 °C.

Fig. 6
figure 6

Bisection of equine embryos at 8 days post-ovulation into an ICMandTE and a TE_part. a Step 1: Identification of the ICM and the TE; b Step 2: Definition of the cutting plane to isolate the ICM in one of the parts; c Step 3: Cutting of the embryo and separation of the two parts. ICMandTE: inner cell mass + trophoblast; TE part: pure trophoblast. NB: Here the trophoblast represents trophectoderm + endoderm and ICM is composed of epiblast cells

Full size image

RNA was extracted from each hemi-embryo following the manufacturer’s instructions using PicoPure RNA isolation kit (PicoPure RNA isolation kit, Applied Biosystems, France), which included a DNAse treatment. RNA quality and quantity were assessed with the 2100 Bioanalyzer system using RNA 6000 Pico or nano kit (Agilent Technologies, France) according to the manufacturer’s instructions (Additional Figure 2).

RNA sequencing

Five nanograms of total RNA were mixed (Additional Figure 2) with ERCC spike-in mix (Thermofisher Scientific, France) according to manufacturer recommendations. Messenger RNAs were reverse transcribed and amplified using the SMART-Seq V4 ultra low input RNA kit (Clontech, France) according to the manufacturer recommendations. Nine PCR cycles were performed for each hemi-embryo. cDNA quality was assessed on an Agilent Bioanalyzer 2100, using an Agilent High Sensitivity DNA Kit (Agilent Technologies, France). Libraries were prepared from 0.15 ng cDNA using the Nextera XT Illumina library preparation kit (Illumina, France). Libraries were pooled in equimolar proportions and sequenced (Paired end 50–34 pb) on NextSeq500 instrument, using a NextSeq 500 High Output 75 cycles kit (Illumina, France) (Additional Figure 2). Demultiplexing was performed with bcl2fastq2 version 2.2.18.12 (Illumina, France) and adapters were trimmed with Cutadapt version 1.15 [97]. Only reads longer than 10pb were kept.

Embryo sexing

Embryo sex was determined for each embryo in order to take in consideration embryo sex in the statistical analysis as sexual dimorphism was observed in blastocysts in other species [98].

The expression of X Inactive Specific Transcript (XIST) was analyzed. XIST is a long noncoding RNA from chromosome X responsible for X inactivation when two X chromosomes are present, meaning that only female embryos express XIST. XIST was shown to be expressed in equine female blastocysts on Day 8 post ovulation [37, 44]. XIST is not annotated in the available equine Ensembl 3.0.99 or NCBI annotation. Thus, a de novo transcript assembly was made using StringTie version 2.1.2 [99] following manufacturer instructions to identify XIST in equine embryo samples (Additional Figure 2). The genome annotation file output was aligned on Ensembl 99 EquCab3.0 assembly using IGV version 2.7.2 [100]. As XIST position and length on X chromosome has been previously determined [37] and as XIST is annotated in humans, the most similar contig was selected. The obtained XIST sequence was 22,600 bp long and composed of 8 exons and 7 introns. XIST sequence was added to FASTA files from Ensembl 99 EquCab3.0 assembly and to the genome annotation file from Ensembl 99 EquCab3.0 assembly as a unique fictitious chromosome with one gene composed of one transcript.

XIST expression differed between individuals but was relatively close in quantification between related hemi-embryos. Therefore, seven embryos were determined as females (4 in the YM group and 3 in the OM group) while 4 were considered as males (1 in the YM group, and 3 in the OM group).

RNA mapping and counting

Alignment was performed using STAR version 2.6 [101] (2-pass mapping with the genome annotation file using the following parameters: --alignSoftClipAtReferenceEnds No --alignEndsType EndToEnd --alignIntronMax 86,545 --alignMatesGapMax 86,545) on previously modified Ensembl 99 EquCab3.0 assembly and annotation. Genes were then counted with featureCounts [102] from Subreads package version 1.6.1 (parameters: -p -B -C -M) (Additional Figure 2).

Data analysis

All statistical analyses were performed by comparing OM to YM (set as reference group) using R version 4.0.2 [103] on Rstudio software version 1.3.1056 [104].

Embryo recovery rate, embryo diameter and total RNA content analysis

Embryo recovery rate (ERR) per mare was calculated as number of attempts with at least one embryo collected/total number of attempts. ERR per ovulation was determined as number of collected embryos/total number of observed ovulations. They were analyzed using the Exact Fisher test to determine if age influences ERR.

For total RNA content analysis, as embryos were bisected without strict equality for each hemi-embryo, a separate analysis of ICMandTE and TE_part RNA quantities would not be meaningful. Thus, ICMandTE and TE_part RNA quantities were initially summed up and analyzed using a mixed linear model from nlme package version 3.1–148 [105], followed by 1000 permutations using PermTest function from pgirmess package version 1.6.9 [106]. Maternal age and embryo sex were considered as random effects.

Embryo diameter was analyzed with a fixed effects linear model of nlme package version 3.1–148 [105] including maternal age and embryo sex, followed by 1000 permutations using PermTest function from pgirmess package version 1.6.9 [106]. Variables were kept in models when statistically significant differences were observed. Differences were considered as significant for p < 0.05.

Moreover, the relation between embryo size and total RNA quantity/embryo was studied and represented using the nonlinear regression analysis (fitness to exponential growth equation) in GraphPad Prism software 8.0.1 for Windows (Graphpad Software, San Diego, California USA, www.graphpad.com).

Deconvolution of gene expression in ICMandTE using DeMixT

A filter on counts obtained with featureCounts was applied: all genes with more than 3 non null count values in at least one group (OM or YM) per hemi-embryo (ICMandTE or TE_part) were conserved. The ICMandTE hemi-embryo was composed of both trophoblast and inner cell mass in unknown relative proportions. In order to estimate the relative gene expression of both cell types, the DeMixT R package version 1.4.0 was used [38, 39]. Starting from a dataset obtained from one cell type, DeMixT estimates the proportion of cells in the mixed samples and performs a deconvolution algorithm on gene expression. TE_part (reference cells) and ICMandTE (mixed cells type) datasets were thus analyzed using the following options: niter = 10; nbin = 50; if.filter = TRUE; ngene.selected.for.pi = 1500; ngene. Profile.selected = 1500; nspikein = 0. Since the deconvolution algorithm cannot be performed in the presence of null values, the value “1” was added to all count values in ICMandTE and TE_part data tables. Moreover, variance in reference cell dataset (TE_part) must not be null and gene with null variance had been removed prior to the deconvolution. Output datasets were DeMixT_ICM_cells and DeMixT_TE_cells, corresponding to the deconvoluted gene expression in ICM cells and TE cells of ICMandTE, respectively.

Quality checks for deconvolution

At the end of deconvolution, a quality check is automatically performed by the analysis: for every gene, if the difference of average of deconvoluted expression of reference cells in mixed samples (here, gene expression of TE cells in ICMandTE dataset, i.e., DeMixT_TE_cells) and average of expression of reference cells (here, gene expression of TE cells in TE_part) is > 4, then, the deconvolution for concerned genes is not reliable. Therefore, genes in this case are filtered out.

Moreover, outputs of DeMixT_ICM_cells vs DeMixT_TE_cells, DeMixT_ICM_cells vs TE_part and ICMandTE vs TE_part were compared with Deseq2 version 1.28.1 [107] to confirm that the deconvolution was effective at separating gene expression. Moreover, several genes proper to ICM and TE cells in equine embryos were identified using literature search [37]. Their raw and deconvoluted expressions were analyzed to check the reliability of deconvolution. Results of these analyses were represented through manually drawn Venn diagrams as well as principal component analysis (PCA) graphics of individuals, using FactoMineR version 2.4 [108], ggplot2 version 3.3.3 [109] and factoextra version 1.0.7 [110].

Maternal age comparison for gene expression

For the comparison according to maternal age, ICMandTE, DeMixT_ICM_cells and TE_part datasets were used and a filter was applied on datasets: all genes with an average expression equal or above 10 counts in at least one group (OM or YM) per hemi-embryo (ICM or TE) were conserved. Differential analyses were performed with Deseq2 version 1.28.1 [107] with YM group as reference, without independent filtering. Genes were considered differentially expressed (DEG) for FDR < 0.05 after Benjamini-Hochberg correction (also known as false discovery rate, FDR). As ovulation was checked only every 48 h and because embryos growth is exponential in the uterus, embryo diameter was considered as a cofactor in the model as well as embryo sex. Moreover, PCAs using factominer version 2.4 [108] on R software followed with a Partial redundancy analysis using rda function from vegan package version 2.5–7 [111] were performed. For these analyses, embryo size was considered as a factor of two modalities: embryos < 1100 μm were considered as small and > 1100 μm were considered as large. Results were represented with a biplot graphic using the package factoextra version 1.0.7 [110]. The combined analyses were used to represent and determine the specific contribution of each of the two factors in the total variability observed in DEGs. The objective was to verify that the age of the mare was the main factor explaining the difference between the DEGs and not the differences in embryo size. On biplot graphics, only the 20 variables, i.e., genes, that contributed the most in axis construction, were represented.

Equine Ensembl IDs were converted into Human Ensembl IDs and Entrez Gene names using gorth function in gprofiler2 package version 0.1.9 [112]. Genes without Entrez Gene names using gprofiler2 were manually converted when Entrez Gene names were available, using Ensembl web search function [113].

Moreover, DEG were analyzed for classification, over- and under- representation (FDR < 0.05) using Fisher’s exact test corrected with FDR and using the database Gene Ontology (GO) Biological Process complete on PANTHER version 16.0 web-based software [114]. Only the most specific subclass pathways, provided by PANTHER’s automated hierarchy sorting, were conserved and represented.

Gene set enrichment analyses (GSEA)

The count values were log transformed using RLOG function of DESeq2 version 1.28.1. Gene set enrichment analyses (GSEA) were performed on expressed genes using GSEA software version 4.0.3 (Broad Institute, Inc., Massachusetts Institute of Technology, and Regents of the University of California) [115, 116] to identify biological gene sets disturbed by age. GSEA was performed using the following parameters: 1000 gene set permutations, weighted enrichment analysis, gene set size between 15 and 500, signal-to-noise metrics. Molecular Signatures Database [117] version 7.1 (C2: KEGG: Kyoto Encyclopedia of Genes and Genomes, C5: BP: GO biological process) were used to identify most perturbed pathways. Pathways were considered significantly enriched for FDR < 0.05. YM group was considered as the reference group, which indicates that when the normalized enrichment score (NES) was positive, the gene set was enriched in the OM group while when NES was negative, the gene set was enriched in the YM group.

Enriched terms from GO BP and KEGG database were represented using SUMER analysis from SUMER R package version 1.1.5 and using FDR q-values [40]. First, SUMER analysis reduces redundancy using weighted set cover algorithm and then clusters the results using affinity propagation algorithm. The weighted set cover algorithm performs a selection of the fewest number of gene sets that include all genes associated with the enriched sets, giving priority for the most significant terms. Then, the affinity propagation algorithm groups similar gene sets and identifies one representative gene set. This algorithm will repeat this clustering to obtain the most representative gene sets. Results were represented with graph modified using Cytoscape version 3.8.2 [41]. Gene sets were represented by nodes and the gene set size was represented by the size of the node. Node shape represented the gene set database (GO BP or KEGG). Blue nodes represented gene sets enriched in YM (NES < 0) while red nodes represented gene sets enriched in OM (NES > 0). Edge width represents the level of connection between representative gene sets (thinner edges represent the first clustering while thicker edges represent the second clustering of the affinity propagation algorithm).

Data availability statement

The RNA sequençing data supporting the conclusions of this article are available in the repository [accession: GSE162893; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE162893] from the GEO SuperSeries [accession: GSE193676; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE193676].

Availability of data and materials

The RNA sequençing data supporting the conclusions of this article are available in the repository [accession: GSE162893; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE162893] from the GEO SuperSeries [accession: GSE193676; https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE193676].

References

  1. Scoggin CF. Not just a number: effect of age on fertility, pregnancy and offspring vigour in thoroughbred brood-mares. Reprod Fertil Dev. 2015;27(6):872.

    Article  PubMed  Google Scholar 

  2. Allen WR, Brown L, Wright M, Wilsher S. Reproductive efficiency of Flatrace and National Hunt Thoroughbred mares and stallions in England. Equine Veterinary J. 2007;39(5):438–45.

    Article  CAS  Google Scholar 

  3. Langlois B, Blouin C. Statistical analysis of some factors affecting the number of horse births in France. Reprod Nutr Dev. 2004;44(6):583–95.

    Article  PubMed  Google Scholar 

  4. Malheiros de Souza JR, Dias Gonçalves PB, Bertolin K, Ferreira R, Sardinha Ribeiro AS, Ribeiro DB, et al. Age-dependent effect of foal heat breeding on pregnancy and embryo mortality rates in Thoroughbred mares. J Equine Veterinary Sci. 2020;5:102982.

    Article  Google Scholar 

  5. Chevalier-Clément F. Pregnancy loss in the mare. Anim Reprod Sci. 1989;20(3):231–44.

    Article  Google Scholar 

  6. Ball BA, Little TV, Hillman RB, Woods GL. Pregnancy rates at days 2 and 14 and estimated embryonic loss rates prior to day 14 in normal and subfertile mares. Theriogenology. 1986;26(5):611–9.

    Article  CAS  PubMed  Google Scholar 

  7. Doig PA, McKnight JD, Miller RB. The use of endometrial biopsy in the infertile Mare. Can Vet J. 1981;22:72–6.

    CAS  PubMed  PubMed Central  Google Scholar 

  8. Siemieniuch MJ, Gajos K, Kozdrowski R, Nowak M. Advanced age in mares affects endometrial secretion of arachidonic acid metabolites during equine subclinical endometritis. Theriogenology. 2017;103:191–6.

    Article  CAS  PubMed  Google Scholar 

  9. Carnevale EM, Ginther OJ. Relationships of age to uterine function and reproductive efficiency in mares. Theriogenology. 1992;37(5):1101–15.

    Article  CAS  PubMed  Google Scholar 

  10. Ball BA, Little TV, Weber JA, Woods GL. Survival of Day-4 embryos from young, normal mares and aged, subfertile mares after transfer to normal recipient mares. Reproduction. 1989;85(1):187–94.

    Article  CAS  Google Scholar 

  11. Brinsko SP, Ball BA, Ellington JE. In vitro maturation of equine oocytes obtained from different age groups of sexually mature mares. Theriogenology. 1995;44(4):461–9.

    Article  CAS  PubMed  Google Scholar 

  12. Carnevale E, Uson M, Bozzola J, King S, Schmitt S, Gates H. Comparison of oocytes from young and old mares with light and electron microscopy. Theriogenology. 1999;51(1):299.

    Article  Google Scholar 

  13. Frank BL, Doddman CD, Stokes JE, Carnevale EM. Association of equine oocyte and cleavage stage embryo morphology with maternal age and pregnancy after intracytoplasmic sperm injection. Reprod Fertil Dev. 2019;31(12):1812.

    Article  CAS  PubMed  Google Scholar 

  14. Carnevale EM, Coutinho da Silva MA, Panzani D, Stokes JE, Squires EL. Factors affecting the success of oocyte transfer in a clinical program for subfertile mares. Theriogenology. 2005;64(3):519–27.

    Article  CAS  PubMed  Google Scholar 

  15. Altermatt JL, Suh TK, Stokes JE, Carnevale EM. Effects of age and equine follicle-stimulating hormone (eFSH) on collection and viability of equine oocytes assessed by morphology and developmental competency after intracytoplasmic sperm injection (ICSI). Reprod Fertil Dev. 2009;21(4):615.

    Article  CAS  PubMed  Google Scholar 

  16. Campos-Chillon F, Farmerie TA, Bouma GJ, Clay CM, Carnevale EM. Effects of aging on gene expression and mitochondrial DNA in the equine oocyte and follicle cells. Reprod Fertil Dev. 2015;27(6):925.

    Article  CAS  PubMed  Google Scholar 

  17. Cox L, Vanderwall DK, Parkinson KC, Sweat A, Isom SC. Expression profiles of select genes in cumulus–oocyte complexes from young and aged mares. Reprod Fertil Dev. 2015;27(6):914.

    Article  CAS  PubMed  Google Scholar 

  18. Rizzo M, Stout TAE, Cristarella S, Quartuccio M, Kops GJPL, De Ruijter-Villani M. Compromised MPS1 Activity Induces Multipolar Spindle Formation in Oocytes From Aged Mares: Establishing the Horse as a Natural Animal Model to Study Age-Induced Oocyte Meiotic Spindle Instability. Front Cell Dev Biol. 2021;9:657366.

    Article  PubMed  PubMed Central  Google Scholar 

  19. Rizzo M, du Preez N, Ducheyne KD, Deelen C, Beitsma MM, Stout TAE, et al. The horse as a natural model to study reproductive aging-induced aneuploidy and weakened centromeric cohesion in oocytes. Aging. 2020; [cité 9 nov 2020]; Disponible sur: http://www.aging-us.com/article/104159/text.

  20. Rizzo M, Ducheyne KD, Deelen C, Beitsma M, Cristarella S, Quartuccio M, et al. Advanced mare age impairs the ability of in vitro-matured oocytes to correctly align chromosomes on the metaphase plate. Equine Vet J. 2019;51(2):252–7.

    Article  CAS  PubMed  Google Scholar 

  21. Rambags BPB, van Boxtel DCJ, Tharasanit T, Lenstra JA, Colenbrander B, Stout TAE. Advancing maternal age predisposes to mitochondrial damage and loss during maturation of equine oocytes in vitro. Theriogenology. 2014;81(7):959–65.

    Article  CAS  PubMed  Google Scholar 

  22. Catandi GD, Obeidat YM, Broeckling CD, Chen TW, Chicco AJ, Carnevale EM. Equine maternal aging affects oocyte lipid content, metabolic function and developmental potential. Reproduction. 2021;161(4):399–409.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Hendriks WK, Colleoni S, Galli C, Paris DBBP, Colenbrander B, Roelen BAJ, et al. Maternal age and in vitro culture affect mitochondrial number and function in equine oocytes and embryos. Reprod Fertil Dev. 2015;27(6):957.

    Article  CAS  PubMed  Google Scholar 

  24. Woods GL, Baker CB, Hillman RB, Schlafer DH. Recent studies relating to early embryonic death in the mare. Equine Vet J. 1985;17(S3):104–7.

    Article  Google Scholar 

  25. Carnevale EM, Griffin PG, Ginther OJ. Age-associated subfertility before entry of embryos into the uterus in mares. Equine Vet J. 1993;25(S15):31–5.

    Article  Google Scholar 

  26. Brinsko SP, Ball BA, Miller PG, Thomas PGA, Ellington JE. In vitro development of day 2 embryos obtained from young, fertile mares and aged, subfertile mares. Reproduction. 1994;102(2):371–8.

    Article  CAS  Google Scholar 

  27. Willmann C, Schuler G, Hoffmann B, Parvizi N, Aurich C. Effects of age and altrenogest treatment on conceptus development and secretion of LH, progesterone and eCG in early-pregnant mares. Theriogenology. 2011;75(3):421–8.

    Article  CAS  PubMed  Google Scholar 

  28. Camargo Ferreira J, Linhares Boakari Y, Sousa Rocha N, Saules Ignácio F, Barbosa da Costa G, de Meira C. Luteal vascularity and embryo dynamics in mares during early gestation: effect of age and endometrial degeneration. Reprod Dom Anim. 2019;54(3):571–9.

    Article  CAS  Google Scholar 

  29. Morel MCGD, Newcombe JR, Swindlehurst JC. The effect of age on multiple ovulation rates, multiple pregnancy rates and embryonic vesicle diameter in the mare. Theriogenology. 2005;63(9):2482–93.

    Article  PubMed  Google Scholar 

  30. Vanderwall DK. Early embryonic development and evaluation of equine embryo viability. Vet Clin N Am Equine Pract. 1996;12(1):61–83.

    Article  CAS  Google Scholar 

  31. Panzani D, Rota A, Marmorini P, Vannozzi I, Camillo F. Retrospective study of factors affecting multiple ovulations, embryo recovery, quality, and diameter in a commercial equine embryo transfer program. Theriogenology. 2014;82(6):807–14.

    Article  PubMed  Google Scholar 

  32. Goszczynski DE, Tinetti PS, Choi YH, Hinrichs K, Ross PJ. Genome activation in equine in vitro–produced embryos. Biol Reprod. 2022;106(1):66–82.

    Article  CAS  PubMed  Google Scholar 

  33. Grondahl C, Hyttel P. Nucleologenesis and ribonucleic acid synthesis in preimplantation equine embryos. Biol Reprod. 1996;55(4):769–74.

    Article  CAS  PubMed  Google Scholar 

  34. Carnevale EM, Bergfelt DR, Ginther OJ. Aging effects on follicular activity and concentrations of FSH, LH, and progesterone in mares. Anim Reprod Sci. 1993;31(3–4):287–99.

    Article  CAS  Google Scholar 

  35. Chavatte-Palmer P, Velazquez MA, Jammes H, Duranthon V. Review: epigenetics, developmental programming and nutrition in herbivores. Animal. 2018;12(s2):s363–71.

    Article  CAS  PubMed  Google Scholar 

  36. McKinnon AO, Squires EL. Morphologic assessment of the equine embryo. J Am Vet Med Assoc. 1988;192(3):401–6.

    CAS  PubMed  Google Scholar 

  37. Iqbal K, Chitwood JL, Meyers-Brown GA, Roser JF, Ross PJ. RNA-Seq Transcriptome Profiling of Equine Inner Cell Mass and Trophectoderm. Biol Reprod. 2014;90(3) [cité 2 août 2019] Disponible sur: https://academic.oup.com/biolreprod/article-lookup/doi/10.1095/biolreprod.113.113928.

  38. Cao S, Wang JR, Ji S, Yang P, Chen J, Shen JP, et al. Differing total mRNA expression shapes the molecular and clinical phenotype of cancer. bioRxiv. 2020;57. https://doi.org/10.21203/rs.3.rs-600171/v1.

  39. Wang Z, Cao S, Morris JS, Ahn J, Liu R, Tyekucheva S, et al. Transcriptome Deconvolution of Heterogeneous Tumor Samples with Immune Infiltration. iScience. 2018;9:451–60.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. Savage SR, Shi Z, Liao Y, Zhang B. Graph algorithms for condensing and consolidating gene set Analysis results. Mol Cell Proteomics. 2019;18(8):S141–52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Shannon P. Cytoscape: A Software Environment for Integrated Models of Biomolecular Interaction Networks. Genome Res. 2003;13(11):2498–504.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Marinone AI, Losinno L, Fumuso E, Rodríguez EM, Redolatti C, Cantatore S, et al. The effect of mare’s age on multiple ovulation rate, embryo recovery, post-transfer pregnancy rate, and interovulatory interval in a commercial embryo transfer program in Argentina. Anim Reprod Sci. 2015;158:53–9.

    Article  CAS  PubMed  Google Scholar 

  43. Betteridge KJ, Eaglesome MD, Mitchellt D, Flood PF, Beriault R. Development of horse embryos up to twenty two days after ovulation: observations on fresh specimens. J Anat. 1982;135(Pt 1):191–209.

    CAS  PubMed  PubMed Central  Google Scholar 

  44. Beckelmann J, Budik S, Bartel C, Aurich C. Evaluation of Xist expression in preattachment equine embryos. Theriogenology. 2012;78(7):1429–36.

    Article  CAS  PubMed  Google Scholar 

  45. Colchen S, Battut I, Fiéni F, Tainturier D, Siliart B, Bruyas JF. Quantitative histological analysis of equine embryos at exactly 156 and 168 h after ovulation. J Reprod Fertil Suppl. 2000;56:527–37.

    Google Scholar 

  46. Leisinger CA, Medina V, Markle ML, Paccamonti DL, Pinto CRF. Morphological evaluation of day 8 embryos developed during induced aluteal cycles in the mare. Theriogenology. 2018;105:178–83.

    Article  CAS  PubMed  Google Scholar 

  47. Hisey E, Ross PJ, Meyers SA. A review of OCT4 functions and applications to equine embryos. J Equine Veterinary Sci. 2021;98:103364.

    Article  Google Scholar 

  48. Enders AC, Schlafke S, Lantz KC, Liu IKM. Endoderm cells of the equine yolk sac from day 7 until formation of the definitive yolk sac placenta. Equine Vet J. 1993;25(S15):3–9.

    Article  Google Scholar 

  49. Enders AC, Lantz KC, Liu IKM, Schlafke S. Loss of polar trophoblast during differentiation of the blastocyst of the horse. J Reprod Fertil. 1988;83(1):447–60.

    Article  CAS  PubMed  Google Scholar 

  50. Grøndahl ML, Andersen CY, Bogstad J, Nielsen FC, Meinertz H, Borup R. Gene expression profiles of single human mature oocytes in relation to age. Hum Reprod. 2010;25(4):957–68.

    Article  PubMed  CAS  Google Scholar 

  51. Reyes JM, Silva E, Chitwood JL, Schoolcraft WB, Krisher RL, Ross PJ. Differing molecular response of young and advanced maternal age human oocytes to IVM. Human Reprod. 2017;32(11):2199–208.

    Article  CAS  Google Scholar 

  52. Hamatani T, Falco G, Carter MG, Akutsu H, Stagg CA, Sharov AA, et al. Age-associated alteration of gene expression patterns in mouse oocytes. Human Mol Genet. 2004;13(19):2263–78.

    Article  CAS  Google Scholar 

  53. Kawai K, Harada T, Ishikawa T, Sugiyama R, Kawamura T, Yoshida A, et al. Parental age and gene expression profiles in individual human blastocysts. Sci Rep. 2018;8(1):2380.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  54. McCallie BR, Parks JC, Trahan GD, Jones KL, Coate BD, Griffin DK, et al. Compromised global embryonic transcriptome associated with advanced maternal age. J Assist Reprod Genet. 2019;36(5):915–24.

    Article  PubMed  PubMed Central  Google Scholar 

  55. Hauck CR, Hsia DA, Schlaepfer DD. The Focal Adhesion Kinase--A Regulator of Cell Migration and Invasion. IUBMB Life. 2002;53(2):115–9.

    Article  CAS  PubMed  Google Scholar 

  56. Provenzano PP, Keely PJ. Mechanical signaling through the cytoskeleton regulates cell proliferation by coordinated focal adhesion and Rho GTPase signaling. J Cell Sci. 2011;124(8):1195–205.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Zhao X, Guan J-L. Focal adhesion kinase and its signaling pathways in cell migration and angiogenesis. Adv Drug Delivery Rev. 2011;63(8):610–5.

    Article  CAS  Google Scholar 

  58. Hynes RO. Integrins: bidirectional, allosteric signaling machines. Cells. 2002;110(6):673–87.

    Article  CAS  Google Scholar 

  59. Cheuk IW-Y, Siu MT, Ho JC-W, Chen J, Shin VY. ITGAV targeting as a therapeutic approach for treatment of metastatic breast cancer. Am J Cancer Res. 2020;10(1):211–23.

    CAS  PubMed  PubMed Central  Google Scholar 

  60. Riento K, Ridley AJ. ROCKs: multifunctional kinases in cell behaviour. Nat Rev Mol Cell Biol. 2003;4(6):446–56.

    Article  CAS  PubMed  Google Scholar 

  61. Brouns MR. p190 RhoGAP is required for neural development. Development. 2000;127:4891–903.

    Article  CAS  PubMed  Google Scholar 

  62. Rane CK, Minden A. P21 activated kinases: structure, regulation, and functions. Small GTPases. 2014;5(1):e28003.

    Article  PubMed  PubMed Central  Google Scholar 

  63. Chan PM, Lim L, Manser E. PAK Is Regulated by PI3K, PIX, CDC42, and PP2Cα and Mediates Focal Adhesion Turnover in the Hyperosmotic Stress-induced p38 Pathway. J Biol Chem. 2008;283(36):24949–61.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  64. Bacon C, Endris V, Rappold GA. The cellular function of srGAP3 and its role in neuronal morphogenesis. Mechanisms Dev. 2013;130(6–8):391–5.

    Article  CAS  PubMed  Google Scholar 

  65. Wu MY, Hill CS. TGF-β superfamily signaling in embryonic development and homeostasis. Dev Cell. 2009;16(3):329–43.

    Article  CAS  PubMed  Google Scholar 

  66. Kemp C, Willems E, Abdo S, Lambiv L, Leyns L. Expression of all Wnt genes and their secreted antagonists during mouse blastocyst and postimplantation development. Dev Dyn. 2005;233(3):1064–75.

    Article  CAS  PubMed  Google Scholar 

  67. Smith CA, Hurley TM, McClive PJ, Sinclair AH. Restricted expression of DMRT3 in chicken and mouse embryos. Mech Dev. 2002;119:S73–6.

    Article  PubMed  Google Scholar 

  68. Zhang Y, Huang L, Fu H, Smith OK, Lin CM, Utani K, et al. A replicator-specific binding protein essential for site-specific initiation of DNA replication in mammalian cells. Nat Commun. 2016;7(1):11748.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  69. Yuan K, O’Farrell PH. Cyclin B3 is a mitotic cyclin that promotes the metaphase-anaphase transition. Curr Biology. 2015;25(6):811–6.

    Article  CAS  Google Scholar 

  70. Vinopal S, Černohorská M, Sulimenko V, Sulimenko T, Vosecká V, Flemr M, et al. γ-Tubulin 2 Nucleates Microtubules and Is Downregulated in Mouse Early Embryogenesis. PLoS ONE. 2012;7(1):e29919.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Fry AM, O’Regan L, Sabir SR, Bayliss R. Cell cycle regulation by the NEK family of protein kinases. J Cell Sci. 2012;125(19):4423–33.

    CAS  PubMed  PubMed Central  Google Scholar 

  72. Losada A, Hirano T. Shaping the metaphase chromosome: coordination of cohesion and condensation. Bioessays. 2001;23(10):924–35.

    Article  CAS  PubMed  Google Scholar 

  73. Schmidt TI, Kleylein-Sohn J, Westendorf J, Le Clech M, Lavoie SB, Stierhof Y-D, et al. Control of centriole length by CPAP and CP110. Curr Biol. 2009;19(12):1005–11.

    Article  CAS  PubMed  Google Scholar 

  74. Watanabe S, Shioi G, Furuta Y, Goshima G. Intra-spindle microtubule assembly regulates clustering of microtubule-organizing centers during early mouse development. Cell Reports. 2016;15(1):54–60.

    Article  CAS  PubMed  Google Scholar 

  75. Ladurner R, Kreidl E, Ivanov MP, Ekker H, Idarraga-Amado MH, Busslinger GA, et al. Sororin actively maintains sister chromatid cohesion. EMBO J. 2016;35(6):635–53.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  76. Rivera T, Losada A. Shugoshin and PP2A, shared duties at the centromere. Bioessays. 2006;28(8):775–9.

    Article  CAS  PubMed  Google Scholar 

  77. Welburn JPI. The molecular basis for kinesin functional specificity during mitosis. Cytoskeleton. 2013;70(9):476–93.

    Article  CAS  PubMed  Google Scholar 

  78. Zhou C-J, Wang X-Y, Han Z, Wang D-H, Ma Y-Z, Liang C-G. Loss of CENPF leads to developmental failure in mouse embryos. Cell Cycle. 2019;18(20):2784–99.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  79. Monzo P, Gauthier NC, Keslair F, Loubat A, Field CM, Marchand-Brustel YL, et al. Clues to CD2-associated protein involvement in cytokinesisϒD ϒV. Mol Biol Cell. 2005;16:12.

    Article  Google Scholar 

  80. Babst M, Davies BA. Regulation of Vps4 during MVB sorting and cytokinesis. Traffic. 2011;12:1298–305.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  81. Hales BF, Grenier L, Lalancette C, Robaire B. Epigenetic programming: from gametes to blastocyst. Birth Defects Res Part A Clin Mol Teratol. 2011;91(8):652–65.

    Article  CAS  Google Scholar 

  82. Marcho C, Cui W, Mager J. Epigenetic dynamics during preimplantation development. Reproduction. 2015;150(3):R109–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  83. Suo L, Meng Q-G, Pei Y, Yan C-L, Fu X-W, Bunch TD, et al. Changes in acetylation on lysine 12 of histone H4 (acH4K12) of murine oocytes during maternal aging may affect fertilization and subsequent embryo development. Fertil Steril. 2010;93(3):945–51.

    Article  CAS  PubMed  Google Scholar 

  84. Dumollard R, Duchen M, Carroll J. The role of mitochondrial function in the oocyte and embryo. In: The mitochondrion in the germline and early development: Elsevier; 2007. p. 21–49. [cité 18 sept 2020] (Current topics in developmental biology; vol. 77). Disponible sur: https://linkinghub.elsevier.com/retrieve/pii/S0070215306770028.

    Chapter  Google Scholar 

  85. Balaban RS, Nemoto S, Finkel T. Mitochondria, oxidants, and aging. Cell. 2005;120(4):483–95.

    Article  CAS  PubMed  Google Scholar 

  86. Wilding M, Dale B, Marino M, di Matteo L, Alviggi C, Pisaturo ML, et al. Mitochondrial aggregation patterns and activity in human oocytes and preimplantation embryos. Human Reproduction. 2001;16(5):909–17.

    Article  CAS  PubMed  Google Scholar 

  87. Gardner DK, Harvey AJ. Blastocyst metabolism. Reprod Fertil Dev. 2015;27(4):638.

    Article  CAS  PubMed  Google Scholar 

  88. Lane M, O’Donovan MK, Squires EL, Seidel GE, Gardner DK. Assessment of metabolism of equine morulae and blastocysts. Mol Reprod Dev. 2001;59(1):33–7.

    Article  CAS  PubMed  Google Scholar 

  89. Hutchison CA, Newbold JE, Potter SS, Edgell MH. Maternal inheritance of mammalian mitochondrial DNA. Nature. 1974;251(5475):536–8.

    Article  CAS  PubMed  Google Scholar 

  90. Krisher RL. Maternal age affects oocyte developmental potential at both ends of the age spectrum. Reprod Fertil Dev. 2019;31(1):1.

    Article  Google Scholar 

  91. Hendriks WK, Colleoni S, Galli C, Paris DBBP, Colenbrander B, Stout TAE. Mitochondrial DNA replication is initiated at blastocyst formation in equine embryos. Reprod Fertil Dev. 2019;31(3):570.

    Article  CAS  PubMed  Google Scholar 

  92. Waelchli RO, MacPhee DJ, Kidder GM, Betteridge KJ. Evidence for the presence of sodium- and potassium-dependent adenosine Triphosphatase alpha 1 and beta 1 subunit isoforms and their probable role in blastocyst expansion in the Preattachment horse conceptus. Biol Reprod. 1997;57:630–40.

    Article  CAS  PubMed  Google Scholar 

  93. do Nascimento AD, JCC M, ARR C, Batista AM, Kastelic JP, Câmara DR. Inhibition of Na+, K+ −ATPase with ouabain is detrimental to equine blastocysts. Anim Reprod. 2020;17(1):e20190079.

    Article  PubMed  Google Scholar 

  94. Arnaud MG, Dubroeucq MH, Rivot MD. Notation de l’état corporel des chevaux de selle et de sport: guide pratique. 1ère ed. Paris: Institut de l’élevage; 1997. p. 40.

    Google Scholar 

  95. Duchamp G, Bour B, Combarnous Y, Palmer E. Alternative solutions to hCG induction of ovulation in the mare. J Reprod Fertil Suppl. 1987;35:221–8.

    CAS  PubMed  Google Scholar 

  96. Bucca S, Carli A. Efficacy of human chorionic gonadotropin to induce ovulation in the mare, when associated with a single dose of dexamethasone administered at breeding time: efficacy of human chorionic gonadotropin to induce ovulation when associated with dexamethasone. Equine Vet J. 2011;43:32–4.

    Article  Google Scholar 

  97. Martin M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet J. 2011;17(1):10–2.

    Article  Google Scholar 

  98. Gutiérrez-Adán A, Perez-Crespo M, Fernandez-Gonzalez R, Ramirez M, Moreira P, Pintado B, et al. Developmental consequences of sexual dimorphism during pre-implantation embryonic development. Reprod Domest Anim. 2006;41(s2):54–62.

    Article  PubMed  Google Scholar 

  99. Pertea M, Pertea GM, Antonescu CM, Chang T-C, Mendell JT, Salzberg SL. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat Biotechnol. 2015;33(3):290–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  100. Robinson JT, Thorvaldsdóttir H, Winckler W, Guttman M, Lander ES, Getz G, et al. Integrative genomics viewer. Nat Biotechnol. 2011;29(1):24–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  101. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29(1):15–21.

    Article  CAS  PubMed  Google Scholar 

  102. Liao Y, Smyth GK, Shi W, et al. Bioinformatics. 2014;30(7):923–30.

    Article  CAS  PubMed  Google Scholar 

  103. R Core Team. A Language and Environment for Statistical Computing. Vienna: R foundation for Statistical Computing; 2020. Disponible sur: https://www.R-project.org/

    Google Scholar 

  104. Rstudio Team. RStudio: integrated development for R. RStudio. Boston: PBC; 2020. Disponible sur: http://www.rstudio.com/

    Google Scholar 

  105. Pinheiro J, Bates D, Debroy S, Sarkar D, R Core Team. nlme: Linear and Nonlinear Mixed Effects Models. 2020. Disponible sur: https://CRAN.R-project.org/package=nlme

    Google Scholar 

  106. Giraudoux P. pgirmess: Spatial Analysis and Data Mining for Field Ecologists. 2018. Disponible sur: https://CRAN.R-project.org/package=pgirmess

    Google Scholar 

  107. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  108. Lê S, Josse J, Husson F. FactoMineR: An R package for multivariate analysis. J Stat Soft. 2008;25(1):18.

    Article  Google Scholar 

  109. Gómez-Rubio V. ggplot2 - Elegant Graphics for Data Analysis (2nd Edition). J Stat Soft. 2017;77 [cité 9 août 2020](Book Review 2). Disponible sur: http://www.jstatsoft.org/v77/b02/.

  110. Kassambara A, Mundt F. factoextra: Extract and Visualize the Results of Multivariate Data Analyses. 2020. Disponible sur: https://CRAN.R-project.org/package=factoextra

    Google Scholar 

  111. Oksanen J, Blanchet FG, Friendly M, Kindt R, Legendre P, McGlinn D, et al. vegan: Community Ecology Package. 2020. Disponible sur: https://CRAN.R-project.org/package=vegan

    Google Scholar 

  112. Kolberg L, Raudvere U. gprofiler2: Interface to the “ g:Profiler ” Toolset. 2020. Disponible sur: https://CRAN.R-project.org/package=gprofiler2

    Google Scholar 

  113. Yates AD, Achutan P, Akanni W, Allen J, Allen J, Alvarez-Jarreta J, et al. Ensembl 2020.

    Google Scholar 

  114. Mi H, Ebert D, Muruganujan A, Mills C, Albou L-P, Mushayamaha T, et al. PANTHER version 16: a revised family classification, tree-based classification tool, enhancer regions and extensive API. Nucleic Acids Res. 2021;49(D1):D394–403.

    Article  CAS  PubMed  Google Scholar 

  115. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene set enrichment analysis: A knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci. 2005;102(43):15545–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  116. Mootha VK, Lindgren CM, Eriksson K-F, Subramanian A, Sihag S, Lehar J, et al. PGC-1α-responsive genes involved in oxidative phosphorylation are coordinately downregulated in human diabetes. Nat Genet. 2003;34(3):267–73.

    Article  CAS  PubMed  Google Scholar 

  117. Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P. The molecular signatures database Hallmark gene set collection. Cell Systems. 2015;1(6):417–25.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgments

The authors are grateful to the staff of the Institut Français du Cheval et de l’Equitation (IFCE) experimental farm (Domaine de la Valade, Chamberet, France) for care and management of animals. We acknowledge the High-throughput sequencing facility of I2BC for its sequencing and bioinformatics expertise. The bioinformatics analyses were performed thanks to Core Cluster of the Institut Français de Bioinformatique (IFB) (ANR-11-INBS-0013). Many thanks to Matthias Zytnicki and Christophe Klopp for their advice on RNA-seq de novo analysis. Many thanks to Pablo Ross who kindly provided the coordinates for the XIST gene. The authors thank Shavahn Loux for her information on PLAC8A. The authors are grateful to Alice Jouneau and Sophie Calderari who kindly revised this article.

Funding

This work was supported by the “Institut Français du Cheval et de l’Equitation” (grant numbers CS_2018_23, 2018). The National Research Institute for Agriculture, Food and Environment (INRAE) department Animal Physiology and Breeding Systems also supported this research.

Author information

Authors and Affiliations

  1. Université Paris-Saclay, UVSQ, INRAE, BREED, 78350, Jouy-en-Josas, France

    Emilie Derisoud, Luc Jouneau, Catherine Archilla, Nathalie Daniel, Nathalie Peynot, Michèle Dahirel, Véronique Duranthon & Pascale Chavatte-Palmer

  2. Ecole Nationale Vétérinaire d’Alfort, BREED, 94700, Maisons-Alfort, France

    Emilie Derisoud, Luc Jouneau, Catherine Archilla, Nathalie Daniel, Nathalie Peynot, Michèle Dahirel, Véronique Duranthon & Pascale Chavatte-Palmer

  3. IFCE, Plateau technique de Chamberet, 19370, Chamberet, France

    Cédric Dubois, Juliette Auclair-Ronzaud & Laurence Wimel

  4. Institute for Integrative Biology of the Cell (I2BC), UMR 9198 CNRS, CEA, Paris-Sud University F-91198, Gif-sur-Yvette, France

    Yan Jaszczyszyn

  5. Institut Pasteur—Bioinformatics and Biostatistics Hub—Department of Computational Biology, Paris, France

    Rachel Legendre

Authors
  1. Emilie Derisoud
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Luc Jouneau
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Cédric Dubois
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Catherine Archilla
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Yan Jaszczyszyn
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Rachel Legendre
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Nathalie Daniel
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Nathalie Peynot
    View author publications

    You can also search for this author in PubMed Google Scholar

  9. Michèle Dahirel
    View author publications

    You can also search for this author in PubMed Google Scholar

  10. Juliette Auclair-Ronzaud
    View author publications

    You can also search for this author in PubMed Google Scholar

  11. Laurence Wimel
    View author publications

    You can also search for this author in PubMed Google Scholar

  12. Véronique Duranthon
    View author publications

    You can also search for this author in PubMed Google Scholar

  13. Pascale Chavatte-Palmer
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

PCP obtained the funding. PCP and VD conceived the project. LW, PCP, and VD supervised the study. ED, CD, CA, ND, NP, LW, VD and PCP adapted the methodology for the project. ED, CD, CA, YJ, JAR, and LW performed the experiments. CD, CA, ND, NP, MD and LW provided the resources. ED, LJ, YJ and RL performed data curation. ED and LJ analysed the data. ED wrote the original draft. All authors read, revised, and approved the submitted manuscript.

Corresponding authors

Correspondence to Emilie Derisoud or Pascale Chavatte-Palmer.

Ethics declarations

Ethics approval and consent to participate

The experiment was performed at the experimental farm of IFCE (research agreement C1903602 valid until March 22, 2023). The protocol was approved by the local animal care and use committee (“Comité des Utilisateurs de la Station Expérimentale de Chamberet”) and by the regional ethical committee (“Comité Régional d’Ethique pour l’Expérimentation Animale du Limousin”, approved under N° C2EA - 33 in the National Registry of French Ethical Committees for animal experimentation) under protocol number APAFIS#14963–2018050316037888 v2. All experiments were performed in accordance with the European Union Directive 2010/63EU.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Supplementary Figure 1.

Analysis of ICM enriched and TE part before gene expression deconvolution using DeMixT. On first page, the analysis of differential expressed genes is represented. The second page represent the results of the gene set enrichment analysis of ICMandTE gene expression. The third page presents the clustering of gene sets altered by maternal age from SUMER analysis of the GO BP and KEGG terms of mixed ICMandTE and TE_part.

Additional file 2: Supplementary Figure 2.

Pipeline of biostatistical analysis. Once extracted, total RNA qualification and quantification were obtained using a Bioanalyzer. For sequencing, 5 ng/sample of total RNA were used for paired end RNA sequencing (I2BC platform) using Illumina NextSeq 500 high technology. Sequences were trimmed using Cutadapt. To determine XIST sequence, a de novo assembly was made using a mapping with STAR and StringTie. The GTF observed was aligned in IGV and XIST sequence was found at the known position. XIST sequence was included in EquCab 3.0 and a new mapping was performed using STAR. Counting was performed using featureCounts. The differential analysis was performed including embryo sex and diameter using Deseq2. A functional enrichment analysis of DEGs was performed using PANTHER web software. In another time, counts of all genes were normalized using Deseq2 and a gene set analysis (GSEA) was performed using the Gene Ontology (GO) biological process (BP) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) databases with GSEA software form the Broad Institute. To reduce redundancy between terms, SUMER analysis was performed on enriched gene sets.

Additional file 3: Supplementary Table 1.

Sex, diameter, RNA yield and mare age group of each equine embryo at Day 8 post-ovulation. Embryo ID, mare age group, embryo sex and diameter and RNA yield for each embryo are represented in columns while each line represents one individual embryo. OM: old mares; YM: young mares.

Additional file 4: Supplementary Table 2.

Differential gene analysis using DeSeq2 in ICMandTE of equine embryo at Day 8 post-ovulation. Equine ensemble ID, orthologue human Ensembl ID, Orthologue human Entrez Gene ID, gene description, normalized counts for each sample of ICM and parameters obtained after Deseq2 analysis (log2FoldChange, pvalue and padj (after FDR correction)) of genes expressed in ICM enriched part (before gene expression deconvolution using DeMixT) of OM and YM embryos. Young group was the reference group. ICM: Inner cell mass; OM: old mares; YM: young mares.

Additional file 5: Supplementary Table 3.

Gene set enrichment analysis results on gene expression comparing ICMandTE of embryos from old mares to young mares. Gene Set Enrichment Analysis results (pathway name, GO accession when possible and size, Normalized Enrichment Score, p-value and FDR corrected q-value) for GO biological process and KEGG databases on ICMandTE gene expression table. Young group was the reference group.

Additional file 6: Supplementary Table 4.

Differential gene analysis using DeSeq2 in DeMixT_ICM_cells of equine embryo at Day 8 post-ovulation. Equine ensemble ID, orthologue human Ensembl ID, Orthologue human Entrez Gene ID, gene description, normalized counts for each sample of ICM and parameters obtained after Deseq2 analysis (log2FoldChange, pvalue and padj (after FDR correction)) of genes expressed in ICM (after gene expression deconvolution of ICMandTE using DeMixT) of OM and YM embryos. Young group was the reference group. ICM: Inner cell mass; OM: old mares; YM: young mares.

Additional file 7: Supplementary Table 5.

Differential gene analysis using DeSeq2 in TE_part of equine embryo at Day 8 post-ovulation. Equine ensemble ID, orthologue human Ensembl ID, Orthologue human Entrez Gene ID, gene description, normalized counts for each sample of ICM and parameters obtained after Deseq2 analysis (log2FoldChange, pvalue and padj (after FDR correction)) of genes expressed in TE of OM and YM embryos. Young group was the reference group. TE: Trophoblast; OM: old mares; YM: young mares.

Additional file 8: Supplementary Table 6.

Gene set enrichment analysis results on gene expression of DeMixT_ICM_cells of embryos from young and old mares. Gene Set Enrichment Analysis results (pathway name, GO accession when possible and size, Normalized Enrichment Score, p-value and FDR corrected q-value) for GO biological process and KEGG databases on DeMixT_ICM_cells gene expression table (after gene expression deconvolution on ICMandTE using DeMixT). Young group was the reference group.

Additional file 9: Supplementary Table 7.

Gene set enrichment analysis results on gene expression of TE_part of embryos from young and old mares. Gene Set Enrichment Analysis results (pathway name, GO accession when possible and size, Normalized Enrichment Score, p-value and FDR corrected q-value) for GO biological process and KEGG databases on TE_part gene expression table. Young group was the reference group.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Derisoud, E., Jouneau, L., Dubois, C. et al. Maternal age affects equine day 8 embryo gene expression both in trophoblast and inner cell mass. BMC Genomics 23, 443 (2022). https://doi.org/10.1186/s12864-022-08593-7

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12864-022-08593-7

Keywords

BMC Genomics

ISSN: 1471-2164

Contact us