hsa‑miR‑155‑5p concentration as a determinant of leptinemia in umbilical cord plasma

Introduction

The prevalence of overweight and obesity among children and adolescents aged 5-19 has increased dramatically, rising from 8% in 1990 to 20% in 2022(1), and now reaching alarming proportions worldwide. In 2022, the World Health Organization (WHO) estimated that 37 million children under the age of 5 years were overweight globally. Overweight and obesity in childhood and adolescence are associated with an increased risk of developing type 2 diabetes, cardiovascular diseases and mental health disorders early in life (2). Children affected by obesity are also five times more likely to maintain excess weight into adulthood, placing them at higher risk for metabolic diseases and certain cancers (3), among other conditions. Therefore, it is essential to investigate the mechanisms underlying weight gain and adiposity regulation early in life to better understand and prevent obesity and its related comorbidities.

Leptin is a satiety hormone primarily synthesized by white adipose tissue and secreted in circulation proportionally to an individual's fat mass (4). Leptin binds to its soluble leptin receptor (LEPR) isoform to circulate in the bloodstream and interacts with LEPR transmembrane isoforms to trigger specific physiological responses or cross the blood-brain barrier to reach the central nervous system (CNS). In the hypothalamus, leptin suppresses appetite and increases energy expenditure (5). In peripheral tissues such as the liver, skeletal muscles, pancreas, adipose tissue, and within the cardiovascular and immune systems, leptin contributes to the regulation of lipolysis, angiogenesis, inflammation and immune responses, among other pathways (6-11). During pregnancy, leptin also plays a role in programming hypothalamic neural circuits responsible for food intake regulation, with long-term effects on offspring's adiposity (12,13). Indeed, leptinemia rises dramatically after 34 weeks of gestation, coinciding with a significant increase in fetal fat mass during this period. This represents a critical window of leptin sensitivity, essential for neuronal programming, as previously demonstrated in mice (14,15). Interestingly, leptinemia is elevated in individuals with obesity without suppressing food intake, suggesting central leptin resistance (16).

The mechanisms regulating leptin production, secretion and sensitivity are not yet fully understood (5,17), but several studies suggest that microRNAs (miRNAs) may play a role (18-20). MiRNAs are small non-coding RNAs, typically ~19-25 nucleotides in length, that are secreted and stable in the bloodstream. They generally bind the 3' untranslated region (UTR) of messenger RNA (mRNA) targets, inhibiting protein synthesis through mRNA degradation or translational inhibition. A single miRNA often targets multiple mRNAs, thereby influencing numerous metabolic pathways (21).

To date, no study has investigated the association between miRNA profiles and leptinemia in umbilical cord plasma (UCP) to identify potential regulators of leptin during this critical developmental window. We hypothesized that circulating miRNAs are involved in the regulation of leptinemia in UCP. Our objective was to assess the association between UCP leptinemia and miRNAs using Next Generation Sequencing (NGS) and validate these associations with specific miRNAs by re-measuring using an alternative quantification technology, digital PCR (dPCR).

Materials and methods

Participant selection in Gen3G Cohort

The participants were selected from the Genetics of Glucose regulation on Gestational and Growth (Gen3G) prospective pregnancy and early-life cohort, which consists of mother-child pairs recruited between January 2010 and June 2013(22). Briefly, pregnant women aged 18 years and older, having singleton pregnancy and free of pre-gestational diabetes or medication affecting glucose tolerance were recruited at the Centre Hospitalier Universitaire de Sherbrooke (CHUS) between the 5th and 16th weeks of pregnancy. Gestational Diabetes Mellitus (GDM) was diagnosed between the 24th and 28th weeks of pregnancy using a 75 g oral glucose tolerance test (OGTT) according to the criteria of the International Association of Diabetes and Pregnancy Study Groups (IADPSG). Among the 854 women followed until delivery with complete OGTT data, additional inclusion criteria were applied for this study. These included the availability of child follow-up at 5 years with adiposity measurements by Dual-energy X-ray Absorptiometry (DXA), the availability of UCP samples that had never been thawed to minimize RNA degradation, and the availability of relevant maternal and child health (e.g., GDM status, UCP leptin measurement at birth). Applying these criteria resulted in a final sample size of 273 participants for the sequencing analysis, of which 214 participants with non-thawed UCP available were selected for dPCR analysis. All participants provided written informed consent for the collection of samples and data for research purposes, for both themselves and their newborns. The study protocol was approved by the CHUS ethics committee under reference number 2010-198, 07-027-A1 (original approval on 2009-10-27, renewed annually since).

Anthropometric and biochemical measurements

UCP samples and anthropometric measurements were collected, and health and lifestyle questionnaires were completed for each participant at first (between 5th and 16th weeks) and second visit (between 24th and 28th weeks) during pregnancy (22). At delivery, relevant obstetric and neonate data were obtained from electronic medical records (e.g., mode of delivery, gestational age at delivery, newborn sex, height and weight).

UCP sampling

Venous UCP samples were collected in tubes containing EDTA within the first 30 min after the birth by trained research staff following a standardized procedure. After centrifugation (2,500 g for 10 min at 4˚C), UCP samples were aliquoted and stored at -80˚C until biochemical measurements and RNA extraction. UCP leptin concentrations were quantified by multiplexed particle-based flow cytometric assay (Luminex technology; EMD Millipore, Burlington, MA, USA catalog #HMH3-34k) according to the manufacturer's protocol.

RNA extraction

The total RNA was extracted from up to 500 µl of UCP samples (130 to 500 µl, based on availability) with the MirVana PARIS Kit (ThermoFisher Scientific, catalog #AM1556) according to the manufacturer's recommendations. Samples were first randomized, and hemolysis status was assessed using a visual analogue scale ranging from 0 (non-hemolytic) to 5 (highly hemolytic), as recommended by the Center for Disease Control (23). To control for the plasma volume used for RNA extraction and for technical variability in candidate miRNA quantification (for dPCR analyses only), and as no universally accepted endogenous control with stable expression across individuals and disease states has been identified in human plasma (24-28), Caenorhabditis elegans (cel)-miR-39-3p mimics spike-in control (Qiagen, ON, Canada, catalog # 219610) were added as recommended to the UCP samples at the beginning of the RNA extraction procedure, at a concentration of 6.4x105 copies/µl of plasma. Total RNA was eluted in 75 µl of nuclease-free water. For NGS analysis, 70 µl of RNA eluate was precipitated, dried and resuspended in 5 µl of nuclease-free water, as previously described (29). The entire volume was used for library preparation. A second RNA extraction was performed for dPCR analyses following the same protocol.

Libraries preparation and quality control

Sequencing libraries were prepared using the TruSeq Small RNA Sample Prep kit (Illumina, BC, Canada, catalog #RS-200-0012) following the standard procedure adapted for low RNA concentration samples by Burgos et al (30). After randomization, the process included RNA ligation with adaptors at the 3' and 5' ends, reverse transcription, indexing (1-48, one index per sample) and PCR amplification (15 cycles), using half the reagents volumes to maintain optimal reaction ratios. Libraries were then separated on Novex polyacrylamide TBE Gel, 6% (ThermoFisher Scientific, Waltham, United States, catalog #EC265BOX) and fragments between 145 and 160 bp (including adaptors) were excised and purified. Libraries were eluted in 300 µl of nuclease-free water and incubated overnight at room temperature at 500 RPM on a microplate shaker (VWR, ON, Canada, catalog #12620-930). Lastly, libraries were concentrated by precipitation following the manufacturer's instructions, including a 30 min incubation at -80˚C, and resuspended in 25 µl of 10 mM Tris-HCl pH 8.5 buffer.

Sequencing analysis

Sequencing was performed at the Centre d'expertise et de services Génome Québec (Montréal, Canada) using the NovaSeq 6000 platform. Library quantification was carried out with the Kapa Illumina GA with Revised Primers-SYBR Fast Universal kit (Kapa Biosystems, Wilmington, United States). Libraries were normalized, pooled equimolarly (48 libraries with distinct indexes per lane at a molarity of 200 pM), denatured, and clustered on an Illumina NovaSeq S1 flow cell following the Xp protocol according to the manufacturer's recommendations. Sequencing runs were conducted for 100 cycles in single-end mode.

Bioinformatics processing of the sequencing data

Bioinformatics processing of miRNA sequencing data was performed using the extracellular RNA processing toolkit (exceRpt) pipeline (version 4.6) (31). Adapter sequences and low-quality reads (Phred score <20 for ≥80% of the read) were removed using FASTX-Toolkit and FastQC, integrated within the exceRpt pipeline. High-quality sequences were then aligned to the human genome (GRCh37) and miRBase (version 21) (32) using the STAR alignment algorithm (version 2.4.2a) (33). Post-sequencing quality control identified two samples with no detectable RNA, which were excluded from further analysis.

Validation with digital PCR

Digital PCR was used to validate the sequencing results for two miRNAs, providing absolute quantification (copies/µl of plasma). RNA samples were randomized, and reverse transcription (RT) reaction was performed using the miRCURY LNA RT Kit (QIAGEN, ON, Canada, cat #339340) following the manufacturer's procedure for detection of miRNAs from biofluids, adapted for low-abundance plasma miRNAs. Each RT reaction contained 12.0 µl of total RNA (corresponding to 80 µl of original plasma), 2.0 µl of nuclease-free water, 4.0 µl of 5x miRCURY RT SYBR Green Reaction Buffer and 2.0 µl 10x miRCURY RT Enzyme, for a total volume of 20 µl. Reverse transcription was performed for 60 min at 42 ˚C, followed by enzyme inactivation for 5 min at 95 ˚C. Complementary DNA (cDNA) was diluted 1:2 in nuclease-free water and amplified using the miRCURY LNA miRNA PCR Assays (Qiagen, ON, Canada, cat #339306; assay ID #YP02119311 (hsa-miR-155-5p), #YP00206030 (hsa-miR-183-5p) and #YP00203952 (cel-miR-39-3p)) (Table SI for manufacturer informations) and QIAcuity EG PCR Kit (Qiagen, ON, Canada, cat # 250112), on a Qiacuity One 5-Plex instrument (QIAGEN). The primers sequences used in these pre-designed miRCURY LNA™ miRNA PCR assays are proprietary and not disclosed by the manufacturer (QIAGEN). Briefly, 3 µl of diluted cDNA was combined with 3.8 µl of nuclease-free water, 1.2 µl of 10x miRCURY LNA miRNA PCR Assay and 4 µl of 3x EvaGreen PCR master mix, and 11.0 µl of this reaction were loaded in a QIAcuity 8.5k-partition 96-well nanoplate (Qiagen, ON, Canada, cat # 250021). Amplification was performed for 40 cycles with an annealing/extension temperature of 60˚C (Table SII), and imaging was conducted using the green channel with an exposure time of 200 ms. Thresholds were manually set above the highest signal in the no-template control wells and applied uniformly across all samples. Analyses were performed in duplicate, and data from both replicates were combined using the Hyperwell function, yielding approximately 17k partitions per sample. Digital PCR data were excluded for 13 samples, including 12 from the same RNA extraction batch, due to poor quality, resulting in a final dataset of 201 participants.

The concentration of target miRNAs in UCP (copies/µl) was calculated as described in the Qiagen's protocol: concentration of miRNA measured in the dPCR reaction mix (copies/µl) * [(total dPCR reaction volume(µl))/(volume of cDNA added to dPCR reaction (µl))] * dilution factor of cDNA * [(total RT reaction volume (µl))/volume of plasma equivalent added to RT reaction (µl))]; where the volume of plasma equivalent (µl)=volume of RNA in the RT reaction (µl) * [(volume of plasma used for RNA extraction (µl))/(RNA elution volume (µl))]. For statistical analysis, the miRNA concentrations were normalized for the concentration of the cel-miR-39-3p spike-in control as follows: (concentration of target miRNA/concentration of spike-in) * 10,000 (arbitrary number added to avoid low values) and expressed in relative units.

Statistical analyses

The miRNA-sequencing data analysis was performed using 273 biological replicates (independent UCP samples) from mother-child dyads of the Gen3G cohort for which a complete dataset was available for UCP leptinemia and potential covariates. From these, 201 samples were included to the dPCR analysis of the validation step. Each biological replicate was quantified by dPCR in duplicate (two technical replicates) and their data were combined prior to statistical analysis. The characteristics of the mother-child dyads were presented as mean ± standard deviation and a range of minimum to maximum values for continuous variables or number (percentage) for categorical variables in Table I (miRNA-sequencing analysis sample set; n=273) and Table SIII (dPCR analysis sample set; n=201). Data distribution was assessed using the Shapiro-Wilk test and visual inspection of histograms and Q-Q plots. Because UCP leptin concentrations were not normally distributed, values were natural log-transformed prior to analysis. The mathematical transformation of non-normally distributed variables, using transformation like logarithm or square root, is a commonly used strategy to normalize the distribution of a variable, in order to satisfy the assumption of normality needed to apply parametric tests, which are more powerful than non-parametric tests (34). The normalized concentrations of the two candidate miRNAs selected for the validation step, measured by NGS and dPCR, were also natural log-transformed to meet the normality assumptions required for partial Pearson correlation.

Table I Characteristics of the study participants from the Gen3G cohort (n=273).

Table I

Characteristics of the study participants from the Gen3G cohort (n=273).

Characteristics	Values	Range
Mothers
Age at first trimester, years	28.5±4.1	18.0-47.0
No. of pregnancies carried to term	0.7±0.8	0.0-4.0
Gestational age at delivery, weeks	39.3±1.3	34.2-41.5
Mode of delivery, n (%) vaginal birth	225 (82.4%)
Newborn
Sex, n (%) female	129 (47.3%)
Height, cm	51.0±2.3	43.0-58.5
Weight, g	3,419.2±448.2	2,070.0-5,130.0
UCP leptin level, pg/ml	15,946.9±15,479.2	1,110.0-87,394.0

[i] Data are presented as Mean ± SD or n (%). Gen3G, Genetics of Glucose regulation in Gestation and Growth birth cohort; SD, standard deviation; UCP, umbilical cord plasma.

To identify miRNAs associated with leptinemia in UCP (continuous variable), we used the DESeq2 R package (35), which applies a negative binomial generalized linear model to estimate variance and fold change (FC) in sequencing counts. Read counts correspond to the number of sequencing reads mapped to each miRNA. These counts were normalized using the DESeq2 median-of-ratios method, which adjusts for differences in sequencing depth and RNA composition across samples (35). Briefly, all 2,218 miRNAs detected in our samples were included in the analysis. The statistical model was adjusted for batch effect (sequencing run and region), hemolysis degree, maternal age at first trimester, parity, mode of delivery, child's sex and gestational age at delivery. Default parameters were applied, including a Wald test to identify miRNAs associated with leptinemia, a P-value threshold of P<0.05 and a False Discovery Rate (FDR)-adjusted q-value threshold of q≤0.1. Since no statistically significant association was detected at q-value≤0.1, miRNAs were considered nominally associated with leptinemia at P<0.05. A FC represents the variation in miRNA abundance (mean normalized counts) for a one-unit increase in natural log-transformed leptinemia (pg/ml) in UCP.

Stepwise linear regression analysis was then performed to identify miRNAs independently associated with UCP leptinemia. To limit the number of miRNAs included, we selected those associated with leptinemia at P<0.01 and with a mean normalized read counts >20. The model also included covariates known to influence circulating miRNA profiles and leptinemia, including maternal age at first trimester, mode of delivery, child's sex, the gestational age at delivery and birth weight and height. The full stepwise linear regression model is presented in Table SIV.

Biological pathway analysis was conducted using the DIANA-miRPath v3.0 software (DIANA Tools, University of Thessaly, Larissa, Greece) (36). Kyoto Encyclopedia of Genes and Genomes (KEGG) functional annotations and pathway union analysis mode were applied to identify signaling pathways enriched in genes targeted by the most significant miRNAs associated with UCP leptinemia, with DESeq2 (P<0.01) and with a mean normalized read counts >20. Default settings were applied, including the FDR-adjusted q-value threshold of <0.05. The Tarbase 7.0 database (https://dianalab.e-ce.uth.gr/html/universe/index.php?r=tarbase/index) was selected as the source of experimentally validated miRNA:mRNA interactions (37). We also queried the DIANA microT 2023 webserver (https://dianalab.e-ce.uth.gr/microt_webserver/#/), which provides miRNA target predictions based on in silico microT-CDS algorithm (38).

Two candidate miRNAs were selected for validation by dPCR based on the strength of their association with leptinemia, their abundance in UCP and their independent contribution to explaining variability in leptinemia in UCP. Partial Pearson correlations, adjusted for sequencing batch effect (run and region; NGS data only) hemolysis degree, maternal age at first trimester, parity, mode of delivery, child's sex and gestational age at birth, were calculated to confirm associations between natural log-transformed normalized miRNA concentrations measured by NGS and natural log-transformed leptinemia in UCP. The same model was applied to test associations between natural log-transformed normalized miRNA concentrations measured by dPCR and natural log-transformed leptinemia in UCP. Correlations were considered statistically significant at P<0.05. All analyses were performed using R software (version 4.3.2; R Foundation for Statistical Computing, Vienna, Austria) within RStudio (version 2023.12.0+369; Posit Software, Boston, MA, USA).

Results

Participant's characteristics

The characteristics of the 273 mother-child pairs from the Gen3G prospective pregnancy and early life cohort selected for this study are shown in Table I. Briefly, mothers were on average aged 28.5±4.1 years old [mean ± standard deviation (SD)] at the time of enrolment in the study (1st trimester of pregnancy). The mean gestational age at birth was 39.3±1.3 weeks, and most were born by vaginal delivery (82.4%). The newborns (47.3% girls) had a mean height of 51.0±2.3 cm and a mean weight of 3419±448 g. The characteristics of the 201 participants included in the dPCR analyses were shown in Table SIII and were similar to the characteristics of the 273 participants.

Association between miRNA levels measured by NGS and leptinemia in UCP

A total of 73 miRNAs (P<0.05) were nominally associated with leptinemia in UCP, out of the 2,218 miRNAs detected (Fig. 1). Among these, 47 showed positive associations and 26 showed negative associations with UCP leptinemia. Table II shows the 16 most significant (0.01>P>0.0002) and abundant (normalized read counts >20) miRNAs identified by NGS, along with their corresponding fold change (FC; representing the variation in normalized miRNA counts for a one-unit increase in natural log-transformed leptin concentration) and their p-value. The complete list of 73 miRNAs nominally associated with leptinemia in UCP is shown in Table SV. None of these associations reached statistical significance after FDR adjustment for multiple testing.

Figure 1 Circulating miRNAs associated with leptinemia in UCP. The volcano plot shows the differentially expressed (47 positively and 26 negatively) miRNAs according to leptinemia in UCP. Each point represents a single miRNA. The horizontal dotted line represents the P<0.05 threshold. The model was adjusted for the sequencing batch effect (run and region), the hemolysis degree of UCP samples, the maternal age at the first trimester, the number of pregnancies carried to term, the mode of delivery, the child's sex and the gestational age. Fold-change represents the change in miRNA abundance (NGS data normalized using DESeq2 median-of-ratios) per one unit increase in ln transformed leptinemia (pg/ml). NS, non-significant; miRNA, microRNA; UCP, umbilical cord plasma; NGS, next generation sequencing.

Table II Top 16 miRNAs identified by NGS significantly associated with leptin levels in umbilical cord plasma (P<0.01) and a mean normalized counts >20 determined by DESeq2.

Table II

Top 16 miRNAs identified by NGS significantly associated with leptin levels in umbilical cord plasma (P<0.01) and a mean normalized counts >20 determined by DESeq2.

miRNAs	Normalized read counts (mean ± SD)	Fold-change	P-value
hsa-miR-3909	123.7±114.2	1.14	0.0002
hsa-miR-183-5p	251.4±316.9	1.16	0.0007
hsa-miR-185-5p	1,713.0±2,848.4	1.21	0.0008
hsa-miR-34c-5p	24.7±177.7	1.28	0.0008
hsa-miR-515-3p	23.3±46.0	1.25	0.0010
hsa-miR-518a-5p|hsa-miR-527	74.1±256.6	1.27	0.0010
hsa-miR-424-3p	2,839.9±1,601.2	1.10	0.0030
hsa-miR-182-5p	5,638.9±7,408.9	1.12	0.0040
hsa-miR-561-5p	132.0±305.9	1.15	0.0040
hsa-miR-1294	391.9±533.9	1.15	0.0040
hsa-miR-155-5p	876.9±438.9	0.91	0.0050
hsa-miR-192-5p	97,037.7±100,649.3	1.11	0.0060
hsa-miR-365a-3p|hsa-miR-365b-3p	35.3±58.5	1.16	0.0060
hsa-miR-339-3p	635.9±234.8	0.94	0.0080
hsa-miR-3187-3p	20.7±13.7	1.13	0.0100
hsa-miR-371a-5p	61.5±130.0	1.22	0.0100

[i] Associations between miRNA abundance and leptin levels in UCP were assessed using the Wald test integrated into the DESeq2 R package. The model was adjusted for the sequencing batch effect (run and region), the hemolysis degree of UCP samples, the maternal age at the first trimester, the number of pregnancies carried to term, the mode of delivery, the child's sex and the gestational age. Normalized read counts correspond to sequencing counts normalized using the DESeq2 median-of-ratios method. Effect sizes of the associations are represented by the fold-change values; these are defined as the multiplicative change in normalized miRNA counts associated with a one-unit increase in natural log-transformed leptin levels. Values >1 indicate a positive association, while values <1 indicate a negative association. NGS, Next Generation Sequencing; SD, Standard deviation; miR, microRNA.

Identification of variables independently associated with leptinemia in UCP

To identify variables that are independently associated with leptinemia in UCP and contribute to its variability, we performed stepwise linear regression including the 16 most relevant miRNAs (mean of normalized read counts >20 and P<0.01) along with maternal and newborn's characteristics. The analysis revealed that hsa-miR-515-3p, hsa-miR-155-5p, hsa-miR-424-3p, hsa-miR-183-5p, hsa-miR-339-3p and hsa-miR-3187-3p, together with newborn's weight and sex, explained 35.2% of the variability of UCP leptinemia (p=2.2x10-16; Table III).

Table III Independent predictors of leptin levels in UCP identified by stepwise regression.

Table III

Independent predictors of leptin levels in UCP identified by stepwise regression.

Variables	Estimate	Standard error	T-value	P-value
Weight at birth, g	8.55x10-4	9.85x10-5	8.68	4.06x10-16
Sex of the child (reference: Male)	-6.82x10-1	9.00x10-2	-7.58	5.76x10-13
hsa-miR-515-3p (normalized read counts)	3.03x10-3	1.18x10-3	2.56	1.10x10-2
hsa-miR-155-5p (normalized read counts)	-3.11x10-4	1.19x10-4	-2.62	9.40x10-3
hsa-miR-424-3p (normalized read counts)	-7.37x10-5	4.39x10-5	-1.68	9.40x10-2
hsa-miR-183-5p (normalized read counts)	4.15x10-4	1.77x10-4	2.35	2.00x10-2
hsa-miR-339-3p (normalized read counts)	-7.16x10-4	1.99x10-4	-3.59	3.90x10-4
hsa-miR-3187-3p (normalized read counts)	7.17x10-3	3.40x10-3	2.11	3.60x10-2

[i] A stepwise linear regression analysis was applied to identify miRNAs independently associated with leptin levels in UCP. The dependent variable was natural log-transformed leptinemia. The model included the 16 most relevant miRNAs (mean normalized read counts (sequencing counts normalized using the DESeq2 median-of-ratios method) >20 and P<0.01) associated with leptinemia in UCP, as well as maternal age in the first trimester, mode of delivery, child's sex (male/female), weight (g), height (cm), and gestational age (weeks) at birth. Only variables with P<0.10 were retained in the final model. Regression coefficients (estimates) represent the expected change in natural log-transformed leptin levels associated with a one-unit increase in each predictor variable while holding all other variables constant. Statistical significance (P-value) of each individual regression coefficient within the multivariable linear model was assessed by a one-sample t-test. Model statistics: Residual standard error =0.72 (df=264); Multiple R²=0.3708; Adjusted R²=0.3517; F-statistic =19.45 (df=8, 264), P<2.2x10^-16.

Biological pathway analysis

To better understand the potential pathophysiological roles of the 16 most relevant miRNAs associated with UCP leptinemia, we performed a KEGG pathway analysis. A total of 18 pathways were significantly enriched with the mRNA targets of these miRNAs (P<0.05; Fig. 2). These include the fatty acid biosynthesis (P<1x10-325; hsa-miR-34c-5p and hsa-miR-182-5p), fatty acid metabolism (P=1.21x10-7; hsa-miR-34c-5p, hsa-miR-1294 and hsa-miR-155-5p), and biosynthesis of unsaturated fatty acids (P=4.13x10-4; hsa-miR-1294, hsa-miR-155-5p and hsa-miR-365a-3p) pathways, among others. Based on the DIANA microT 2023 webserver, one of our miRNAs of interest, hsa-miR-155-5p, is predicted to target three miRNA response elements (MRE) in the 3'UTR region of the longest leptin receptor (LEPR) isoform's mRNA.

Figure 2 KEGG pathways targeted by circulating miRNAs associated with leptinemia in UCP. This bar graph shows the KEGG pathways enriched with targets of miRNAs associated with leptinemia in UCP. Each bar represents one pathway, ranked according to its level of significance (FDR adjusted q-value). The number of miRNAs involved in the pathway is shown on the bar. ECM, extracellular matrix; FDR, false discovery rate; KEGG, Kyoto Encyclopedia of Genes and Genomes; UCP, umbilical cord plasma; miRNA, microRNA.

Validation of two candidate miRNAs associated with leptinemia using digital PCR

Hsa-miR-155-5p and hsa-miR-183-5p were selected for dPCR quantification based on the strength of their association with leptinemia (DESeq2: log2(FC)=-0.14, P=0.005 and log2(FC)=0.22, P=0.0007, respectively), their relatively high mean normalized read counts (876.86±438.86 and 251.40±316.94, respectively) and their significant and independent contribution to explain variability in UCP leptinemia in stepwise linear regression analysis. Partial Pearson correlations confirmed the negative association for hsa-miR-155-5p (Fig. 3A; r=-0.13, P=0.03) and positive associations for hsa-miR-183-5p (Fig. 3B; r=0.20, P=0.0009). According to dPCR results, hsa-miR-155-5p and hsa-miR-183-5p had mean concentrations of 412.47±536.27 copies/µl of plasma and 192.82±612.92 copies/µl of plasma, respectively. Consistent with NGS findings, miRNA concentrations measured by dPCR were negatively correlated with leptinemia for hsa-miR-155-5p (r=-0.12, p (one sided)=0.04) and positively correlated for hsa-miR-183-5p (r=0.07, p (one sided)=0.17), although the latter was weak and not statistically significant (Fig. 3C and D).

Figure 3 Pearson partial correlation between circulating levels of candidate miRNAs measured by NGS and dPCR and leptinemia in UCP. Correlations between (A) the abundance of hsa-miR-155-5p measured by NGS, (B) the abundance of hsa-miR-183-5p measured by NGS, (C) the concentration of hsa-miR-155-5p measured by dPCR and (D) the concentration of hsa-miR-183-5p measured by dPCR and leptinemia in UCP. To meet the normality assumption when applying Pearson correlation, the leptinemia and miRNA abundance variables were natural log-transformed (34). NGS data were normalized using DESeq2 median-of-ratios method and dPCR data were normalized for the measured concentration of the cel-miR-39-3p spike-in control. The model was adjusted for the sequencing batch effect (run and region), the hemolysis degree of UCP samples, the maternal age in the first trimester, the number of pregnancies carried to term, the mode of delivery, the child's sex and the gestational age at birth. dPCR, digital PCR; NGS, Next Generation Sequencing; UCP, umbilical cord plasma; miRNA/miR, microRNA.

Discussion

In this study, we aimed to identify miRNAs potentially involved in the regulation of leptinemia in UCP. We identified a total of 73 miRNAs nominally associated (P<0.05) with leptinemia. Although none remained statistically significant after FDR correction, exploratory analyses were conducted to identify the variables that are independently associated with UCP leptinemia. Stepwise linear regression identified six miRNAs (hsa-miR-515-3p, hsa-miR-155-5p, hsa-miR-424-3p, hsa-miR-183-5p, hsa-miR-339-3p and hsa-miR-3187-3p) independently associated with leptinemia. Among these, hsa-miR-155-5p and hsa-miR-183-5p were selected for validation using dPCR quantification based on their significant association with leptinemia, independent contribution to its variability in UCP and abundance in UCP. The biological relevance of these candidates is supported by their known roles in adipogenesis. MiR-155 has been identified as a regulator of key adipogenic genes such as C/EBPβ and PPARγ (39,40), while miR-183 promotes adipogenesis by interfering with the Wnt/β catenin pathway (41). Although dPCR quantification results were generally consistent with NGS findings for both miRNAs, significant association with UCP leptinemia was fully confirmed only for hsa-miR-155-5p. As expected for complex traits and considering the typically subtle regulatory effects of miRNAs on multiple targets rather than being strong regulators of individual genes (42), the observed partial Pearson correlation between hsa-miR-155-5p and leptinemia in UCP was of low effect size (r<0.3), suggesting that other uncontrolled factors also contribute to its variance.

Our study reports that hsa-miR-155-5p concentrations are negatively correlated with leptinemia in UCP, which may suggest that higher miRNA concentration could reduce leptin production and secretion or increase its clearance. However, we cannot exclude the possibility that leptin itself regulates hsa-miR-155-5p expression or secretion. Previous studies have reported the potential involvement of miR-155 in inflammation and adiposity regulation (43-45) and metabolic syndrome (46), as well as its association with cholesterol efflux capacity (47) and triglycerides levels (48). Its low expression was also associated with obesity in Type 2 Diabetes (49). Interestingly, animal and in vitro studies suggest that miR-155 acts as a repressor of adipogenesis. Specifically, miR-155 has been shown to inhibit adipogenesis through downregulation of C/EBPβ and CREB genes expression, thereby blocking TNFa-induced adipogenesis in 3T3-L1 cells (50). Furthermore, miR-155, which is highly expressed in preadipocytes and decreases with maturation, was shown to repress the differentiation of brown and beige adipocytes through a double-negative feedback loop with C/EBPβ (40). Importantly, C/EBPβ has been reported as a transcriptional activator of LEP, and C/EBPβ-/- mice exhibit reduced weight, fat mass and serum leptinemia (51,52). Although this hypothesis requires confirmation, hsa-miR-155-5p may indirectly contribute to attenuating leptinemia through its downregulation of C/EBPβ in adipose tissue.

Our pathway analysis with the microT-CDS database revealed that hsa-miR-155-5p is predicted to target the mRNA of the longest leptin receptor (LEPR) isoform, which is primarily expressed in the hypothalamus and plays a critical role in leptin-dependent central regulation of food intake and energy expenditure (53). Interestingly, previous studies have shown that UCP leptinemia is positively associated with adiposity at birth, but this relationship reverses around 3 years of age (54-56), suggesting the existence of a critical period of leptin sensitivity in early life followed by a phase of resistance to its action. This early-life window of higher leptin sensitivity has been demonstrated in animal studies to be essential for programming of the hypothalamic neural circuits involved in the regulation of food intake and energy expenditure, with long-term effects on fat mass deposition and adiposity in mice (13,57,58). In the absence of leptin resistance during early life, a decrease in LEPR level could be accompanied by a parallel decrease in leptin concentrations, potentially due to a positive feedback loop regulating LEP gene expression or a shorter half-life of the hormone, aimed at maintaining free and bioactive leptin levels in UCP (59). Although further mechanistic studies are needed to clarify the regulation of circulating UCP leptin and the potential role of hsa-miR-155-5p in this process through the LEPR gene and its different isoforms, our findings are consistent with existing literature.

The discrepancy between the associations of hsa-miR-183-5p with leptinemia in UCP observed using NGS and dPCR likely reflects platform-specific biases, which explain the need for validation. NGS libraries are subject to sequence- and structure-dependent ligation and sequencing biases, meaning that read counts do not always accurately represent the absolute abundance of individual miRNAs (60). In contrast, dPCR relies on reverse transcription and highly specific primers and probes; therefore, suboptimal assay design or reduced analytical sensitivity near the detection limit can compromise power (61). Additionally, the two approaches use distinct normalization strategies: NGS normalization is based on DESeq2's median-of-ratios method, whereas dPCR normalization relies on an exogenous spike-in control. Finally, the lack of validation for hsa-miR-183-5p may partly result from differences in the detection of miRNA isoforms (isomiRs) between these technologies. Previous studies have reported that the strength of correlation between NGS and qPCR data varies depending on the isomiR (62,63). Further investigations incorporating detailed isomiR characterization are needed to clarify the role of this miRNA in regulating UCP leptinemia.

Other studies have reported the involvement of circulating miRNAs in the regulation of metabolic homeostasis (64), obesity (18,64-67), adiposity programming (68), and even leptin production and sensitivity (20,69). Our work advances this literature by shifting the focus to the perinatal window and providing the first evidence that circulating miRNAs are potentially associated with leptinemia in UCP at birth. Thus, compared with previous reports, our findings highlight a previously unexplored potential regulatory role of miRNAs in the regulation of leptin levels at birth that may contribute to adiposity programming.

Strengths and limitations

To the best of our knowledge, this is the largest study investigating the circulating miRNAs associated with leptinemia in UCP using NGS. Among its strengths, we included 273 participants from the prospective longitudinal Gen3G cohort, which allowed consideration of potential confounding factors during pregnancy and at birth. Another strength is the use of NGS to comprehensively profile and quantify plasma miRNAs. NGS is a robust technology for miRNA quantification, enabling the identification and measurement of all miRNAs in a sample with high sensitivity and accuracy (70). Furthermore, we validated the NGS results using dPCR, which confirmed the association between leptin and hsa-miR-155-5p levels and provided absolute concentration estimates in UCP.

Our study also has some limitations. Firstly, the sample size remains relatively modest and all participants were of European descent, limiting the generalizability of our results to other population. Secondly, the cross-sectional design precludes establishing causality or the directionality of the association between circulating hsa-miR-155-5p and leptinemia. Thirdly, there is no established gold standard strategy for miRNA data normalization in plasma and no universally accepted endogenous control with stable expression across individuals and disease states (24,25), including U6 small nuclear RNA which may display interindividual variability in certain conditions (71-73). Accordingly, the miRNA quantification data measured by dPCR in this study were normalized for the exogenous cel-miR-39-3p spike-in control, which controls for the plasma volume used for RNA extraction and for technical variability in miRNA extraction and quantification. Finally, although in silico prediction of LEPR as a target of hsa-miR-155-5p is useful for hypothesis generation, it remains a limitation. Further experiments are needed to confirm the regulatory role of hsa-miR-155-5p on the LEPR gene. For example, an interesting perspective could be to perform an in vitro experiment using adipocytes treated with hsa-miR-155-5p agomiRs and antagomiRs to assess its effect on LEPR mRNA and protein levels, as well as luciferase reporter assays to validate predicted MREs (74).

In conclusion, using NGS, we identified 73 circulating miRNAs nominally associated with UCP leptinemia at birth. The negative association between hsa-miR-155-5p plasmatic concentration and UCP leptinemia was confirmed using dPCR quantification. Our biological pathways analysis also indicated that hsa-miR-155-5p is predicted to directly target the LEPR mRNA, potentially influencing leptinemia in UCP. Additional studies are required to elucidate the functions of hsa-miR-155-5p in LEP/LEPR regulation, the directionality of this relationship, and its broader role in leptin-dependent central regulation of food intake, energy expenditure programming, adipogenesis and the development of obesity.

Supplementary Material

Manufacturer information on the miRCURY LNA™ miRNA PCR assays used to quantify candidate miRNAs by digital PCR.

Thermocycling conditions used for digital PCR quantifications.

Characteristics of the study participants from the Gen3G cohort selected for the dPCR validation analysis (n=201).

Variables included in the stepwise linear regression model.

73 miRNAs identified as nominally associated with UCP leptinemia with DESeq2.

Acknowledgements

We would like to thank the Génome Québec Centre of Expertise and Services (CHU Sainte-Justine, Montréal, Canada) for the next generation sequencing service performed as part of this study.

Funding

Funding: LG was supported by a master research award from Fonds de la recherche du Québec en santé (FRQS). KT and AAC were supported by a doctoral research award from FRQS. LB is a Senior Research Scholar from the FRQS and member of the FRQS-Centre Hospitalier Universitaire de Sherbrooke Research Center. Gen3G was supported over the years by an FRQS operating grant (grant no. 20697); Canadian Institutes of Health Research (CIHR) operating grants (grant nos. MOP 115071, PJT-152989, IGH-155183 and PJT-190076) and PP Diabète Québec and La Fondation de ma vie. This study was supported by the Réseau de recherche en santé cardiométabolique, diabète et obésité of the FRQS in collaboration with the Fondation des étoiles.

Availability of data and materials

The raw sequencing data generated in the present study are accessible through the public platform Gene Expression Omnibus database under accession number GSE334144 or at the following URL: https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE334144. The additional data generated in the present study may be requested from the corresponding author.

Authors' contributions

LG performed the digital PCR data collection, performed the statistical analyses, and wrote the manuscript. VD contributed to data interpretation and critically revised the manuscript. LB contributed to study conception and design, data interpretation, manuscript revision, and supervised the project. KT and VD performed the NGS libraries preparation, coordination of sequencing services, as well as the management of the sequencing data generated by the service provider. KT and AAC confirmed the authenticity of all the raw data. AAC contributed to bioinformatics analysis. PP and MFH participated in the conception of the study design and thoroughly revised the manuscript. All authors read and approved the final version of the manuscript.

Ethics approval and consent to participate

Our study was approved by the ethical committee from the Centre Hospitalier Universitaire de Sherbrooke (approval no. 2010-198, 07-027-A1). Written informed consent was obtained from all participants at enrollment in the Gen3G prospective pregnancy and early-life cohort.

Patient consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Use of artificial intelligence tools

During the preparation of this work, artificial intelligence tools were used to improve the readability and language of the manuscript, and subsequently, the authors revised and edited the content produced by the artificial intelligence tools as necessary, taking full responsibility for the ultimate content of the present manuscript.

References