Kinetics of the chromosome 14 microRNA cluster ortholog and its potential role during placental development in the pregnant mare

Background The human chromosome 14 microRNA cluster (C14MC) is a conserved microRNA (miRNA) cluster across eutherian mammals, reported to play an important role in placental development. However, the expression kinetics and function of this cluster in the mammalian placenta are poorly understood. Here, we evaluated the expression kinetics of the equine C24MC, ortholog to the human C14MC, in the chorioallantoic membrane during the course of gestation. Results We demonstrated that C24MC-associated miRNAs presented a higher expression level during early stages of pregnancy, followed by a decline later in gestation. Evaluation of one member of C24MC (miR-409-3p) by in situ hybridization demonstrated that its cellular localization predominantly involved the chorion and allantoic epithelium and vascular endothelium. Additionally, expression of predicted target transcripts for C24MC-associated miRNAs was evaluated by RNA sequencing. Expression analysis of a subset of predicted mRNA targets showed a negative correlation with C24MC-associated miRNAs expression levels during gestation, suggesting the reciprocal control of these target transcripts by this miRNA cluster. Predicted functional analysis of these target mRNAs indicated enrichment of biological pathways related to embryonic development, endothelial cell migration and angiogenesis. Expression patterns of selected target mRNAs involved in angiogenesis were confirmed by RT-qPCR. Conclusion This is the first report evaluating C24MC kinetics during pregnancy. The findings presented herein suggest that the C24MC may modulate angiogenic transcriptional profiles during placental development in the horse. Electronic supplementary material The online version of this article (10.1186/s12864-018-5341-2) contains supplementary material, which is available to authorized users.


Background
Small non-coding RNAs (ncRNAs) constitute a group of RNAs which do not encode for proteins, but instead modify the expression of protein-coding mRNAs [1]. In plants and animals, microRNAs (miRNA) were one of the first identified families of regulatory ncRNAs [2,3]. MiRNAs are single-stranded small RNAs, approximately 22 nucleotides in length, that post-transcriptionally modulate protein-coding genes [4]. These ncRNAs regulate numerous cellular processes such as metabolism, cell proliferation, apoptosis, and cell differentiation in almost all cell types. The full or partial complementarity of base pairs between miRNA and their respective target mRNA lead to degradation or blocking the translation of the target mRNAs [4][5][6].
One of the unique features of miRNAs is that they can be transcribed as clusters. A miRNA cluster is a group of two or more miRNAs which are transcribed from physically adjacent miRNA genes in the same orientation [7,8]. These miRNA clusters are exclusively or preferentially expressed in a tissue-specific manner [9]. Recent reports have described several miRNA clusters expressed in human placenta, including the chromosome 14 miRNA cluster (C14MC), the chromosome 19 miRNA cluster (C19MC) and the miR-371-3 cluster also located in chromosome 19 [9][10][11][12]. The C19MC is exclusively found in primates while C14MC appears to be conserved among all eutherian species studied so far [12]. To date, the expression of miRNA clusters in equine placenta has not been reported.
The chromosome 14 miRNA cluster is one of the largest miRNA clusters transcribed in human and mouse (C12MC) placenta and consists of over 50 miRNAs. This cluster was originally described as two families located within two closely neighboring segments spanning about 40 kb [12][13][14]. The expression of C14MC has a tissue-specific pattern in mice and in humans [9,14]. In mice, C12MC is mostly expressed in the head and trunk of the developing embryo and it is restricted to the brain and placental tissues in the adult mice [14]. In healthy humans, expression of the C14MC is predominantly limited to the placenta and other tissues of epithelial origin [9]. The expression of C14MC was also observed in abnormal tissues in humans and mice, as potential tumor suppressor factors in gastric and liver cancers, mainly by inhibiting cell migration and invasion [15][16][17]. Currently, the tissue-specific expression of this cluster has not been characterized in species other than humans and mice.
While eutherian mammals have a wide diversity in the morphology of their placenta, the orthologous miRNA clusters are conserved without significant structural changes in all sequenced genomes and appear to be maintained in evolution for approximately 100 million years [13]. It has been suggested that the emergence of the C14MC was one of the factors that facilitated evolution of the mammalian placenta and may play important biological roles in placental and embryonic development [13,18]. For instance, the importance of the murine C12MC, orthologous to C14MC, has been experimentally verified in a mouse model by deletion of this cluster, which led to neonatal lethality [19]. The C12MC-associated miRNA knockout mice showed placental over growth and a defective feto-maternal interface [20,21]. In contrast, overexpression of C12MC-associated miRNAs was associated with growth retardation and postnatal lethality [22]. Although the C14MC has been reported to play an important role in placental development, its expression kinetics and its functions are poorly understood. Morales-Prieto et al. described the expression of miRNAs within C14MC in human trophoblast cells and showed a significant decrease in the expression levels of this cluster between first to third trimester of pregnancy, suggesting a strong decline in the expression of C14MC with advancing gestational age in humans [11]. Gu et al. also demonstrated that expression of miR-NAs within C14MC was significantly upregulated in the first trimester human placenta in comparison to postpartum placental samples [23]. These studies investigated the expression of C14MC at only two time-points due to limited sample availability, so there is minimal information on the kinetics of this cluster throughout pregnancy. Characterizing the expression pattern of miRNAs during physiological and pathological conditions is an important first step to investigate the function of this critical miRNA cluster.
We hypothesized that equine C24MC-associated miRNAs, the ortholog of human C14MC, have a differential expression pattern in chorioallantoic membrane (CAM) during different stages of gestation. In the present study, we determined the expression pattern of equine C24MC-associated miRNAs in CAM during pregnancy. We further analyzed the predicted functions of this cluster by generating a corresponding RNA-sequencing dataset from corresponding CAM. Additionally, we analyzed the cellular localization of one member of this cluster (eca-miR-409-3p) throughout gestation in equine CAM by in situ hybridization (ISH). This study identified that the expression of C24MC-associated miRNAs is downregulated during the course of gestation and determined the cellular localization of one of its members. Furthermore, we showed the presence of a negative correlation between the expression profile of these miRNAs and that of their putative mRNA targets, predicted to be functionally involved in embryonic development, endothelial cell migration and angiogenesis. Further studies to characterize the functions of each member of this miRNA cluster during placental development are warranted.

Animal use
The horses (Equus caballus) used in this study were mixed-breed ranging between 250 and 550 kg and four to sixteen years of age. Mares were born and raised at University of Kentucky research farm and housed on pasture with ad libitum grass hay. Gestational age was determined based on the day of ovulation (day 0).

CAM collection and preparation
Samples of CAM (n = 29) were collected from pregnant mares at 45 days (45d, n = 9), four months (4mo, n = 7), six months (6mo, n = 4) and ten months of gestation (10mo, n = 6), as well as from the CAM after normal parturition (n = 3, postpartum). For collection of CAM samples (except for CAM at 45d and postpartum), the uterus of pregnant mares was recovered immediately after euthanasia (using a barbiturate overdose following the American Veterinary Medical Association [AVMA] guidelines for the euthanasia of animals), the CAM was carefully separated from the endometrium, and full thickness CAM were collected from the body of the placenta. The CAM samples at 45d were retrieved from conceptuses collected by uterine lavage. The postpartum samples were collected immediately after normal parturition.
Collected CAM samples were stored overnight at 4°C in RNAlater™ (Life Technologies, Carlsbad, CA) and then stored at − 80°C until further processed. The remaining tissue samples were fixed in 10% formaldehyde for 24 h, transferred to 70% methanol and paraffin embedded following standard histological procedures [24]. Histological sections were stained with hematoxylin and eosin following routine procedures and assessed to confirm normal CAM without any contamination from endometrium. The study overview is presented in Fig. 1.

MiRNA extraction, cDNA synthesis and qPCR
CAM samples were thawed on ice and total RNA was extracted with Trizol™ (Life Technologies, Carlsbad, CA) following the manufacturer's recommendations. Total RNA concentration was measured using a Nano-Drop DP-1000 spectrophotometer (ThermoFisher Scientific, Waltham, MA).
Complementary DNA synthesis was carried out using the miScript II RT Kit (Qiagen, Hilden, Germany). The reaction carried out in 20 μL reverse transcription reaction mixture consisted of 2 μL of extracted RNA (containing 400 ng of RNA), 4 μL of miScript HiSpec buffer, 2 μL of miScript Nucleics mix, 10 μL of RNase-free water and 2 μL of miScript Reverse Transcriptase mix. Three μL of cDNA product from each sample were combined to make the pooled CAM sample. Pooled samples were only used as RT-qPCR inter-plate controls as previously described [25].
The expression of mature miRNAs was determined by RT-qPCR using the miScript SYBR Green PCR kit (Qiagen, Hilden, Germany) with the miScript Universal Primer along with miRNA-specific primers according to the manufacturer's guidelines. The reaction consisted of 1 ng (1 μl) of cDNA added to a reaction volume (9 μl) containing 2X QuantiTect SYBR® Green PCR Master Mix (5 μl), 10X miScript Universal Primer (1 μl), assay specific primers (1 μl of a 5 μM stock, final concentration 0.5 μM), and RNase free water (2 μl). The cycling conditions included an initial PCR activation step (95°C for 15 min), followed by 40 cycles of denaturation at 94°C for 15 s, annealing at 55°C for 30 s, and extension at 70°C for 30 s. A DNA melting curve was generated to discriminate between specific and non-specific amplification products. Primers were designed using miRprimer software (version 2.0) for 30 candidate miRNAs [26]. The candidate miRNAs were selected in three steps: (i) Based on prediction of C24MC in miRBase.org (miRbase.v.21) (EquCab2.0), 46 distinct miRNAs are included in this cluster with an inter-miRNA distance of less than 10 kb (EquCab2.0). The comparative list of C14MC in domestic animals is presented in the Additional file 1: Table S1. (ii) Based on our previously generated Illumina miRNAsequencing dataset of equine CAM, we narrowed the predicted miRNA list (Gene Expression Omnibus [GEO, NCBI, NIH] Accession GSE113142). Presence of all the predicted C24MC-associated miRNAs (based on miRBase.org) were confirmed in the CAM with the miRNA-sequencing dataset from 4mo (n = 3), 10mo (n = 3) and postpartum (n = 3). Information about miRNA-sequencing, mapping, normalization and statistical analysis are provided in the supplementary section (Additional file 2 and Additional file 3: Table S3 and Additional file 4: Table  S4). Primers were designed for the 11 miRNAs that showed significant differences between the examined time points (P < 0.05) and for the 17 most abundant miRNAs (primers for a total of 28 miRNA targets were designed). (iii)Based on a literature search, two additional miRNAs (hsa-miR-154 and hsa-miR-1247) were added to the list of candidate miRNAs [11,12,23]. Hsa-miR-154 belongs to the human C14MC and the equine orthologous miRNA has not been reported in miRBase database even though it was expressed in our miRNA-sequencing dataset. Similarly, hsa-miR-1247 was added to the candidate miRNAs since it has been jointly studied with the C14MC even though it is expressed from the opposite DNA strand [14]. It has been shown that the expression of these miRNAs change significantly during both normal and abnormal pregnancy conditions in human [27,28].
The list of miRNA candidates and primer sequences are provided in Table 1. Primer efficiency was verified on the pooled samples. Primers with C T values < 35 that did not produce primer-dimers were used for further experimentation. Otherwise, primers were re-designed and re-tested. Eca-miR-106a, eca-miR-8908a and eca-miR-369-5p were used as reference genes [25]. Real-Time qPCR was performed in triplicate for all samples, and geometric mean for the triplicate was calculated and used for further analysis [25,[29][30][31]. PCR efficiencies were calculated using LinRegPCR (version 2012.0) (http://www.hartfaalcentrum.nl) to ensure that all primers resulted in PCR efficiencies between 1.8 and 2.1.
Total RNA isolation, library preparation and RNA-sequencing (RNA-seq) A subset of CAM samples from 45d, 4mo, 6mo, and 10mo (n = 4, at each time point) were thawed on ice and total RNA was extracted using the RNeasy Mini Kit (Qiagen) per manufacturer's instructions. After extraction, RNA was analyzed by NanoDrop® (Thermo Fisher Scientific) and Bioanalyzer® (Agilent, Santa Clara, CA, USA) to evaluate concentration, purity and integrity of recovered RNA. All samples had a 230/260 ratio > 1.8, a 260/280 ratio > 2.0 and an RNA integrity number > 8.0 (mean ± STD = 9.2 ± 0.4).
Library preparation was performed using the TruSeq Stranded mRNA Sample Prep Kit (Illumina, San Diego, CA, USA) per manufacturer's instructions. The adapter for Read 1 was AGATCGGAAGAGCACACGTCTGAA CTCCAGTCACNNNNNNATCTCGTATGCCGTCTTC TGCTTG, with NNNNNN being the index sequence. The read 2 adapter was AGATCGGAAGAGCGTCGTG TAGGGAAAGAGTGTAGATCTCGGTGGTCGCCGTA TCAT. The libraries were quantified by qPCR. Sequencing was performed on a HiSeq 4000 (Illumina) using a HiSeq 4000 sequencing kit version 1, generating 150 bp paired-end reads. FASTQ files were generated and demultiplexed using bcl2fastq v2.17.1.14 Conversion Software (Illumina).

MiRNA target prediction
To assess pathways that might be affected by the C24MC-associated miRNAs, predicted targets of differentially expressed miRNAs were assessed using Ingenuity Pathway Analysis (IPA; Qiagen, Redwood City, CA, USA; Release December 2017). The miRNA and mRNA pairing was performed using the miRNA expression values (evaluated by RT-qPCR) and normalized read counts (RNA-seq dataset) using IPA's microRNA Target Filter, respectively. This allows predicting the function for the miRNAs according to the provided datasets (RNA-seq) (www.ingenuity.com/products/IPA/microRNA.html). Only target mRNAs that were previously (experimentally) confirmed or predicted with high confidence (IPA scoring system) were selected as putative targets. In addition to IPA predicted targets, putative targets were also identified using miRDB (http://mirdb.org), and mRNAs with a target prediction score (TPS) of ≥80 were added to the putative target list [35].
Since miRNAs mainly down regulate their target mRNA expression [36,37], those mRNAs in the RNA-seq data set which exhibited the reverse pattern of expression when compared to the miRNA expression pattern were selected for further analysis. The negative correlation (r ≤ − 0.5, P < 0.05) between the expression pattern of putative target mRNAs and that of their corresponding miRNA was analyzed during gestation (45d, 4mo, 6mo, 10mo; significance set at P-value < 0.05) and those target mRNAs presenting a significant negative correlation with the expression pattern of their corresponding miRNA were selected for canonical pathway and functional analysis [38]. The canonical pathways and functional analysis for the target mRNAs were performed using IPA and biological functions of the target mRNA were analyzed by Protein ANalysis THrough Evolutionary Relationships Classification System (PANTHER; Release 13.1) ontology classification system [39].

Analysis of putative miRNA:mRNA binding sites for selected angiogenic genes
For further analysis of potential miRNA:mRNA interactions, identification of the miRNA binding site at the 3' UTR of the selected equine angiogenic genes was performed. 3' UTR sequences of the genes which were shown to be involved in angiogenesis in our dataset were retrieved from NCBI and the binding site was identified using IntaRNA 2.0 (Freiburg RNA Tools, Universität Freiburg, Germany) [42,43]. The free energy (ΔG) of the predicted miRNA:mRNA interaction was estimated using the RNAhybrid tool (Universität Bielefeld, Germany) [44].

MiRNA localization by in situ hybridization
The expression pattern of eca-miR-409-3p was further investigated by chromogenic ISH using a dual digoxigenin (DIG)-labeled LNA™ probe specific to hsa-miR-409-3p (610701-360; miRCURY LNA™, Exiqon, Vedbaek, Denmark). In addition, a dual DIG-labeled LNA™ probe specific to U6 snRNA (#699002-360; Exiqon) and a scrambled miRNA probe (#699003-360; Exiqon) were used as positive and negative controls, respectively. The sample set evaluated included CAM from 45d (n = 2), 4 mo (n = 2), 6 mo (n = 2), 10 mo (n = 2), and postpartum (n = 2). The ISH assay was performed as recommended by the manufacturer. Briefly, five-micrometer sections of formalin-fixed, paraffin-embedded tissues were mounted on positively charged Superfrost® Plus slides (Fisher Scientific, Pittsburgh, PA), dried overnight at 37°C, and subjected to deparaffinization in xylenes followed by rehydration through an ethanol dilution series and 1X phosphate buffered saline. Automated ISH was performed in a Ventana Discovery Ultra (Ventana Medical Systems, Inc. AZ, USA). The tissue sections were digested with 10 μg/ml of proteinase K for 16 min at room temperature and hybridized with the respective probes (diluted to 40 nM for hsa-miR-409-3p and scrambled miRNA probes, or 0.1 nM for U6 snRNA probe in DIS-COVERY ISH diluent (Ventana)) for two hours at 50°C. Stringency washes were performed with RiboWash (Ventana). The bound dual-labeled probe was detected with an alkaline phosphatase-conjugated polyclonal anti-DIG antibody (DISCOVERY anti-DIG AP Multimer, Ventana) after a 15 min blocking step using the DISCOVERY Antibody Block (Ventana), and the signal was detected using nitro blue tetrazolium chloride/5-bromo-4-chloro-3-indolyl-phosphate (NBT-BCIP; Chromo-Map Blue Kit, Ventana) as the substrate. The sections were incubated with the substrate for two hours and counterstained using the Red Counterstain II (Ventana) for 8 min. Sections were finally dehydrated in successive alcohols and mounted using permanent mounting medium (Eukitt®, Sigma-Aldrich, St. Louis, MO). ISH signal was evaluated and quantified by two different investigators blinded to gestational age with a Zeiss operative microscope equipped with MicroPublisher 5.0 RTV camera (Q-Imaging, Burnaby, BC, Canada). For quantification purposes, sections were scored as follows: 0, no staining; 1, weak or minimal staining; 2, mild staining; 3, moderate staining; 4, intense or abundant staining.

Statistical analysis
Delta C T (ΔC T ) values for CAM samples were calculated where ΔC T = (C T values of miRNA of interest -Geometric mean of the C T values of all three reference miRNAs) [31,45]. CAM results are presented as -ΔC T (negative ΔC T is more intuitive than ΔC T ). Graphs were made in JMP (SAS Institute, version 12.1.0) and "d3heatmap" or "plot" packages in R (version 1.0.136) [46].
Statistical analyses were performed in SAS version 9.2 (SAS Institute Inc., Cary, NC, USA, 2010). Normality of the data (-ΔC T , microRNA-read counts, and FPKM) was verified using "PROC UNIVARIATE". Normal quantile transformation was performed to normalize the data [47]. MiRNA expression levels (-ΔC T ) were compared between pregnancy stages using the mixed model by "PROC MIXED", where stage of pregnancy was set as the fixed effect and animal was set as the random effect [48]. Post hoc analysis was performed using Tukey's Test with significance set at P < 0.05. The correlation between the RT-qPCR results and miRNA-sequencing results was determined using Pearson correlation coefficient performed in the "Hmisc" package in R with significance set to P < 0.05 [49]. The correlation between the expression patterns of each miRNA and its putative target mRNAs was determined using Pearson correlation coefficient during gestation with significance set to P < 0.05 (the false discover rate [FDR] was calculated in the "p.adjust" function in R) [50]. Negative ΔC T was calculated for the predicted genes involved in angiogenesis where ΔC T = (C T values of predicted mRNA -Geometric mean of the C T values of all three housekeeping genes). The correlation between target gene's FPKM and -ΔC T was evaluated using Pearson correlation coefficient during gestation due to the normality of the data. Also, the correlation between miRNAs`-ΔC T and their target-genes`-ΔC T was measured using Pearson correlation coefficient during gestation.
Mean ISH scores between the two independent observers were determined, and mean scores were compared between the different stages of pregnancy using the logistic regression model for overall significance and Chi square test was used for inter-group (stages of pregnancy) comparison. Significance was set at P < 0.05. Due to the non-normal distribution of the ISH scores, the correlation between the expression levels of eca-miR-409-3p determined by RT-qPCR and ISH mean scores was calculated by Spearman correlation (ρ) with significance set at P < 0.05.

Results
Differential expression of C24MC-associated miRNAs in CAM during pregnancy The expression of 26 mature miRNAs belonging to the C24MC was quantified in normal equine CAM (n = 29) at five different stages of pregnancy (45d, 4mo, 6mo, 10mo, and postpartum) by RT-qPCR (Fig. 2), and the differentially expressed miRNAs (n = 23) during the different stages of pregnancy are presented in Table 2 (P < 0.05). A total of 21 miRNAs presented a significantly higher expression level at 45d followed by a downregulation later in gestation ( Fig. 2 and Table 2). In contrast, eca-miR-370-3p had a lower expression at 45d of gestation in comparison to CAM collected at later stages ( Fig. 2 and Table 2). Furthermore, eca-miR-379-5p, eca-miR-432-5p and eca-miR-487b-3p showed a higher expression in postpartum samples in comparison to CAM collected at 6mo. Moreover, eca-miR-3958-3p expression was significantly higher in postpartum compared to CAM from 6mo and 10mo of gestation (Table 2). A significant positive correlation was observed between the normalized read counts derived from miRNA sequencing and -ΔC T values for each of the samples (r ranging from 0.71 to 0.9, P < 0.002) and the overall correlation coefficient between normalized read counts and -ΔC T was 0.77 (P < 0.001). Four miRNAs out of 30 tested miR-NAs (eca-miR-136-3p, eca-miR-136-5p, eca-miR-494-3p, and eca-miR-495-3p) were not expressed in > 50% of the samples, with some C T values higher than the accepted range (C T > 35), therefore these miRNAs were excluded from the dataset.

C24MC-associated miRNA target prediction
Computational target prediction for miRNAs was performed to identify putative mRNA targets and the expression patterns of predicted target mRNAs present in CAM were evaluated in an Illumina RNA-sequencing dataset generated from a subset of samples (n = 16) as indicated under Materials and Methods, subsection 2.4. We generated a total of 723 million reads for 16 samples, and an average of 91.23% (range 89.83 -92.72%) of reads were uniquely mapped to the reference genome (information about individual sample read count and mapping quality can be found in Additional file 6: Table S6).
A total of 796 target mRNAs were predicted using IPA's microRNA Target Filter function for the differentially expressed miRNAs that presented a significantly higher expression level at 45d followed by a decline later in the gestation (n = 21; Table 3). Additionally, 1179 target genes were predicted for this subset of miRNAs using mirDB. However, 229/1179 predicted genes were not expressed in the associated CAM RNAseq dataset and were, therefore, excluded from further analysis (Table 3) (the list of putative target mRNAs for each miRNA is available in Additional file 7: Table S7). There was a low agreement between the target lists within the two programs, which is due to the different algorithms used by these programs. Since miRNAs mainly down regulate their target mRNA expression [36,37], those mRNAs in the RNA-seq data set which exhibited the reverse pattern of expression when compared to the miRNA expression pattern (negative correlation, r ≤ − 0.5, P < 0.05) were selected for further analysis (Table 4) (n = 130 mRNAs, among which seven mRNAs have been predicted for more than one miRNA). High expression levels of C24MC-associated miRNAs in equine CAM collected at 45d are associated with reduced expression levels of the predicted mRNA targets with an inverse relationship as gestational age increases ( Fig. 3a and b).
The physiological functions and canonical pathways for the target mRNAs were predicted using Ingenuity Pathway Analysis. The predicted physiological function analysis included involvement in cardiovascular development, embryonic development, tissue morphology and nervous system development and function (Table 5). Specific functions such as endothelial cell movement/ migration, angiogenesis, differentiation of angioblasts, formation of vascular plexus, vascularization and morphology of blood vessels were among the annotated cardiovascular development processes (Fig. 4a). The expression pattern of the mRNAs and miRNAs involved in angiogenesis and vascular formation during gestation are presented in Fig. 4b and c, respectively. Interestingly, 31% (8/26) of the miRNAs evaluated were putatively associated with this specific biological function. In addition, immunosuppression and hypersensitivity/immediate hypersensitivity reactions were among the immunological functions assigned to the predicted targets. The complete list of functions and diseases associated with this subset of miRNAs can be found in Additional file 8: Table S8. The canonical pathways predicted for the target mRNAs primarily involved signaling pathways, including axonal guidance signaling, protein kinase A signaling, cAMP-mediated signaling, insulin-like growth factor 1 (IGF-1) and insulin receptor signaling, IL-4 signaling, and MAPK signaling, among others. The complete list of associated canonical pathways can be found in Additional file 9: Table S9. Gene ontology analysis of predicted mRNA targets was performed using PANTHER, and the biological processes mainly involved cellular (tyrosine kinase signaling pathway, MAPK cascade, cell proliferation, cycle and growth), metabolic (hydrolase activity and protein kinase activity), developmental (anatomical structure morphogenesis, cell differentiation, embryo development, ectoderm and mesoderm

Expression analysis of putative mRNA targets involved in angiogenesis
Since the top physiological function predicted for targeted mRNAs constituted cardiovascular system development (Table 5) and due to the importance of vessel formation in placenta development and function, we specifically selected a subset of genes (n = 10) participating from vascular morphology and cell movement ( Fig. 4a and Additional file 8: Table S8) to confirm their expression by RT-qPCR during pregnancy in equine CAM (n = 29). These included ACVRL1, AXL, CXCL10, DVL1, EDNRA, EPN2, G6PD, JAM2, NRP1 and TYRO3. There was a significant positive correlation between the expression pattern evaluated by RT-qPCR and RNA-seq data (i.e. -ΔC T vs FPKM, 0.6 < r < 0.8 and P < 0.001). There was also a significant negative correlation between the expression pattern of miRNAs and their putative target mRNA ( Table 6). The lowest expression level for these angiogenic-related genes was observed at 45d with a gradual increase in their expression level toward the end of gestation (Fig. 4d). Interestingly, their respective miRNAs presented the highest expression at 45d, decreasing thereafter (Fig. 4c). In addition, we evaluated the putative binding sites between the miR-NAs and their putative mRNA targets (3' UTR) and determined that these occur with a low free energy and, therefore, these miRNA:mRNA interactions are feasible (ΔG between − 12.6 and − 32.2; Table 7).

Localization of eca-miR-409 by in situ hybridization
In order to study the cellular and subcellular localization of one member of C24MC (eca-miR-409-3p) in CAM, we performed miRNA in situ hybridization (ISH) on CAM sections during five different stages of pregnancy. Since eca-miR-409-3p had a differential expression pattern during different stages of pregnancy as determined by RT-qPCR analyses and presented a 100% homology Data are presented as Mean ± Standard Error of -ΔC T The expression pattern for each miRNA were compared between timepoints using the mixed model and significant differences between timepoints were determined by subsequent pairwise comparisons. Thus, differing superscripts indicate significant differences between stages of pregnancy with significance set to P ≤ 0.05 to hsa-miR-409-3p, this miRNA was chosen as the best candidate for ISH analysis. MiR-409-3p signal was significantly higher at 45d (mean score of 3.25, P = 0.002), and its expression showed a decrease towards the end of pregnancy (Fig. 6c-h). The lowest miR-409-3p signal was observed at 6mo and 10mo (mean scores of 1, Fig.  6e-f and h). There was a significant correlation between -ΔC T values and mean ISH scores (ρ = 0.9, P < 0.01). At 45d, miR-409-3p signal was predominantly observed in the cytoplasm of the chorionic and allantoic epithelium (Fig. 6c). MiR-409-3p was also abundantly expressed in mesenchymal cells and vascular endothelium (Fig. 6c). As indicated above, a significant reduction in miR-409-3p expression was observed at 4mo, 6mo and 10mo as well as in postpartum samples. At 4mo, positive signal was observed in scattered chorionic epithelial cells, which presented a mild intensity staining (Fig. 6d). Low signal was observed in the vascular endothelium, allantoic epithelium and scattered mesenchymal cells. At 6mo and 10mo, a low expression of miR-409-3p was observed in the chorionic epithelium, with basolateral expression in scattered epithelial cells specifically at 10mo (Fig. 6e, f ). Finally, low miR-409-3p expression was observed in the vascular endothelium of scattered blood vessels, while no signal was detected in mesenchymal cells after 6mo.

Discussion
In the present study, for the first time, we comprehensively and comparatively evaluated the kinetics of the equine C24MC (26 miRNAs), the ortholog of human C14MC, in the normal chorioallantoic membrane during the course of gestation. The chromosome 14 miRNA cluster, also referred to as the Mirg cluster [19], the miR-379/miR-410 cluster [51] or the miR-379/miR-656 cluster [13], is one of the largest identified miRNA clusters [14]. This cluster is located at the imprinted domain on the human 14q32 chromosomal interval (mouse chromosome 12 and equine chromosome 24) and abundantly expressed in the normal placenta [14]. By evaluating all sequenced genomes of mammals, it has been shown that this cluster is preserved without significant structural changes [13]. For example, all human C14MC miRNAs are conserved in the equids with the exception of five (hsa-miR-154, hsa-mir-300, mir-654, hsa-mir-665 and hsa-mir-668). Although human and murine C14MC   and C12MC do not encode miR-3958 and miR-3959, these miRNAs have been predicted in equine C24MC and ovine C18MC orthologous clusters. This is the first report to unequivocally confirm the expression of the predicted eca-miR-3958-3p and eca-miR-3958-5p in equine CAM tissues. Here, we determined a significantly higher expression level of these miRNAs at 45d followed by a decline with increasing gestational age, demonstrating that these miRNAs have a similar expression pattern as other members of C24MC. Since current search databases for the identification of mRNA targets are based solely on human and murine miRNA annotation, prediction of eca-miR-3958-5p specific mRNA targets could not be performed. Currently, their predicted mRNA targets remain unknown and, therefore, identification of these targets is required to gain further insight into their functional role during gestation.
In previous studies, it has been shown that C14MC in humans is expressed predominantly during the first trimester of pregnancy and that its expression decreases as pregnancy advances [11,23]. Morales-Prieto et al. analyzed the expression of C14MC in human cytotrophoblast cells collected from healthy first trimester placenta and term placenta (collected from postpartum samples) and showed that the expression of miRNAs within C14MC decreases significantly from first to third trimester [11]. This conclusion was based on a decline in expression levels of 16 miRNAs at only two time points. A recent study also observed the same expression pattern for six miRNAs within this cluster in villous tissue from human placenta obtained at the first trimester of pregnancy and postpartum placenta [23]. In the present study, we examined the expression kinetics of C24MC-associated miRNAs at five gestational time points in the entire CAM. Similar to human C14MC, we Table 4 Correlation analysis between C24MC-associated miRNAs presenting a significantly higher expression level during early stage of pregnancy (45d), followed by a decline later in the gestation, and that of their putative mRNA targets (n = 130) (Continued)  The miRNA expression was evaluated by RT-qPCR (− ΔC T ) and the respective mRNA expression was determined by RNA-sequencing (FPKM) Since miRNAs down regulate their target mRNA/s expression, those target mRNAs in the RNA-seq dataset that exhibited the reverse pattern of expression compared to the miRNA expression pattern (negative correlation, r ≤ −0.5, P < 0.05) were selected. There was no significant correlation between the expression pattern of eca-miR-487b-3p and that of its predicted targets determined that the majority of the evaluated miRNAs (21/26) were highly expressed at the beginning of pregnancy followed by a decline towards the end of pregnancy. Accordingly, we observed a similar pattern of expression for miR-409-3p, a C24MC member, by ISH and we determined the cellular localization of this miRNA within CAM during the course of gestation. It is important to notice that, although some miRNA sequences from this cluster overlap with specific genes in the opposite DNA strand, the high RIN values in the samples evaluated determine the reliability of the miRNA expression analysis performed herein. This is further supported by the fact that ongoing studies in our laboratory to elucidate the interplay between eca-miR-127 and RTL1, encoded in opposite DNA strands, demonstrated that the expression of these exhibit a negative correlation during pregnancy whereby expression of eca-miR-127 decreases throughout gestation while that of RTL1 increases (Dini and Ball, 2018, unpublished). Several studies have been published on the importance and functions of miRNAs on placental performance [52,53]. For example, by restricting genes encoding critical miRNA biosynthetic enzymes such as Dicer [54] or Ago2 [55,56], normal placental function was impaired, leading to early post-implantation embryonic lethality in mice. Others have used more targeted murine gene knockout experiments to elucidate miRNA functions in the placenta [57]. Even though these experimental functional analysis are available, they are limited to small number of miRNAs and target genes [35]. The findings presented herein suggest that the differential expression of several miRNAs within C24MC may be associated with placental development and performance during the course of gestation. In that regard, experimentally evaluating their functional role in vivo or in vitro constitutes a challenging task, particularly for a cluster of miRNAs. Thus, computational prediction of miRNA targets is the main tool for identifying possible miRNA:mRNA interactions and miRNA functions when a group of miRNAs are studied simultaneously [35,58]. Computational algorithms can generate a list of putative genes targeted for the miRNA of interest; however, the presence of these genes in the specific cell or tissue of interest is not guaranteed. Another drawback of computational target predictions involves the occurrence of false putative targets (false positives). It has been shown that most of the target prediction algorithms have false positive rates of several tens of percent and, thus, may compromise functional inferences [59][60][61]. In the current study, to assess the predicted targets and have an insight into the functional role of C24MC-associated miRNAs in CAM, we generated an RNA-seq dataset from the same sample set to characterize the specific transcriptome profile of equine CAM throughout pregnancy as an alternative approach to specifically identify our predicted miRNA-targeted genes in the tissue of interest (CAM) with a higher confidence. Initially, we confirmed the presence of predicted target mRNAs in CAM tissue, followed by extensive correlation analysis between the expression profile of each miRNA with that of its predicted target mRNAs throughout pregnancy. It is well-known that miRNAs downregulate their target mRNA expression and, therefore, the expression pattern of specific miRNAs and that of their targets tend to be negatively correlated [36,37]. We analyzed the negative correlation for the predicted targets and were able to narrow the list of putative targets to 130 genes (out of 1975 initially predicted). We believe this approach provided us with a more representative putative target list for C24MC in equine CAM.
Several functions have been proposed for this cluster including immune suppressive functions, participation in the anti-inflammatory response, protection from the ischemic/hypoxia injury, cell motility and participation in angiogenesis [10,11,51,62]. It has been shown that C14MC has a pivotal effect on the feto-maternal interface of the placenta specifically on development of the endothelial and trophoblast cell interaction [19,63,64]. Mice with knockout C14MC had severe abnormalities associated with fetal capillaries in the placenta [63]. This capillary abnormality leads to overgrowth of the placenta, due to the large number of vacuoles in endothelial cells and the surrounding trophoblast layer of the fetal capillaries in the labyrinth zone causing neonatal lethality [63]. Interestingly, it has been shown that disruptions in the expression of this cluster in tissues other than placenta are related to vascular invasion [65]. In line with this, our computational target prediction analysis suggested that miRNAs from the C24MC target specific mRNAs involved in angiogenesis, angiogenesis of tissue of epithelial origin, vascularization, and development and migration of endothelial cells. Also, in another knockout model in mice, it has been shown that murine C12MC plays a number of roles required for neonatal survival [19]. Additionally, in vitro studies have shown that several miRNAs regulate trophoblast cell proliferation, invasion, and angiogenesis [66]. Consequently, the changes in the expression of C24MC-associated miRNAs during the course of pregnancy suggest that this may lead to modifications in angiogenesis within the placenta as pregnancy advances. For instance, it has been shown that CXCL10, which is a target gene for miR-411, can modify the development of newly formed vasculature [67]. We hypothesize that the reduction in the expression of some of the miRNAs in this cluster triggers the upregulation of angiogenesis-related genes that allow the process of angiogenesis to occur as the placenta develops. Here, we investigated the expression pattern of 10 selected putative mRNAs including CXCL10 during the course of gestation and demonstrated a negative correlation with the respective miRNAs targeting this set of genes. For instance, the low expression of these mRNAs Table 5 The physiological function analysis prediction for the putative mRNA targets of C24MC-associated miRNAs with a significantly higher expression level during early stages of pregnancy (45d), followed by a decline toward the end of the gestation at 45d was accompanied by a high expression of their respective miRNAs. Furthermore, estimation of the miRNA binding sites at the 3' UTR of their respective (predicted) mRNA targets demonstrated that all of the putative miR-NA:mRNA interactions occur with a low free energy, which is indicative that these RNA:RNA interactions are feasible. Even though these results suggest the potential influence of these miRNAs in the regulation of this set of angiogenic genes, further in vivo or in vitro studies to elucidate the effect of these miRNAs on angiogenesis during placental development are warranted.
Finally, we localized eca-miR-409-3p in CAM and compared its expression during different stages of pregnancy. Our data suggest that the expression level and localization of this C24MC-associated miRNA is dependent upon the developmental stage of the CAM. While eca-miR-409 expression was high in the chorion and allantoic epithelia, mesenchymal cells and vascular endothelium in early pregnancy (45d), its expression was significantly reduced after 4mo with gradual loss of expression in these cell types. The downregulation of this miRNA in vascular endothelial cells along with an increase in the expression (See figure on previous page.) Fig. 4 a. Pathway analysis of negatively correlated target mRNAs (n = 130) for differentially expressed C24MC-associated miRNAs presenting a decline in their expression during the course of gestation. Target mRNAs were mainly associated with cardiovascular development and function, and included endothelial cell movement/migration, angiogenesis, differentiation of angioblasts, formation of vascular plexus, and vascularization and morphology of blood vessels. Specific C24MC-associated miRNAs and direct relationships are shown in blue. Relevant physiological functions/pathways involved and their direct relationships are depicted with varying sizes according to the associated number of target mRNAs (ie. larger size indicative of higher number of mRNAs involved). The interactive figure can be accessed here. (http://rpubs.com/pouyadini/431871) b. Expression patterns of mRNAs (n = 29 putative targets), which presented a significant negative correlation to their miRNA expression during gestation and are predicted to be involved in angiogenesis. The heatmap was constructed using FPKM values and scaled for each mRNA. The interactive heatmap can be accessed here. (http://rpubs.com/pouyadini/339221) c. Expression patterns of miRNAs (n = 8), which were predicted to be involved in angiogenesis in CAM during different stages of pregnancy (n = 29 samples). The heatmap was constructed using -ΔC T values and scaled for each miRNA. The interactive heatmap can be accessed here. of angiogenic genes further suggests that this miRNA cluster is likely involved in the regulation of angiogenesis during placental development. Moreover, there was a significant correlation between expression levels of this specific miRNA determined by RT-qPCR and ISH. This is the first report on the expression of one of the C24MC members at the cellular and subcellular level in CAM throughout the pregnancy. The differential expression pattern of this miRNA suggests its role during pregnancy is dynamic and dependent on the stage of gestation.

Conclusion
In conclusion, for the first time, we have comprehensively evaluated the kinetics of the C24MC in CAM, orthologous to the human C14MC, during equine pregnancy. We demonstrated that the expression of C24MC-asssociated miRNAs decreased during the course of gestation. We generated an RNA-seq dataset from equine CAM during pregnancy and predicted the functions of C24MC-associated miRNAs in this tissue.
We have correlated the expression profile of C24MC-associated miRNAs and targeted mRNAs to evaluate the role of this cluster in the modulation of gene expression during placental development in vivo. We observed an  Binding site and free energy estimation not provided due to the lack of annotation of the 3′ UTR in the current equine genome assembly h The bar graph shows the mean value for each ISH score with standard error bars in CAM during the course of gestation. Mean miR-409 ISH scores between the two independent observers were determined, and mean scores were compared between the different stages of pregnancy