Genome-wide microRNA expression analysis in human placenta reveals sex-specific patterns: an ENVIRONAGE birth cohort study.

There is an increasing interest in microRNAs (miRNAs) as they are of utmost importance in gene regulation at the posttranscriptional level. Sex-related susceptibility for non-communicable diseases later in life could originate in early life. Until now, no data on sex-specific miRNA expression are available for the placenta. Therefore, we investigated the difference by sex of newborn's miRNA expression in human placental tissue. Within the ENVIRONAGE birth cohort, miRNA and mRNA expression profiling was performed in 60 placentae (50% boys) using Agilent (8 × 60 K) microarrays. The distribution of chromosome locations was studied and pathway analysis of the identified sex-specific miRNAs in the placenta was carried out. Of the total 2558 miRNAs on the array, 597 miRNAs were expressed in over 70% of the samples and were included for further analyses. A total of 142 miRNAs were significantly (FDR<0.05) associated with the newborn's sex. In newborn girls, 76 miRNAs had higher expression (hsa-miR-361-5p as most significant) and 66 miRNAs had lower expression (hsa-miR-4646-5p as most significant) than in newborn boys. In the same study population, placental differentially expressed genes by sex were also identified using a whole genome approach. The placental gene expression revealed 27 differentially expressed genes by comparing girls to boys. Ultimately, we studied the miRNA-RNA interactome and identified 14 miRNA-mRNA interactions as sex-specific. Sex differences in placental m(i)RNA expression may reveal sex-specific patterns already present during pregnancy, which may influence physiological conditions in early or later life. These molecular processes might play a role in sex-specific disease susceptibility in later life.

Introduction

Until now, a total of 2588 mature miRNAs have been described in humans [1]. Many of these miRNAs have been suggested as biomarkers for human diseases [2] and can be involved in a wide range of developmental and physiological processes during intrauterine life, suggesting a persistent impact throughout the prenatal and postnatal periods [3].

Although men and women share most of their genome, they differ in the presence of the Y chromosome in men, in which a relatively small number of genes are located [4,5]. Furthermore, they can differ in the transcriptional regulation of autosomal genes [6]. Sexual dimorphism has been described in cardiovascular diseases, immune-related diseases, and asthma, and has been attributed to differential gene regulation of sex steroid-responsive genes [7]. These sexual variations in the expression levels of sex hormones already begin in utero and persist throughout life. Interestingly, human miRNA expression has been shown to be different in each sex in both development and disease state [8,9]. Sex-specific miRNA expression has been observed in cord blood of healthy newborns [10], in the developing neonatal brain upon prenatal maternal stress, and is linked to neurodevelopmental disorders later in life [11]. Additionally, sex-related differences in miRNA expression have been reported in several diseases [12,13], including metabolic syndrome [14], cerebral ischaemia [15], hepatocellular [16], and squamous cell carcinomas [17]. Sex-specificity in placental miRNA expression has previously been demonstrated in association with maternal pre-pregnancy BMI [18] and placental telomere length [19], suggesting more pronounced effects in placentas from newborn girls. However, little is known on sex differences in human miRNA expression in early life.

David Barker was the first to show that the intrauterine environment is an important predictor of later-in-life disease risk [20]. The placenta plays a major role in foetal programming as it is the primary interface for nutrient, oxygen, waste transfer, and communication between foetus and mother, and may affect the lifetime health of the offspring. Sex-specific changes in the anatomy and the growth of placenta have already been reported from the beginning of pregnancy, with male placentas being smaller in size and having less reserve capacity towards external insults than females [21].

We enrolled 60 mother–newborn pairs to investigate the newborns’ sex-specific placental miRNA expression by whole miRNA profiling. Global gene expression in placental samples from the same study population was compared between girls and boys, in order to confirm the identified sex-specific miRNA patterns. Lastly, we studied the miRNA-RNA interactome by overlaying the outcome from the whole genome and whole miRNome studies with each other. It is of great importance to unravel possible miRNA mediated mechanisms of sex-specific processes, starting from birth, which could further link early life conditions to normal development or diseases later in life.

Results

Characteristics of the study population

In Table 1, detailed information on newborns’ and parental characteristics are provided. Mother’s and father’s age, gestational age, and season at delivery did not differ significantly between newborn girls and boys. Newborn girls and boys combined had an average (±SD) gestational age of 39.3 weeks (± 1.3). The average maternal and paternal age was 30.2 years (±4.4) and 32.1 years (±4.5), respectively.

Table 1.

Characteristics of the study population. Data are presented as mean (SD) or n (%)

Characteristics	Both sexes (n = 60)	Boys (n = 30)	Girls (n = 30)	P-value†
Maternal age, years	30.2 (4.4)	29.7 (3.1)	30.3 (4.9)	0.8
Paternal age, years	32.1 (4.5)	31.1 (3.3)	33.0 (5.3)	0.09
Newborns
Gestational age, weeks	39.3 (1.3)	39.3 (1.2)	39.3 (1.5)	0.9
Season of delivery
Winter	24 (40%)	10 (33%)	14 (46%)	0.57
Spring	11 (18%)	5 (17%)	6 (20%)
Summer	11 (18%)	6 (20%)	5 (17%)
Fall	14 (24%)	9 (30%)	5 (17%)

†P-values are determined by either two samples t-test (continuous data) or chi-square test (categorical data), between girls and boys.

Sex-specific miRNAs in newborns

Out of the 597 miRNAs which are expressed in more than 70% of the samples, 142 miRNAs were differentially expressed in placental tissue between girls and boys after FDR correction. In newborn girls, 76 miRNAs had higher expression and 66 miRNAs had lower expression compared to newborn boys (Tables 2 and 3). Among the identified differentially expressed miRNAs, miR-361-5p and miR-4646-5p had the highest difference in expression, with miR-361-5p predominantly being expressed in girls, and miR-4646-5p predominantly expressed in boys.

After additional adjustments for maternal and paternal age, gestational age, and date of delivery, 150 miRNAs were identified to be significantly differentially expressed between girls and boys, after FDR correction. Placental expression of miR-361-5p remained the most significantly different between girls and boys. Based on FDR <0.05, 129 miRNAs were identified as differentially expressed by sex in both analyses (Supplementary Table S1).

Sex-specific miRNAs and their chromosome locations

In Figure 1 the distribution of the sex-associated miRNAs is illustrated. 91.1% (n = 544) of the identified miRNAs were located on autosomes, 8.9% (n = 53) of them were found on the X chromosome (chrX), while none were present on the Y chromosome (chrY). The lowest number (0.7%) of identified miRNAs were located on chromosome 18 (chr18). Chromosome 19 (chr19) and 14 (chr14) were the chromosomes harbouring the most sex-associated miRNAs, with, respectively, 96 and 61 differentially expressed miRNAs.

Figure 1.

Manhattan plot of – log10 P-values for the identified differentially expressed miRNAs by sex (n = 597) in placental tissue distributed across different chromosomes, using Adobe Illustrator. Red line indicates the FDR threshold (FDR<0.05)

Based on the FDR-significantly differentially expressed miRNAs (n = 142) by sex in the first analysis, chr19 (16.2%), chrX (14.1%) and chromosome 1 (chr1) (9.2%) were significantly enriched for the sex-associated miRNAs. The chromosome locations of the identified sex-specific miRNAs are listed in Supplementary Table S2. Out of the 76 significant female-specific miRNAs, 16 (21.1%) miRNAs were located on chrX. Overall, mature miRNAs are equally spread over the different chromosomes, except for the Y chromosome that only harbours two mature miRNAs [22].

Pathway analysis

Overrepresentation analysis revealed 69 pathways (Supplementary Table S3) and 75 GO terms (Supplementary Table S4) which are significantly enriched by the sex-associated miRNAs. Metabolic pathway was the most significant pathway, and RNA binding was the top significant GO process. Metabolism-, hormone-, immune-, apoptosis- and neural-related pathways were identified among the top significant pathways. GO analysis identified RNA transcription, metabolism, and neurogenesis processes to be enriched by the significantly expressed sex-specific miRNAs.

Sex-specific gene expression in placenta

The gene expression profiling in the same study population (n = 60) revealed that 27 genes were significantly (FDR<5%) differentially expressed by sex. Out of these 27 genes (Table 4), 14 genes (11 protein coding, mRNAs) were found to be upregulated in placenta from girls compared with boys, while 13 other genes (9 protein coding, mRNAs) showed decreased expression. Among the protein-coding genes, EIF1AX (eukaryotic translation initiation factor 1A, X-linked) was identified as the most upregulated in placenta from girls (downregulated in boys), and DDX3Y [DEAD (Asp-Glu-Ala-Asp) box helicase 3, Y-linked] as the most upregulated in placenta from boys (downregulated in girls). Besides protein-coding genes, the non-coding gene, XIST (X inactive-specific transcript) was found as the most significant upregulated gene in placenta from girls compared with boys. With respect to the genomic coordinates of the sex-associated mRNAs, 92.8% of the girls-upregulated mRNAs were located on chrX, and 53.8% and 38.5% of the boys-upregulated (girls-downregulated) mRNAs on chrY and chrX, respectively.

Table 4.

Differentially expressed genes in placental tissue, ranked based on their significance, given by the FDR-corrected p-value (FDR). The β regression coefficient represents the fold changes in their gene expression in newborn girls compared to boys. Their chromosome locations are also provided

Gene symbol	Gene name	β regr. coefficient	FDR	Genomic coordinates
DDX3Y	DEAD (Asp-Glu-Ala-Asp) box helicase 3, Y-linked	−3.53	1.16E-246	chrY:15,027,862–15,027,921
XIST	X inactive specific transcript (non-protein coding)	8.24	8.21E-220	chrX:73,040,565–73,040,506
EIF1AY	eukaryotic translation initiation factor 1A, Y-linked	−1.34	7.94E-104	chrY:22,749,984–22,751,421
GYG2P1	glycogenin 2 pseudogene 1	−2.49	2.10E-55	chrY:14,518,023–14,517,964
KDM5D	lysine (K)-specific demethylase 5D	−2.17	1.74E-47	chrY:21,867,829–21,867,770
TSIX	TSIX transcript, XIST antisense RNA	2.66	5.35E-33	chrX:73,043,015–73,043,074
EIF1AX	eukaryotic translation initiation factor 1A, X-linked	0.54	5.25E-21	chrX:20,146,255–20,146,196
HDHD1	haloacid dehalogenase-like hydrolase domain containing 1	1.05	8.26E-21	chrX:6,967,420–6,967,361
TXLNGY	taxilin gamma pseudogene, Y-linked	−2.43	5.81E-20	chrY:21,766,231–21,766,290
KDM6A	lysine (K)-specific demethylase 6A	0.96	5.10E-16	chrX:44,971,784–44,971,843
CD99	CD99 molecule	−1.13	7.78E-15	chrX:2,638,454–2,640,715
DDX3X	DEAD (Asp-Glu-Ala-Asp) box helicase 3, X-linked	0.53	1.12E-13	chrX:41,208,627–41,208,686
PCDH11X	protocadherin 11 X-linked	−1.29	9.96E-13	chrX:91,134,246–91,134,305
PCDH11Y	protocadherin 11 Y-linked	−1.24	1.08E-11	chrY:4,968,462–4,968,521
VAMP7	vesicle-associated membrane protein 7	−0.48	8.82E-10	chrX:155,172,553–155,172,612
STS	steroid sulfatase (microsomal), isozyme S	1.01	3.66E-08	chrX:7,271,727–7,271,786
NAA10	N(alpha)-acetyltransferase 10, NatA catalytic subunit	0.36	2.07E-07	chrX:153,197,847–153,197,788
UBA1	ubiquitin-like modifier activating enzyme 1	0.41	4.49E-06	chrX:47,073,739–47,073,798
HSD17B10	hydroxysteroid (17-beta) dehydrogenase 10	0.35	7.30E-05	chrX:53,458,440–53,458,381
PRKY	protein kinase, Y-linked, pseudogene	−0.78	1.07E-04	chrY:7,248,885–7,248,944
ASMTL	acetylserotonin O-methyltransferase-like	−0.37	3.43E-04	chrX:1,522,197–1,522,138
SMC1A	structural maintenance of chromosomes 1A	0.53	3.43E-04	chrX:53,426,606–53,426,547
KDM5C	lysine (K)-specific demethylase 5 C	0.52	2.24E-03	chrX:53,253,984–53,253,925
DNM1P46	DNM1 pseudogene 46	−1.00	4.03E-03	chr15:100,331,156–100,331,097
ARMCX1	armadillo repeat containing, X-linked 1	−0.76	5.20E-03	chrX:100,809,416–100,809,475
CDK16	cyclin-dependent kinase 16	0.30	6.86E-03	chrX:47,089,333–47,089,392
TRERNA1	translation regulatory long non-coding RNA 1	0.60	2.30E-02	chr20:48,657,346–48,657,287

miRNA-mRNA interactions

Further, we explored the interactions between significantly differentially expressed miRNAs (n = 142, from the main analysis) and mRNAs (n = 20) using the miRComp tool [23] in R. Figure 2 illustrates the identified miRNA–mRNA interactions (n = 14), in which the red colour indicates upregulation, while the green indicates downregulation in girls compared with boys. More detailed information on the miRNA-RNA interactome is given in Supplementary Table S5. The strongest inverse correlation (r = −0.58; p < 0.001) was found between the miR-361-5p, which was the most significantly upregulated miRNA in girls compared with boys, and PCDH11X (protocadherin 11 X-linked), which was identified as the most significantly downregulated gene in girls (upregulated in boys). Based on the highest score (score = 3.12), hsa-miR-4287 and DDX3Y had the highest opposite deregulation in girls compared with boys.

Figure 2.

The miRNA-mRNA sex interactome. Significantly (FDR<5%) identified miRNA-mRNA interactions (n = 14) obtained by miRComp, using Cytoscape v3.7.2. Rectangle nodes indicate miRNAs (as sources) and diamond nodes indicate mRNA (as targets). The red colour indicates upregulation, while the green downregulation in girls compared with boys. The width of the arrows connecting rectangle and diamond nodes corresponds to the calculated scores, with thicker arrows indicating higher absolute log ratio scores

Discussion

Our key findings suggest that placental miRNA profiling occurs in a sex-specific manner. Epigenetic changes of sex differences, including DNA methylation and histone acetylation, in human samples including blood, brain, and saliva, have already been reported [24-26]. Early life epigenetic mechanisms have been suggested to have a great impact on foetal growth and development of non-communicable disease in later life [3,27]. Although the placenta is neither male nor female, sex-specific placental responses to in utero environmental conditions have been reported in foetal development [28,29], which could explain the differences in gestational vulnerability to external stressors between girls and boys [30] and sex-specific disease susceptibility in later life. Evidence from several animals and human studies has demonstrated sex-specific miRNA expression patterns in normal development or diseases [13,31-34]. Therefore, the identification of sex-specific miRNA mechanisms in healthy placenta is of great importance to reveal physiological patterns of the foetus or the newborn, as these sex-specific differences could be important in disease development later in life.

Additional adjustment for possible covariates including maternal age, paternal age, gestational age, and date of delivery did not alter the outcome. Most of the observed differentially expressed miRNAs by sex (91%) were significantly expressed in both analyses. One of the most differentially expressed miRNA by sex, miR-361-5p (located on chrX), was found to be highly expressed in placenta from girls. This miR-361-5p has earlier been identified as female-biased miRNA in cerebellum [31]. MiR-361-5p was found to be increased in hypoxia-induced pulmonary artery smooth muscle cells [35], and in blood serum from Alzheimer’s disease patients in two independent cohorts [36]. Higher expression of miR-361-5p was also found in plasma from patients with coronary artery disease (n = 49) compared to healthy controls (n = 27) [37]. MiR-361-5p has been reported to play a critical role in several human tumours [38]. A genome-wide sequencing project indicated downregulation of miR-361-5p in plasma of male late-onset hypogonadism (LOH) patients (n = 22) compared to healthy subjects (n = 22), suggesting miR-361-5p as a novel biomarker for diagnosis of deficiency of serum testosterone [39].

Among the top differentially expressed miRNA in placenta by sex, miR-361-3p, has also been shown as female-specific in peripheral blood [31]. In addition, sex-specificity of miR-340-5p was found in the placenta of female prenatally stressed mice, based on DNA global methylation study [40]. While miR-526b (−3p) was indicated to be male-specific in cerebellum [31].

On the other hand, miR-4646-5p (located on chr6) was the most differentially expressed miRNA in placenta between boys and girls, with higher expression in boys. Based on genome-wide DNA methylation, miR-4646-3p was associated with depression in adolescence [41]. In another study, miR-4646-3p was identified to be upregulated in peripheral blood from schizophrenic patients (n = 55) compared to healthy controls (n = 28) [42].

Moreover, within the top significantly upregulated placental miRNAs in boys, miR-1915-3p and miR-575 have also been shown to be male-specific in peripheral blood and cerebellum, respectively [31]. While the miR-7107-5p has been identified in synovial fluid of female patients with osteoarthritis [43].

Remarkably, some of the identified placental sex-specific miRNAs have been previously reported to be sex–associated in other tissues as well (Table 5). Lizarraga and co-authors [10] have reported 94 significantly differentially miRNAs by sex measured in cord blood samples (n = 89), among them, 19 miRNAs were found to be consistent with our significant findings in newborn boys, while miR-452-5p was significantly differentially expressed in opposite sexes (Figure 3). In another study that reported on sex-differential expression of miRNAs in human brain development (n = 18) [9], 40 miRNAs were found significantly differentially expressed in the prefrontal cortex over time, of which 5 miRNAs were identified also in placental tissue (FDR<0.1), but only miR-4787 (in late childhood and adolescence) and miR-29c (in early childhood) were in the same direction in both tissues. In adult life, four out of nine investigated miRNAs were differentially expressed by sex in blood platelets, which were also differentially expressed in the placenta by sex (FDR<0.1) in the current study [44], with miR-145 as girls-upregulated in both placenta and blood platelets.

Table 5.

Reported sex-specific miRNAs overlapped with the findings of our study

Study	Tissue	Commonly expressed miRNAs in our study
Lizarraga et al., 2016	Cord blood	miR-663a, miR-6727-5p, miR-1207-5p, miR-4281, miR-6803-5p, miR-4516, miR-4505, miR-6088, miR-4739, miR-1268b, miR-3656, miR-638, miR-6086, miR-6780b-5p, miR-4534, miR-2861, miR-6786-5p, miR-4442, miR-6821-5p
Ziats et al., 2014	Prefrontal cortex	miR-145-3p, miR-92a-3p, miR-29 c-5p, miR-4787-5p and miR-30 c-1-3p
Simon et al., 2014	Blood platelets	miR-520e, miR-497, miR-29b, and miR-145

Figure 3.

Venn Diagram of common differentially expressed miRNAs by sex in cord blood [10] and in placental tissue. After FDR correction (<0.05), a total of 93 and 142 miRNAs were differentially expressed in cord blood and placental tissue respectively. Twenty sex-specific miRNAs were found in either tissue. In detail, 19 miRNAs including miR-6727-5p, miR-4442, miR-4281, miR-6780b-5p, miR-1207-5p, miR-2861, miR-3656, miR-4505, miR-4516, miR-1268b, miR-4739, miR-6786-5p, miR-6088, miR-638, miR-6803-5p, miR-663a, miR-4543, miR-6821-5 and miR-6086, were significantly lower expressed in girls compared to boys in either tissue. However, only one overlapped miRNA in both studies, namely miR-452-5p, was found to be significantly differentially expressed in opposite sexes

Moreover, the gene expression analysis revealed that all significantly differentially expressed genes by sex in placenta were located on sex-chromosomes, either on chrX or on chrY. Consistent findings in human placenta have been reported earlier by Gonzalez et al. [45], confirming our findings of 7 upregulated genes (EIF1AX, KDM6A, STS, SMC1A, KDM5C, DDX3X, KDM6A, XIST), and 5 downregulated genes (DDX3Y, KDM5D, TXLNGY, PCDH11Y, EIF1AY) by comparing girls to boys. In addition, in another study [46], human placental transcriptome showed sexual dimorphic gene expression including DDX3X, KDM5D, and EIF1AY with higher expression in boys, and KDM6A with higher expression in girls. These common genes were confirmed as sex-linked genes in both late first-trimester placenta [45] and term placenta (our study), emphasizing their sex-specificity already from the early pregnancy.

Further, the integrative analysis of the placental miRNA-mRNA interactome showed 14 significant (FDR<0.05) interactions between the identified sexually differentially expressed miRNAs and mRNAs. The girls-downregulated (boys-upregulated) PCDH11X gene was inversely correlated with three miRNAs, including hsa-miR-361-5p, hsa-miR-139-5p, and hsa-miR-4769-3p, which were found to be upregulated in girls. Although the PCDH11X, the X-linked homolog of PCDH11Y, was not found to be associated with sex in first-trimester placenta [45], it was shown to be boys-upregulated in term placenta based on our findings. Moreover, KDM6A, upregulated in girls, was negatively correlated with three downregulated miRNAs (hsa-miR-1207-5p, hsa-miR-4687-3p and hsa-miR-4530). It has been reported that the girls-upregulated X-linked KDM6A gene can affect chromatin modification and is conserved during pregnancy and into adulthood [47]. Likewise, the girls-downregulated CD99 was inversely correlated with the six girls-upregulated miRNAs (hsa-miR-520e, hsa-miR-518d-3p, hsa-miR-331-3p, hsa-miR-423-3p, hsa-miR-342-5p, and hsa-miR-330-3p). CD99 gene, encoding a transmembrane glycoprotein, is known to be highly expressed in placenta [48], and in consistency with our study, it has been shown to be overexpressed in men, in human cell lines by a large-scale population study [49] and in inflammatory cytokines [50]. The girls-downregulated DDX3Y, DEAD (Asp-Glu-Ala-Asp) box helicase 3 was inversely correlated with the girls-upregulated hsa-miR-4287. As aforementioned, the Y-linked DDX3Y gene has been reported to show significant higher expression in male than female human placentas [46]. Lastly, the girls-upregulated (boys-downregulated) CDK16, cyclin-dependent kinase 16, was negatively correlated with girls-downregulated (boys-upregulated) hsa-miR-4287. However, there are no available data on the sex-specific pattern of CDK16.

Enrichment analysis of the sex-specific significantly expressed miRNAs revealed pathways involved in metabolism, gonadotropin-releasing hormone (GnRH) signalling, depression, apoptosis signalling, activation of B and T cell, MAPK signalling, and axon guidance. Additionally, GO enrichment analysis showed the top processes to be related to RNA transcription, metabolism, and neurogenesis. Taken together, our enrichment analysis revealed a wide range of processes to be related to sex-specific processes in foetal tissues. Twelve common GO biological processes in placenta and cord blood [10] were identified to be regulated by the sex-specific miRNAs, involved mainly in RNA transcription, gene expression, and cellular metabolic processes (data not shown).

Sex differences in foetal developmental programming of metabolism have previously been described, indicating that females are more susceptible to develop increased adiposity and impaired glucose haemostasis in response to exposure to under- or over-nutrition during early life [51]. Moreover, differences in myocardial metabolism between newborn girls and boys have been reported, with girls to exhibit lower energy levels and greater tissue lactic acidosis, both linked to an increased susceptibility to ischaemic injury and weak myocardial function [52]. Sex differences have also been reported in neurodevelopment in either newborns [10] or older children [53]. This suggests an important role in brain development which may explain cognitive processes in early life. Placenta GnRH can influence steroid release and varies according to the gestational age [54].

Furthermore, evident early differences in development showed faster cell divisions in XY than XX embryos in mammalian species, including humans [55]. Sex-specific placental responses have been reported to play an important role in placental functions and pathologies [28]. Clifton et al. [56,57] have shown altered placental glucocorticoid metabolism in female foetuses of normal pregnancies, resulting in decreased foetal growth, besides altered placental cytokine expression, and altered insulin-like growth factor pathways, in response to the presence of maternal asthma. Human placental transcriptomic analyses have also revealed sex differences in gene expression, with females expressing more genes located on autosomal chromosomes, specifically many immune-related genes [58]. Distinct sex-specific placental gene expression is involved in inflammation in association with the development of allergies in childhood (at age 4 years), affecting the foetal immune system in early life [59]. In another study, enrichment of signalling pathways previously reported to mediate graft-versus-host disease and other transcripts involved in immune function and inflammation were identified in the placenta in boys [60]. Lastly, differentially methylated genes have also been observed in a sex-specific manner, resulting in sex-specific differences in the placental transcriptome [61].

This is the first study in assessing sex-specific miRNA patterns in placenta tissue from normal pregnancies. Our findings are supported by possible links of sex-specific processes earlier mentioned in foetal tissues with the identified pathways (Table 6). In addition, placental gene expression confirmed the sex-specificity of miRNA targets. Sex-specific miRNA expression in the placenta was identified by a high-throughput method for miRNA profiling. The obtained results were corrected for possible technical errors such as hybridization effect which could influence our statistical analysis. Most of our observed findings remained consistent across the two performed statistical analyses, which strengthens our analysis. Nevertheless, we cannot exclude the possibility that other factors bias our analysis such as unknown or unmeasured variables.

Table 6.

Different processes in human placenta or cord blood are linked to sex-specific responses and to identified pathways and/or GO terms. The ‘x’ displays the sex of the newborn, in which each process has been identified as specific

Processes/Diseases	Girls	Boys	Potential link to identified pathways or GO terms	Ref.
Placental glucocorticoid metabolism (Placenta)	x (normal pregnancy)		hsa04912_GnRH_signalling_pathwayP04373_5HT1_type_receptor_mediated_signalling_pathwayhsa04720_Long_term_potentiationhsa04914_Progesterone_mediated_oocyte_maturation	[70]
Glucocorticoid receptor a D2 (GRaD2) (preterm Placenta)		x (maternal stress)		[71]
Human chorionic gonadotropin (HCG) (Placenta)	X			[72]
11β-hydroxysteroid dehydrogenase type 2 (11βHSD2) & Cortisol (Placenta)	x (after glucocorticoid treatment)			[73]
Oxidative stress (Placenta)		x (after glucocorticoid treatment)	WP111_Electron_Transport_Chain	[74]
Immune function (Placenta)	X		P00010_B_cell_activationP00053_T_cell_activationhsa04664_Fc_epsilon_RI_signalling_pathwayhsa04662_B_cell_receptor_signalling_pathwayhsa04660_T_cell_receptor_signalling_pathwayP00031_Inflammation_mediated_by_chemokine_and_cytokine_signalling_pathway	[75]
Immunological disease and immune cell trafficking networks (Placenta)	X			[59]
Allergy (early childhood): Dermatological and respiratory disease networks (Placenta)		x
Immune function-Graft vs. host diseases (Placenta)		x		[76]
Cytokine expression- Maternal asthma (Placenta)	x			[77]
Metabolism
Foetal growth-Maternal asthma (Placenta)	x		GO0008286_insulin_receptor_signalling_pathway hsa01100_Metabolic_pathwayshsa00270_Cysteine_and_methionine_metabolismhsa00561_Glycerolipid_metabolismGO0005975_carbohydrate_metabolic_process GO0006629_lipid_metabolic_processhsa00600_Sphingolipid_metabolismhsa00010_Glycolysis__Gluconeogenesishsa01040_Biosynthesis_of_unsaturated_fatty_acids	[56,70]
Maternal obesity (Placenta)	X			[78]
Insulin-like growth factor (Cord blood)	x (maternal asthma)		GO0008286_insulin_receptor_signalling_pathway hsa01100_Metabolic_pathways	[79]
Renin-Angiotensin system (RAS) (Placenta)	x		hsa04370_VEGF_signalling_pathwayhsa04270_Vascular_smooth_muscle_contractionWP536_Calcium_Regulation_in_the_Cardiac_Cell	[80]

In conclusion, these findings can provide more insights into molecular mechanisms underlying sex differences already present in newborns, which could contribute to a better understanding of physiological or pathological conditions in early or later life. As metabolic pathways and processes are regulated by the differentially expressed miRNAs by sex, our observations could indicate a novel mechanism to explain differential disease susceptibility by sex.

Materials and methods

Study population

We selected 60 participants (30 boys) from the ENVIRONAGE (ENVIRomental influence ON early AGEing) birth cohort [62]. The recruitment procedure of the mother–newborn pairs was approved by the Ethics Committees of Hasselt University and East-South-Limburg Hospital and was performed according to the declaration of Helsinki. The participants were recruited before delivery (within the period of January 2014 – April 2017), and provided written informed consent and filled out questionnaires after delivery. In these questionnaires, information about maternal and paternal age, occupation, education, smoking behaviour, alcohol consumption, place of residence, use of medication, newborn’s ethnicity, and parity was given. The newborn’s ethnicity was defined as of European-Caucasian origin when at least two grandparents were European, while as non-European when at least three grandparents were non-European.

Experimental design

In our main study population (n = 60), we only included mother–newborn pairs with placentas of gestational age ≥36 weeks, for which the newborn was of European-Caucasian ethnicity, and the mothers were never-smokers. We had an equal number of newborn girls (n = 30) and boys (n = 30).

Sample collection, m(i)RNA isolation

Fresh placental tissue was collected within 1 h after delivery. Four standardized biopsies were taken at fixed locations across the middle point of the placenta at around 4 cm distance from the umbilical cord, on the foetal side as detailed by Janssen and colleagues [62,63]. The freshly collected biopsies were stored in RNA later (ThermoFisher Scientific, Massachusetts, USA) at 4°C for at least 12 hours, and maximally 1 week, after which the biopsies were stored at −20°C until extraction.

Total RNA and miRNA were extracted from a single fresh placental biopsy using the miRNeasy mini kit (Qiagen, Venlo, The Netherlands) according to the manufacturer’s protocol. RNase-Free DNase treatment was performed on RNA samples according to the manufacturer’s instructions (Qiagen, Venlo, The Netherlands). RNA quantity and sample purity were assessed by spectrophotometry (Nanodrop 1000, Isogen, Life Science, Belgium) and RNA integrity by Agilent 2100 Bioanalyzer (Agilent Technologies, Amstelveen, the Netherlands). All RIN (RNA integrity number) values were above 6, as required for good quality microarray-based analysis [64,65].

miRNA expression profiling by microarrays

After labelling, total RNA was hybridized onto the Sureprint G3 Human V19 miRNA Agilent 8 × 60 K microarrays. Raw data were extracted by Agilent Feature Extraction Software. The corresponding data quality assessment, preprocessing, and quantile normalization were performed by R package AgiMicroRna [66,67]. miRNAs with less than 70% present data over all samples were omitted and the missing values were imputed by K-nearest neighbour imputation (K = 15). After filtering of 2558 miRNAs, expression data of 597 miRNAs remained for further analyses.

Global gene expression profiling by microarrays

Samples of RNA samples were hybridized onto Agilent Whole Human Genome 8 × 60 K microarrays, which were scanned by an Agilent G2505 C DNA Microarray Scanner (Agilent Technologies, Amstelveen, Netherlands). Scan images were converted into pixel intensities using the Agilent Feature Extraction Software (Version 10.7.3.1, Agilent Technologies, Amstelveen, Netherlands). Raw intensity data were pre-processed by an in-house-developed quality control pipeline in R software as follows: local background correction, omission of controls, flagging of bad spots and spots with too low intensity, log2-transformation, and quantile normalization using arrayQC. The R-scripts of the pipeline and information about the flagging are available at https://github.com/BiGCAT-UM/arrayQC_Module. Based on the processed data, probes showing >30% flagged data were excluded, after which replicate probes were merged based on the median and missing values were imputed by means of K-nearest neighbour imputation (K = 15). For genes with multiple probes, only the probe with the largest interquartile range (IQR) was selected, leading to a dataset with 18,847 probes. We further excluded probes that were not annotated by a known gene, leaving 17,672 probes for the following analysis.

Statistical analysis

Using a mixed linear regression model (RStudio, R version 3.6.3), in a first analysis, correcting only for the hybridization batch effect, the association of newborn’s sex with miRNA expression was investigated in placental tissue. In a second analysis, we additionally adjusted our model for maternal and paternal age, gestational age, and date of delivery. The obtained regression coefficients indicate the difference in the placental miRNA expression in newborn girls compared to newborn boys (i.e., the reference in our model). P-values were adjusted for multiple testing by the Benjamini-Hochberg false discovery rate (FDR), considering a level of significance of 5%. Categorical data were presented as frequency and percentage, whereas continuous data were presented as mean (SD).

Similarly, the differential mRNA expression by sex was assessed in the same mother–newborn pairs (n = 60), using linear mixed model (RStudio, R version 3.6.3) correcting for the hybridization effect. Probes with FDR<0.05 were statistically differentially expressed (DE) mRNAs by sex.

Significant differentially expressed miRNAs (n = 142) and mRNAs (n = 20) excluding non-coding RNAs and transcripts of pseudogenes were used to conduct the analysis of miRNA–mRNA interactions with an R package miRComp [23]. Briefly, pairwise Pearson correlation coefficients between all miRNA-mRNA pairs were calculated. After correction for FDR, the significant interaction pairs were matched to target prediction databases [TargetScan [81] and MicroCosm [82,83], with a total of 563,179 interactions]. An interaction was assumed to be present if there was a significantly negative correlation between the pair and the target was found in at least one database. The network of miRNA–mRNA interaction, referred to as the miRNA-mRNA interactome, was plotted by Cytoscape v3.7.2 [84]. The log ratios (fold-changes comparing girls to boys) obtained significantly by DE of mRNAs and miRNAs were used to calculate a score for each interaction pair, indicating whether the mRNA and the miRNA in a pair were concordant or opposite in their deregulation, with the following formula: score = −2(logratio_miRNA×logratio_mRNA).

Genomic coordinates of sex-specific miRNAs

Using miRBase v.22 database [1] of published miRNA sequences and annotation, we retrieved information about the chromosome location of the sex-associated miRNAs in placental tissue. Manhattan plot, using R package qqman [68] and Adobe Illustrator, was used to visualize the distribution of the identified differentially expressed miRNAs by sex (n = 597) in placental tissue across all different chromosomes. The significantly differentially expressed miRNAs by sex are provided after applying an FDR threshold (FDR<0.05).

Enrichment and pathway analysis

MiRNomics from GeneTrail2 (v1.6) [69] was used for enrichment analysis of differentially expressed miRNA (n = 142) by sex in placental tissue, based on the outcome of the first analysis. GeneTrail2 is a web-interface providing access to tools for statistical analysis of molecular signatures. MiRNomics combines information from various databases (e.g., miRWalk, miRTarBase, miRDB) to perform overrepresentation analysis. Significant (FDR<0.05) overrepresented pathways and GO processes were identified.

Click here for additional data file.

Table 2.

Sex-specific miRNAs with higher expression in placental tissue in newborn girls compared to boys, top 35 miRNAs, ranked based on their significance (P-value). Mean log2 miRNA expression in both or separated sexes with variances are provided

miRNA	Mean log2 expression (overall)	Variance (overall)	Mean log2 expression (boys)	Variance (boys)	Mean log2 expression (girls)	Variance (girls)	β regr. coefficient	P-value	FDR
hsa-miR-361-5p	5.54	0.27	5.18	0.13	5.89	0.17	0.68	2.35E-09	1.40E-06
hsa-miR-139-5p	1.31	0.59	0.94	0.76	1.66	0.19	0.76	5.88E-05	1.75E-02
hsa-miR-664a-3p	2.93	0.27	2.66	0.28	3.21	0.12	0.48	1.12E-04	1.76E-02
hsa-miR-361-3p	2.36	0.56	1.99	0.67	2.73	0.18	0.67	1.39E-04	1.76E-02
hsa-miR-340-5p	1.72	0.53	1.25	0.50	2.13	0.20	0.63	1.47E-04	1.76E-02
hsa-miR-374 c-5p	1.07	0.57	0.63	0.59	1.45	0.25	0.64	2.86E-04	2.26E-02
hsa-miR-520e	0.70	0.54	0.40	0.62	1.02	0.27	0.64	2.88E-04	2.26E-02
hsa-miR-526b-3p	1.12	0.58	0.80	0.49	1.46	0.46	0.67	3.04E-04	2.26E-02
hsa-miR-18b-5p	0.81	0.46	0.50	0.55	1.10	0.21	0.58	4.30E-04	2.68E-02
hsa-miR-30 c-2-3p	1.26	0.33	1.02	0.44	1.50	0.10	0.46	5.61E-04	2.68E-02
hsa-miR-1323	10.41	0.09	10.29	0.09	10.53	0.05	0.25	7.48E-04	2.68E-02
hsa-miR-186-5p	3.16	0.22	2.93	0.18	3.38	0.17	0.38	7.66E-04	2.68E-02
hsa-miR-7-5p	1.23	0.56	0.93	0.69	1.55	0.24	0.58	7.96E-04	2.68E-02
hsa-miR-24-1-5p	1.76	0.49	1.49	0.54	2.05	0.29	0.57	9.94E-04	2.84E-02
hsa-miR-518 f-5p	5.94	0.16	5.77	0.14	6.11	0.13	0.33	1.01E-03	2.84E-02
hsa-miR-934	1.00	0.55	0.71	0.57	1.29	0.37	0.60	1.03E-03	2.84E-02
hsa-miR-377-5p	0.74	0.49	0.46	0.64	1.01	0.22	0.56	1.06E-03	2.84E-02
hsa-miR-338-3p	1.21	0.43	0.93	0.42	1.50	0.29	0.54	1.13E-03	2.84E-02
hsa-miR-766-3p	2.75	0.25	2.55	0.23	2.95	0.20	0.41	1.14E-03	2.84E-02
hsa-miR-378i	4.10	0.23	3.90	0.22	4.29	0.17	0.39	1.19E-03	2.85E-02
hsa-miR-141-5p	0.74	0.55	0.42	0.54	1.05	0.37	0.60	1.26E-03	2.86E-02
hsa-miR-30d-3p	0.85	0.65	0.52	0.57	1.18	0.54	0.64	1.36E-03	2.86E-02
hsa-miR-4287	4.48	0.29	4.26	0.28	4.69	0.20	0.42	1.49E-03	2.86E-02
hsa-miR-371a-3p	2.67	1.46	2.23	1.86	3.19	0.54	0.85	1.53E-03	2.86E-02
hsa-miR-378a-3p	4.31	0.23	4.12	0.21	4.50	0.18	0.38	1.61E-03	2.86E-02
hsa-miR-512-5p	7.38	0.22	7.20	0.20	7.57	0.17	0.37	1.63E-03	2.86E-02
hsa-miR-197-3p	3.77	0.26	3.57	0.32	3.98	0.13	0.40	1.71E-03	2.87E-02
hsa-miR-518a-5p	3.72	0.35	3.49	0.30	3.96	0.30	0.46	1.73E-03	2.87E-02
hsa-miR-498	7.26	0.23	7.07	0.25	7.45	0.14	0.38	1.81E-03	2.93E-02
hsa-miR-629-5p	0.81	0.55	0.54	0.67	1.07	0.32	0.54	2.05E-03	3.04E-02
hsa-miR-376a-5p	1.77	0.40	1.54	0.40	2.02	0.30	0.49	2.07E-03	3.04E-02
hsa-miR-4769-3p	1.59	0.48	1.34	0.59	1.85	0.24	0.52	2.11E-03	3.04E-02
hsa-miR-452-5p	3.80	0.34	3.57	0.30	4.02	0.28	0.45	2.14E-03	3.04E-02
hsa-miR-629-3p	0.60	0.50	0.32	0.51	0.90	0.33	0.55	2.19E-03	3.04E-02
hsa-miR-181a-3p	1.85	0.68	1.52	0.80	2.18	0.35	0.62	2.26E-03	3.04E-02

Table 3.

Sex-specific miRNAs with lower expression in placental tissue in newborn girls compared to boys, top 35 miRNAs, ranked based on their significance. Mean log2 miRNA expression in both or separated sexes with variances are provided

miRNA	Mean log2 expression (overall)	Variance (overall)	Mean log2 expression (boys)	Variance (boys)	Mean log2 expression (girls)	Variance (girls)	β regr. coefficient	P-value	FDR
hsa-miR-4646-5p	2.08	0.31	2.30	0.29	1.85	0.22	−0.47	4.65E-04	2.68E-02
hsa-miR-4793-5p	2.46	0.22	2.67	0.22	2.25	0.14	−0.37	4.97E-04	2.68E-02
hsa-miR-6778-5p	2.39	0.58	2.69	0.59	2.03	0.36	−0.62	5.92E-04	2.68E-02
hsa-miR-7107-5p	4.33	0.42	4.61	0.43	4.06	0.27	−0.55	6.78E-04	2.68E-02
hsa-miR-1915-3p	5.87	0.64	6.22	0.77	5.52	0.28	−0.63	8.09E-04	2.68E-02
hsa-miR-663a	2.72	0.89	3.10	0.84	2.34	0.67	−0.76	1.34E-03	2.86E-02
hsa-miR-575	3.69	0.47	3.96	0.51	3.41	0.29	−0.55	1.48E-03	2.86E-02
hsa-miR-6757-5p	1.74	0.51	2.03	0.63	1.42	0.20	−0.54	1.50E-03	2.86E-02
hsa-miR-3141	2.50	0.30	2.71	0.30	2.29	0.23	−0.42	2.25E-03	3.04E-02
hsa-miR-6727-5p	3.99	0.63	4.28	0.85	3.69	0.27	−0.61	2.54E-03	3.04E-02
hsa-miR-1207-5p	5.73	0.58	6.02	0.70	5.44	0.31	−0.58	2.57E-03	3.04E-02
hsa-miR-1225-5p	6.40	0.51	6.68	0.56	6.11	0.32	−0.54	2.58E-03	3.04E-02
hsa-miR-4298	2.69	0.37	2.91	0.47	2.46	0.17	−0.47	2.67E-03	3.04E-02
hsa-miR-4281	6.79	0.62	7.09	0.79	6.49	0.30	−0.60	2.81E-03	3.04E-02
hsa-miR-7847-3p	2.93	0.32	3.14	0.39	2.72	0.17	−0.43	3.02E-03	3.04E-02
hsa-miR-6510-5p	2.98	0.46	3.23	0.49	2.72	0.31	−0.51	3.18E-03	3.04E-02
hsa-miR-6775-5p	3.31	0.43	3.57	0.51	3.05	0.22	−0.47	3.24E-03	3.04E-02
hsa-miR-6800-5p	5.83	0.73	6.16	0.94	5.49	0.32	−0.63	3.31E-03	3.04E-02
hsa-miR-6803-5p	5.42	0.55	5.73	0.58	5.11	0.34	−0.50	3.35E-03	3.04E-02
hsa-miR-4530	6.56	0.69	6.90	0.78	6.22	0.38	−0.55	3.37E-03	3.04E-02
hsa-miR-4741	5.76	0.49	6.01	0.41	5.50	0.45	−0.52	3.51E-03	3.04E-02
hsa-miR-4800-5p	3.02	0.43	3.28	0.45	2.76	0.27	−0.47	3.53E-03	3.04E-02
hsa-miR-6867-5p	2.39	0.29	2.54	0.33	2.21	0.18	−0.40	3.57E-03	3.04E-02
hsa-miR-4516	9.85	1.57	10.32	1.62	9.39	1.13	−0.93	3.64E-03	3.04E-02
hsa-miR-619-5p	2.31	0.69	2.63	0.74	2.00	0.46	−0.61	3.66E-03	3.04E-02
hsa-miR-6087	8.39	1.52	8.86	1.96	7.92	0.67	−0.89	3.81E-03	3.04E-02
hsa-miR-4669	3.30	0.18	3.46	0.19	3.13	0.12	−0.29	4.36E-03	3.25E-02
hsa-miR-4299	7.64	0.74	7.94	0.87	7.34	0.44	−0.62	4.50E-03	3.32E-02
hsa-miR-6768-5p	2.79	1.10	3.16	1.23	2.42	0.72	−0.75	4.82E-03	3.48E-02
hsa-miR-1972	1.92	0.60	2.20	0.68	1.58	0.31	−0.49	5.04E-03	3.50E-02
hsa-miR-371a-5p	3.38	0.59	3.65	0.59	3.12	0.47	−0.55	5.12E-03	3.51E-02
hsa-miR-4505	5.54	0.85	5.88	1.18	5.19	0.31	−0.62	5.28E-03	3.56E-02
hsa-miR-6088	5.54	0.36	5.79	0.40	5.29	0.21	−0.39	5.41E-03	3.59E-02
hsa-miR-6858-5p	3.85	0.97	4.19	0.93	3.50	0.81	−0.69	5.63E-03	3.65E-02
hsa-miR-4739	6.59	1.84	7.06	2.12	6.11	1.16	−0.95	5.99E-03	3.74E-02