Introduction
Uterine fibroids (UFs), or leiomyomas, are among the most common benign tumors of the female reproductive system, with prevalence rates reaching up to 68% among women of reproductive age.1,2 This condition is often associated with symptoms such as heavy menstrual bleeding and pelvic pain, significantly affecting reproductive health and overall quality of life.3 Although UFs affect a substantial proportion of women, the etiological pathways and mechanisms of disease progression remain incompletely understood, warranting further investigation into their pathogenic drivers.4
In recent years, genome-wide association studies (GWAS) have become an important tool for identifying single nucleotide polymorphisms (SNPs) associated with human diseases in general and with UF in particular.5–10 These findings have significantly advanced our understanding of the molecular mechanisms underlying UF pathogenesis, paving the way for the development of targeted therapeutic approaches.11 However, the interplay between GWAS-identified genetic loci and environmental risk factors remains insufficiently explored.
Environmental factors, including pelvic inflammatory disease (PID) and induced abortions, play a notable role in UF development and progression.12–15 Growing interest in the interactions between genetic and environmental factors has highlighted the importance of examining their influences in a holistic context.16,17 On one hand, SNPs identified through GWAS may predispose individuals to UF by modulating inflammatory, hormonal, and tissue-related processes.18 On the other hand, environmental factors such as inflammation or myometrial trauma could amplify the activity of specific genetic loci, thereby modifying disease risk.19 Previously, we examined GWAS loci associated with UF risk, selecting loci from the GWAS Catalog. That work was the first to demonstrate that these genetic associations could be modified by environmental factors such as smoking, intake of fresh vegetables and fruits, and reproductive risk factors.19
This study aimed to investigate the associations between newly selected GWAS loci and reproductive risk factors in the development of uterine UF, with a particular focus on how induced abortions and PID modulate the associations between these loci and UF risk, explored here for the first time. Our study examines how environmental factors interact with genetic risk, either strengthening or weakening the effects of GWAS-identified loci on disease susceptibility. In the future, these data will help identify key mechanisms of pathogenesis and potential targets for preventive and therapeutic interventions.
Materials and methods
Study participants
We analyzed the same Central Russian cohort (n = 1,188) as in our previous work,19 consisting of 737 clinically confirmed UF patients and 451 healthy controls. The studies were conducted in accordance with the guidelines of the Declaration of Helsinki (as revised in 2024), local legislation, and institutional requirements. The Ethics Committee of Kursk State Medical University approved all procedures (Protocol №5, May 11, 2021), with written consent obtained from all participants. Inclusion required self-reported Russian ethnicity and Central Russian birth origin. Table 1 summarizes the baseline characteristics of the participants.
Table 1Baseline and clinical characteristics of the study cohort
| Baseline characteristics of the study cohort | UF group (n = 737) | Control group (n = 451) | p-value |
|---|
| Age, Ме [Q1; Q3] | 48 [43; 52] | 51 [43; 59] | <0.001 |
| Smoking | Yes, N (%) | 108 (14.7%) | 44 (9.7%) | >0.05 |
| No, N (%) | 629 (85.3%) | 407 (90.3%) | |
| Low fruit/vegetable consumption | Yes, N (%) | 603 (81.8%) | ND | – |
| No, N (%) | 134 (18.2%) | ND | |
| Infertility history | Yes, N (%) | 8 (1.1%) | 2 (0.4%) | >0.05 |
| No, N (%) | 566 (76.8%) | 376 (83.4%) | |
| ND, N (%) | 163 (22.1%) | 73 (16.2%) | |
| Pelvic inflammatory diseases (PID) | Yes, N (%) | 108 (14.6%) | 61 (13.5%) | >0.05 |
| No, N (%) | 464 (63.0%) | 318 (70.5%) | |
| ND, N (%) | 165 (22.4%) | 72 (16%) | |
| Family history of UFs | Yes, N (%) | 205 (27.8%) | 36 (8%) | <0.01 |
| No, N (%) | 532 (72.2%) | 415 (92%) | |
| ND, N (%) | – | – | |
| Menarche age (years) | Ме [Q1; Q3] | 12 [12; 14] | 13 [12; 14] | >0.05 |
| Pregnancy history | Yes, N (%) | 622 (84.4%) | 401 (88.9%) | >0.05 |
| No, N (%) | 29 (3.9%) | 13 (2.9%) | |
| ND, N (%) | 86 (11.7%) | 37 (8.2%) | |
| Gravidity (total pregnancies) | Ме [Q1; Q3] | 3 [2; 5] | 3 [2; 4] | >0.05 |
| Parity | Yes, N (%) | 604 (81.9%) | 394 (87.4%) | >0.05 |
| No, N (%) | 44 (6%) | 19 (4.2%) | |
| ND, N (%) | 89 (12.1%) | 38 (8.4%) | |
| Number of deliveries | Ме [Q1; Q3] | 2 [1; 2] | 2 [1; 2] | >0.05 |
| Medical abortion history | Yes, N (%) | 438 (59.4%) | 225 (49.9%) | >0.05 |
| No, N (%) | 201 (27.3%) | 178 (39.5%) | |
| ND, N (%) | 98 (13.3%) | 48 (10.6%) | |
| Medical abortion numbers | Ме [Q1; Q3] | 1 [0; 2] | 1 [0; 2] | >0.05 |
| Miscarriages in anamnesis | Yes, N (%) | 133 (18%) | 82 (18.2%) | >0.05 |
| No, N (%) | 478 (64.9%) | 318 (70.5%) | |
| ND, N (%) | 126 (17.1%) | 51 (11.3%) | |
| Number of miscarriages | Ме [Q1; Q3] | 0 [0; 0] | 0 [0; 0] | >0.05 |
| Periods prolongation (>7 days) | Yes, N (%) | 180 (24.4%) | 83 (18.4%) | >0.05 |
| No, N (%) | 442 (60%) | 305 (67.6%) | |
| ND, N (%) | 115 (15.6%) | 63 (14%) | |
| Periods regularity | Yes, N (%) | 357 (48.4%) | 200 (44.3%) | <0.001 |
| No, N (%) | 159 (21.6%) | 2 (0.5%) | |
| ND, N (%) | 221 (30%) | 249 (55.2%) | |
| Dysmenorrhea | Yes, N (%) | 227 (30.8%) | 23 (5.1%) | <0.001 |
| No, N (%) | 236 (32%) | 179 (39.7%) | |
| ND, N (%) | 274 (37.2%) | 249 (55.2%) | |
| Menorrhagia (heavy menstrual bleeding) | Yes, N (%) | 370 (50.2%) | 36 (8%) | <0.001 |
| No, N (%) | 144 (19.5 %) | 166 (36.8%) | |
| ND, N (%) | 223 (30.3%) | 249 (55.2%) | |
| Age at uterine fibroid diagnosis (years) | Ме [Q1; Q3] | 40 [36; 45] | – | – |
| Uterine dimensions at time of diagnosis (gestational week equivalents) | Ме [Q1; Q3] | 7 [6; 9] | – | – |
| Current uterine size (weeks of gestation) | Ме [Q1; Q3] | 9 [7; 12] | – | – |
| Growth of UFs (weeks of gestation/year) | Ме [Q1; Q3] | 0.33 [0; 0.76] | – | – |
Case participants with ultrasound-verified UF were enrolled from 2021 to 2023 at two tertiary care facilities: the Perinatal Centre and Kursk City Maternity Hospital. The control cohort, comprising individuals with no signs of UFs on clinical or ultrasound evaluation, was assembled through systematic screening during preventive health visits at regional medical centers and occupational health settings.20,21Figure 1 presents a comprehensive flowchart detailing the participant selection process and applied research methodology.
Selection of environment-associated risk factors for UFs
The following environmental risk factors for UF development were analyzed:
Medical (induced) abortions: These procedures may contribute to UF development through mechanisms such as the absence of hormonal regulation typically seen in late pregnancy following early estrogen surges, tissue repair processes resembling keloid formation, and cytokine-mediated inflammation. Collectively, these processes activate genetic and signaling pathways that drive abnormal tissue growth and fibroid formation.13,14
PID: Disorders associated with trauma, infection, or inflammatory processes disrupt immune homeostasis by upregulating T-helper cytokines while suppressing regulatory T cell activity. This dysregulation promotes excessive fibrotic tissue formation and proliferation.12
Genes and polymorphisms selection
Genes and polymorphisms were selected using the online Reproductive System Knowledge Portal (RSKP), which aggregates data from meta-analyses of GWAS of reproductive diseases worldwide (https://reproductive.hugeamp.org/ ) with the search query “uterine fibroids”. SNPs with a minor allele frequency < 0.05 were excluded, as were those presenting technical challenges for allele-specific fluorescent probe genotyping due to low GC content, absence of GC clamps, or homopolymeric nucleotide repeats. A total of ten SNPs were included in the genotyping: rs2235529 (WNT4), rs59760198 (DNM3), rs10929757 (GREB1), rs1812266 (LOC105375949), rs9419958 (STN1), rs11031731 (THEM7P, WT1), rs2553772 (LOC105376626), rs641760 (PITPNM2), rs1986649 (FOXO1), and rs66998222 (LOC102723323).
Genetic analysis
Genotyping was conducted at the Laboratory of Genomic Research, Research Institute for Genetic and Molecular Epidemiology, Kursk State Medical University (Kursk, Russia). Participants provided up to 5 mL of venous blood, collected in ethylenediaminetetraacetic acid (EDTA)-coated tubes and stored at –20°C until analysis. Genomic DNA was extracted using standard protocols, including phenol/chloroform extraction and ethanol precipitation, and its purity, quality, and concentration were assessed using a NanoDrop spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA).
SNP genotyping was conducted via allele-specific fluorescent probes using custom-developed polymerase chain reaction (PCR) protocols. Primers were designed using Primer3 software.22 Real-time PCR was performed in a 25 µL reaction volume containing 1.5 U of Hot Start Taq DNA polymerase (Biolabmix, Novosibirsk, Russia), approximately 10 ng of template DNA, and the following reagent concentrations: 0.25 µM of each primer, 0.1 µM of each probe, 250 µM of each dNTP, and varying MgCl2 concentrations (1.5 mM for rs66998222; 3 mM for rs10929757, rs9419958, rs1812266, rs1986649, rs2553772, and rs11031731; 3.5 mM for rs641760 and rs2235529; 4 mM for rs59760198). PCR buffer (1×) included 67 mM Tris-HCl (pH 8.8), 16.6 mM (NH4)2SO4, and 0.01% Tween-20. The thermal cycling protocol consisted of an initial denaturation at 95°C for 10 m, followed by 39 cycles of 92°C for 30 s and annealing/elongation at various temperatures for 1 m (60°C for rs1812266; 61°C for rs9419958; 64°C for rs1986649, rs2553772, rs11031731, rs2235529, and rs59760198; 65°C for rs10929757 and rs641760; 66°C for rs66998222).
For quality control, 10% of samples underwent blinded duplicate genotyping, demonstrating >99% concordance. The rs9419958 (STN1) SNP, which deviated from Hardy–Weinberg equilibrium in controls, was re-genotyped and showed 100% concordance with initial results, confirming data reliability.
Data analysis methods
Statistical analyses were performed using STATISTICA software (version 13.3, Santa Clara, CA, USA). The normality of data distributions was evaluated with the Shapiro–Wilk test. As most quantitative variables deviated from normality, results are presented as medians with interquartile ranges. Continuous variables were analyzed using the Mann–Whitney nonparametric test, while categorical variables were evaluated using Pearson’s χ2 test, incorporating Yates’ continuity correction when appropriate for small sample sizes.
Hardy–Weinberg equilibrium for genotype distributions was assessed using Fisher’s exact test. The SNPStats web-based platform (https://www.snpstats.net ) was used to perform logistic regression analyses examining potential correlations between genotype distributions and UF susceptibility. Analyses followed an additive genetic model and were adjusted for confounding factors, including age and family history of UFs. To account for the potential influence of risk factors on genetic associations, separate analyses were conducted for patients with and without exposure to these factors.
The model-based (MB) multifactor dimensionality reduction (MDR) method was applied to investigate combinations of genotypes at two, three, and four loci, assessing both gene–gene (G×G) interactions and interactions between genotypes and UF risk factors (gene–environment, G×E).23 Risk factors such as medical abortions and PID were included in the G×E interaction analyses. Empirical p-values (pperm) for each model were calculated using permutation testing with 1,000 iterations—the default procedure for simultaneously testing all potential interactions of a given complexity.24 Associations with permutation-adjusted significance (pperm < 0.05) were considered statistically reliable. All regression models were adjusted for patient age and documented family history of UFs.
Statistical computations were conducted in R software (version 3.6.3, R Foundation for Statistical Computing, Vienna, Austria). For each level of interaction, three to four models with the highest Wald statistics and most significant p-values were included in the final analysis. The MB-MDR method also allowed the identification of specific genotype combinations significantly associated with the studied phenotypes (p < 0.05). Analyses were performed using the MB-MDR program compatible with the R software (version 3.6.3).
To investigate the functional implications of the studied SNPs, several bioinformatics tools were utilized (described in detail in our previous research25–27):
GTEx Portal (http://www.gtexportal.org/ ) was used to analyze SNP associations with expression quantitative trait loci (eQTLs) across various tissues, including reproductive organs, adipose tissue, blood, and endocrine glands.
eQTLGen (https://www.eqtlgen.org/ ) provided additional insights into SNP–eQTL relationships, particularly in peripheral blood samples.
HaploReg v4.2 (https://pubs.broadinstitute.org/mammals/haploreg/haploreg.php ) was used to assess SNP positioning in regulatory elements, including DNase hypersensitive regions and histone modifications.
atSNP Function Prediction (http://atsnp.biostat.wisc.edu/search ) evaluated how SNP variants affected transcription factor (TF) binding affinity based on reference and alternative alleles.
Gene Ontology (http://geneontology.org/ ) was used to identify biological processes enriched among TFs associated with the studied SNPs, linking these processes to UF pathogenesis.
RSKP (https://cd.hugeamp.org/ ) integrated findings from genetic association studies, providing additional data on potential relationships between identified SNPs and UF-related clinical manifestations, including abnormal heavy menstrual bleeding and elevated body mass index (BMI).
Results
Association of GWAS loci with UF risk in Russian women
The distribution of SNP genotypes across study cohorts is presented in Table S1. As genetic associations may influence population equilibrium patterns, we first verified Hardy–Weinberg equilibrium in control subjects. All examined polymorphisms conformed to Hardy–Weinberg equilibrium expectations (p > 0.05), with the exception of rs9419958 in the STN1 gene (Table S1). Technical validation through replicate genotyping confirmed 100% concordance for rs9419958, warranting its inclusion in subsequent analyses.
Association analysis revealed statistically significant relationships with UF risk (Table 2). The A allele of rs66998222 (LOC102723323) demonstrated a protective effect (odds ratio (OR) = 0.81, 95% confidence interval (CI) = 0.66–0.99, p = 0.038). Conversely, rs11031731 (THEM7P, WT1) was associated with increased disease susceptibility (effect allele A, OR = 1.39, 95% CI = 1.08–1.80, p = 0.01).
Table 2Results of the analysis of associations between GWAS SNPs and UFs in the entire group
| Genetic variant | Effect allele | Other allele | N | OR [95% CI]a | pb |
|---|
| Entire group |
| rs66998222 LOC102723323 | A | G | 1056 | 0.81 [0.66–0.99] | 0.038 |
| rs641760 PITPNM2 | T | C | 1050 | 0.81 [0.65–1.01] | 0.06 |
| rs2553772 LOC105376626 | T | G | 1052 | 0.87 [0.73–1.04] | 0.13 |
| rs10929757 GREB1 | A | C | 1054 | 0.97 [0.81–1.16] | 0.71 |
| rs2235529 WNT4 | T | C | 1047 | 1.06 [0.82–1.36] | 0.67 |
| rs59760198 DNM3 | T | C | 1037 | 1.10 [0.91–1.32] | 0.32 |
| rs1812266 LOC105375949 | C | G | 1056 | 1.12 [0.93–1.33] | 0.23 |
| rs1986649 FOXO1 | T | C | 1055 | 1.12 [0.89–1.40] | 0.33 |
| rs9419958 STN1 | T | C | 1056 | 1.26 [1.00–1.60] | 0.05 |
| rs11031731 THEM7P, WT1 | A | G | 1056 | 1.39 [1.08–1.80] | 0.01 |
Gene-gene interactions analysis (MB-MDR, MDR modeling)
Using the MB-MDR approach, we identified eight significant gene–gene interaction models involving GWAS-identified loci associated with UF: four two-locus models, two three-locus models, and two four-locus models (pperm ≤ 0.05) (Table 3). Notably, five polymorphic loci, rs9419958 (STN1), rs10929757 (GREB1), rs66998222 (LOC102723323), rs2553772 (LOC105376626), and rs2235529 (WNT4), were incorporated into two or more of the most statistically significant gene–gene interaction patterns, suggesting their potential role in polygenic disease mechanisms. At the next stage, the interactions of these genetic variants were analyzed using the MDR method (Fig. 2).
Table 3UF-associated gene-gene interactions (MB-MDR modeling)
| Gene-gene interaction models | NH | βH | WH | NL | βL | WL | Wmax | pperm |
|---|
| The best two-locus models of intergenic interactions (for G×G models with pmin < 0.002, 1000 permutations) |
| rs9419958STN1 × rs10929757GREB1 | 2 | 0.1434 | 10.253 | 1 | −0.0800 | 6.156 | 10.253 | 0.028 |
| rs641760 PITPNM2 × rs9419958STN1 | 1 | 0.0730 | 3.121 | 2 | −0.1064 | 10.030 | 10.030 | 0.029 |
| rs2553772LOC105376626 × rs9419958STN1 | 2 | 0.1433 | 10.110 | 1 | −0.1226 | 8.554 | 10.110 | 0.033 |
| rs9419958STN1 × rs66998222LOC102723323 | 1 | 0.1158 | 4.637 | 1 | −0.1108 | 9.856 | 9.856 | 0.033 |
| The best three-locus models of intergenic interactions (for G×G models with pmin < 1×10−4, 1000 permutations) |
| rs2553772LOC105376626 × rs2235529WNT4 × rs9419958STN1 | 4 | 0.2192 | 20.519 | 4 | −0.1415 | 8.741 | 20.52 | 0.002 |
| rs9419958STN1 × rs59760198 DNM3 × rs66998222LOC102723323 | 0 | NA | NA | 4 | −0.1652 | 22.016 | 22.02 | 0.006 |
| The best four-locus models of gene-gene interactions (for G×G models with pmin < 1×10−8, 1000 permutations) |
| rs11031731 THEM7P, WT1 × rs9419958STN1 × rs10929757GREB1 × rs66998222LOC102723323 | 2 | 0.1961 | 6.218 | 6 | −0.2643 | 41.58 | 41.58 | < 0.001 |
| rs2553772LOC105376626 × rs2235529WNT4 × rs10929757GREB1 × rs66998222LOC102723323 | 4 | 0.1868 | 17.256 | 8 | −0.4425 | 40.01 | 40.01 | 0.001 |
The MDR method indicated that the genetic variants included in the best G×G models exhibited multidirectional effects (synergism, antagonism, and additive effects). Alongside interaction analyses, we separately evaluated the individual (mono) effects of the genetic variants comprising the most significant gene–gene interaction models. Their contributions to UF entropy (0.02–0.39%) were comparable with the effects of intergenic interactions (0.01%–0.42%). The most pronounced mono-effect was observed for rs9419958 (STN1; 0.39% entropy), while the strongest interactions occurred between rs2235529 (WNT4) and rs2553772 (LOC105376626; 0.42% entropy; pronounced synergism). Among the strongest associations with UF risk were specific genotype combinations from interacting polymorphic variants, revealing allelic patterns that may influence disease susceptibility: 1) Decrease risk of UF: rs9419958 STN1 C/C × rs10929757 GREB1 C/A (beta = −0.08, p = 0.01), rs641760 PITPNM2 C/T × rs9419958 STN1 C/C (beta = −0.087, p = 0.01), rs9419958 STN1 C/C × rs66998222 LOC102723323 A/G (beta = −0.11, p = 0.002), rs9419958 STN1 C/C × rs59760198 DNM3 C/C × rs66998222 LOC102723323 A/G (beta = −0.198, p = 0.001), rs11031731 THEM7P, WT1 G/G × rs9419958 STN1 C/C × rs10929757 GREB1 C/A × rs66998222 LOC102723323 A/A (beta = −0.45, p = 0.0001), rs2553772 LOC105376626 T/G × rs2235529 WNT4 T/C × rs10929757 GREB1 C/С × rs66998222 LOC102723323 А/G (beta = −0.557, p = 0.003); 2) Increase disease risk: rs2553772 T/T × rs9419958 T/C (beta = 0.186, p = 0.006), rs2553772 LOC105376626 T/T × rs2235529 WNT4 C/C × rs9419958 STN1 T/C (beta = 0.29, p = 0.001) (Table S2).
Environmental risk factor stratification analysis of GWAS SNPs
The investigation extended to interactions with environmental risk factors, with detailed stratified analyses presented in Table S3 and summarized in Table 4. We found that rs66998222 (LOC102723323; OR = 0.70, 95% CI = 0.54–0.91, p = 0.009) and rs641760 (PITPNM2; OR = 0.71, 95% CI = 0.53–0.96, p = 0.026) were associated with a reduced risk of UF in patients with a history of induced abortions. Stratified analysis revealed that rs2553772 (LOC105376626) conferred protection against UF specifically in women without a history of induced abortion (OR = 0.61, 95% CI = 0.45–0.82, p = 0.001). Conversely, rs11031731 (THEM7P/WT1) showed risk-enhancing effects both in women with a history of abortion (OR = 1.60, 95% CI = 1.12–2.28, p = 0.008) and in those without a history of PID (OR = 1.43, 95% CI = 1.04–1.96, p = 0.02).
Table 4Effect modification by prior induced abortions and PID history on GWAS SNP–UF risk associations
| SNP | Effect allele | Other allele | N | OR [95% CI] | p | N | OR [95% CI]a | pb |
|---|
| | | Without induced abortions | Induced abortions |
| rs2553772 LOC105376626 | T | G | 331 | 0.61 (0.45–0.82) | 0.001 | 585 | 1.04 (0.82–1.32) | 0.75 |
| rs66998222 LOC102723323 | A | G | 335 | 0.93 (0.65–1.33) | 0.69 | 585 | 0.70 (0.54–0.91) | 0.009 |
| rs641760 PITPNM2 | T | C | 329 | 1.04 (0.72–1.49) | 0.84 | 585 | 0.71 (0.53–0.96) | 0.026 |
| rs11031731 THEM7P, WT1 | A | G | 335 | 1.27 (0.80–2.01) | 0.3 | 585 | 1.60 (1.12–2.28) | 0.008 |
| | | Without PID | With PID |
| rs11031731 THEM7P, WT1 | A | G | 696 | 1.43 (1.04–1.96) | 0.02 | 136 | 2.01 (0.91–4.43) | 0.07 |
Gene-environment interactions of UF GWAS genes (MB-MDR, MDR modeling)
MB-MDR analysis identified five significant gene-environment interaction models: two two-level, two three-level, and one four-level interaction (Table 5). Remarkably, SNPs rs2553772 (LOC105376626) and rs1986649 (FOXO1), along with induced abortions, participated in two or more of the best gene–environment interaction models, and the interaction of induced abortions × rs2553772 (LOC105376626) contributed to three of these five top G×E interactions. At the next stage, the interactions of these genetic variants and risk factors were analyzed using the MDR method (Fig. 3).
Table 5UF-associated gene–environment interactions (MB-MDR modeling)
| Gene-environment interaction models | NH | βH | WH | NL | βL | WL | Wmax | pperm |
|---|
| The best two-order models of gene-smoking interactions (for G×E models with pmin < 1×10−5, 1000 permutations) |
| Ind_abort × rs2553772LOC105376626 | 1 | 0.0597 | 3.265 | 2 | −0.1653 | 21.60 | 21.60 | < 0.001 |
| Ind_abort × rs1986649FOXO1 | 1 | 0.0950 | 9.321 | 1 | −0.1671 | 21.48 | 21.48 | < 0.001 |
| The best three-order models of gene-interactions (for G×E models with pmin < 1×10−6, 1000 permutations) |
| Ind_abort × rs11031731 THEM7P, WT1 × rs1986649FOXO1 | 3 | 0.0921 | 7.449 | 1 | −0.222 | 30.39 | 30.39 | < 0.001 |
| Ind_abort × rs2553772LOC105376626 × rs1812266 LOC105375949 | 1 | 0.1230 | 5.179 | 3 | −0.245 | 26.31 | 26.31 | 0.002 |
| The best four-order models of gene-interactions (for G×E models with pmin < 1×10−9, 1000 permutations) |
| Ind_abort × rs2553772LOC105376626 × rs59760198 DNM3 × rs66998222 LOC102723323 | 5 | 0.1649 | 17.794 | 10 | −0.2872 | 40.14 | 40.14 | < 0.001 |
First, the MDR methodology quantified induced abortions as having the strongest individual environmental effect (1.15% entropy contribution), exceeding both individual SNP effects (0.11–0.16%) and G×E interaction effects (0.05–0.34%). Second, induced abortions were characterized by pronounced synergism in interaction with rs2553772 (LOC105376626) and additive (independent) effects in interaction with rs1986649 (FOXO1). Third, SNPs characterizing the best G×E models (rs2553772 LOC105376626 and rs1986649 FOXO1) exhibited moderate antagonism in interaction with each other. Fourth, the following gene–environment interactions showed the strongest associations with reduced UF risk: no induced abortions × rs2553772 LOC105376626 T/T (beta = −0.236, p = 5.54×10−5), no induced abortions × rs11031731 THEM7P, WT1 G/G × rs1986649 FOXO1 C/C (beta = −0.222, p = 4.59×10−8), no induced abortions × rs2553772 LOC105376626 T/T × rs1812266 LOC105375949 C/G (beta = −0.345, p = 3.21×10−5), induced abortions × rs2553772 LOC105376626 G/G× rs59760198 DNM3 C/C × rs66998222 LOC102723323 A/G (beta = −0.21, p = 0.016). The following G×E interaction was most strongly associated with increased UF risk: induced abortions × rs1986649 FOXO1 C/C (beta = 0.095, p = 2.33×10−3) (Table S4).
GWAS loci associated with clinical course parameters of UFs
We also assessed the influence of genotype associations on the quantitative characteristics of the UF clinical course (Fig. 4, Table S5). A pronounced association was observed for rs1986649 (FOXO1) with age of disease onset (p = 0.049) and UF growth rate (p = 0.017), so this SNP was included in further study.
Functional annotation of SNPs
QTL effects
According to the GTEx Portal, the C allele of SNP rs641760 (PITPNM2) is associated with increased expression of CDK2AP1 in subcutaneous adipose tissue and whole blood, RP11-282O18.3 in the uterus and subcutaneous adipose tissue, C12orf65 in subcutaneous adipose tissue, and ARL6IP4 in whole blood (Table 6). Conversely, this allele is linked to decreased expression of KMT5A and ABCB9 in whole blood.
Table 6Cis-eQTL–mediated gene expression modulation by UF-linked GWAS SNPs (GTEx Portal and eQTLGen data)
| Genetic variant | Assessed allele | Expressed gene | p | Effect (NES) | Tissue | Assessed allele | Symbol | Z-score | FDR |
|---|
| GTEx portal | eQTLgene |
| rs641760 PITPNM2 (T/C) | C | CDK2AP1 | 2.3×10−32 | ↑(0.40) | Adipose - Subcutaneous | T | SETD8 | ↑(21.4017) | 0 |
| | C12orf65 | 5.0×10−12 | ↑(0.25) | Adipose - Subcutaneous | | MPHOSPH9 | ↑(15.2575) | 0 |
| | RP11-282O18.3 | 3.2×10−7 | ↑(0.22) | Adipose - Subcutaneous | | SBNO1 | ↓(−8.2818) | 0 |
| | RP11-282O18.3 | 2.4×10−5 | ↑(0.57) | Uterus | | PITPNM2 | ↑(6.0034) | 2×10−5 |
| | CDK2AP1 | 8.1×10−12 | ↑(0.29) | Whole blood | | SNRNP35 | ↑(5.9739) | 2.6×10−5 |
| | KMT5A | 1.9×10−9 | ↓(−0.19) | Whole blood | | OGFOD2 | ↑(5.3297) | 0.0003 |
| | ARL6IP4 | 3.4×10−5 | ↑(0.072) | Whole blood | | EIF2B1 | ↑(4.9115) | 0.0026 |
| | ABCB9 | 4.2×10−5 | ↓(−0.13) | Whole blood | | | | |
| rs2553772 LOC105376626 (T/G) | G | CD44 | 8.8×10−8 | ↑(0.13) | Adipose - Subcutaneous | T | RP1-68D18.4 | ↓(−5.2045) | 0.0006 |
| rs1986649 FOXO1 (C/T) | T | RLIMP1 | 0.0001 | ↑(0.24) | Adipose - Subcutaneous | T | MRPS31 | ↓(−11.9938) | 0 |
| | | | | | | FOXO1 | ↑(6.7807) | 0 |
| | | | | | | SLC25A15 | ↓(−5.9359) | 2.6×10−5 |
| | | | | | | ELF1 | ↓(−4.8444) | 0.0036 |
| | | | | | | MIR621 | ↓(−4.6621) | 0.0083 |
| | | | | | | KBTBD7 | ↓(−4.3822) | 0.0304 |
The eQTLGen Browser reports that the T allele of SNP rs641760 (PITPNM2) increases the expression of SETD8, MPHOSPH9, PITPNM2, SNRNP35, OGFOD2, and EIF2B1, while decreasing expression of SBNO1 in blood (Table 6).
For SNP rs2553772 (LOC105376626), the G allele increases CD44 expression in adipose tissue, while the T allele decreases RP1-68D18.4 expression (Table 6).
Finally, the T allele of SNP rs1986649 (FOXO1) correlates with increased FOXO1 expression and decreased levels of MRPS31, SLC25A15, ELF1, MIR621, and KBTBD7 in blood, as well as reduced RLIMP1 expression in subcutaneous adipose tissue (Table 6).
Transcription factors
Analysis of TFs revealed that the risk allele A of rs66998222 (LOC102723323) creates DNA binding sites for 28 TFs, co-regulating negative regulation of cell population proliferation (GO:0008285; false discovery rate (FDR) = 0.0099) and cellular response to cytokine stimulus (GO:0071345; FDR = 0.0107). The protective allele G of rs66998222 (LOC102723323) creates DNA binding sites for 24 TFs, co-regulating interleukin (IL)-1 beta production (GO:0032651; FDR = 0.0229), response to hypoxia (GO:0001666; FDR = 0.0302), and positive regulation of cytokine production (GO:0001819; FDR = 0.0190) (Table S6).
The risk allele A of rs11031731 (THEM7P, WT1) generates DNA binding sites for 33 TFs involved in the regulation of adipose tissue development (GO:1904177; FDR = 0.0318) and cell population proliferation (GO:0008283; FDR = 0.0478). The protective allele G of rs11031731 (THEM7P, WT1) creates DNA binding sites for 35 TFs, regulating the following biological processes: regulation of cellular response to growth factor stimulus (GO:0090287; FDR = 0.0457), negative regulation of cell population proliferation (GO:0008285; FDR = 0.0230), positive regulation of apoptotic process (GO:0043065; FDR = 0.0311), negative regulation of cell differentiation (GO:0045596; FDR = 0.0130), and regulation of growth (GO:0040008; FDR = 0.0495) (Table S7).
The risk allele G of rs2553772 (LOC105376626) creates DNA binding sites for 20 TFs, jointly involved in the SREBP signaling pathway (GO:0032933; FDR = 0.0045), cellular response to transforming growth factor beta stimulus (GO:0071560; FDR = 0.0301), steroid metabolic process (GO:0008202; FDR = 0.0091), and negative regulation of cell population proliferation (GO:0008285; FDR = 0.0037). The protective allele T of rs2553772 (LOC105376626) generates DNA binding sites for 42 TFs, regulating the following biological processes: positive regulation of vascular endothelial growth factor production (GO:0010575; FDR = 0.0419), positive regulation of angiogenesis (GO:0045766; FDR = 0.0063), angiogenesis (GO:0001525; FDR = 0.0008), muscle organ development (GO:0007515; FDR = 0.0436), and regulation of cell population proliferation (GO:0042127; FDR = 0.0175) (Table S8).
The T allele of rs1986649 (FOXO1) creates DNA binding sites for 30 TFs, jointly involved in the SREBP signaling pathway (GO:0032933; FDR = 0.0067), positive regulation of transforming growth factor beta production (GO:0071636; FDR = 0.0206), regulation of IL-5 production (GO:0032674; FDR = 0.0239), SMAD protein signal transduction (GO:0060395; FDR = 0.0340), positive regulation of IL-4 production (GO:0032753; FDR = 0.0380), cellular response to transforming growth factor beta stimulus (GO:0071560; FDR = 0.0042), and regulation of transforming growth factor beta receptor signaling pathway (GO:0017015; FDR = 0.0427). The reference allele C of rs1986649 (FOXO1) creates DNA binding sites for 57 TFs, co-regulating positive regulation of vitamin D receptor signaling pathway (GO:0070564; FDR = 0.0069), intrinsic apoptotic signaling in response to DNA damage by p53-class mediator (GO:0042771; FDR = 0.0035), negative regulation of canonical Wnt signaling pathway (GO:0090090; FDR = 0.0118), transforming growth factor beta receptor superfamily signaling pathway (GO:0141091; FDR = 0.0216), response to decreased oxygen levels (GO:0036293; FDR = 0.0176), and response to growth factor (GO:0070848; FDR = 0.0214) (Table S9).
Histone modifications
Using HaploReg v4.2, UF-associated SNPs were annotated for regulatory potential, identifying characteristic histone modification patterns in blood cells, vessels, adipose tissue, and ovaries (Table 7).
Table 7The impact of UF-associated GWAS SNPs on histone marks in various tissues
| SNP (Ref/Alt allele) | Tissue marks | Blood | Blood vessels | Ovary | Adipose |
|---|
| rs66998222 LOC102723323 (G/A) | H3K4me1 | – | – | – | – |
| H3K4me3 | – | – | – | – |
| H3K27ac | – | – | – | – |
| H3K9ac | – | – | – | – |
| rs11031731 THEM7P, WT1 (G/A) | H3K4me1 | Enh | – | Enh | – |
| H3K4me3 | – | – | – | – |
| H3K27ac | – | – | – | – |
| H3K9ac | Pro | – | – | – |
| rs641760 PITPNM2 (T/C) | H3K4me1 | Enh | – | – | Enh |
| H3K4me3 | Pro | – | – | – |
| H3K27ac | Enh | – | – | Enh |
| H3K9ac | Pro | – | – | Pro |
| rs2553772 LOC105376626 (T/G) | H3K4me1 | Enh | Enh | – | Enh |
| H3K4me3 | Pro | – | – | – |
| H3K27ac | Enh | – | – | Enh |
| H3K9ac | Pro | – | – | – |
| rs1986649 FOXO1 (C/T) | H3K4me1 | Enh | – | Enh | Enh |
| H3K4me3 | Pro | – | – | – |
| H3K27ac | Enh | – | – | Enh |
| H3K9ac | – | – | – | Pro |
The genomic regions of SNPs rs11031731 (THEM7P, WT1), rs641760 (PITPNM2), rs2553772 (LOC105376626), and rs1986649 (FOXO1) were marked by characteristic histone modifications in blood cells, including mono-methylation (H3K4me1) and tri-methylation (H3K4me3) at the fourth lysine residue of histone H3, along with acetylation at the 27th (H3K27ac) and 9th (H3K9ac) lysine residues of histone H3. Similarly, SNP rs2553772 (LOC105376626) is associated with histone marks in blood vessels; SNPs rs11031731 (THEM7P, WT1) and rs1986649 (FOXO1) are associated with histone marks in ovaries; and SNPs rs641760 (PITPNM2), rs2553772 (LOC105376626), and rs1986649 (FOXO1) are linked with histone marks in adipose tissue (Table 7).
Computational analysis of UF GWAS SNPs and associated clinical manifestations
Cross-referencing with the RSKP bioinformatics resource revealed population-wide consistency for the identified UF-associated SNPs. The risk variants rs11031731 (THEM7P/WT1), rs641760 (PITPNM2), rs2553772 (LOC105376626), and rs1986649 (FOXO1) demonstrated conserved risk-increasing effects across ethnic groups, while rs66998222 (LOC102723323) maintained its protective association. These directional effects extended to both UF susceptibility and heavy menstrual bleeding phenotypes in diverse populations. Additionally, SNPs were linked with the following UF-related phenotypes: decreased age of menarche (rs1986649 FOXO1), increased age of menarche (rs66998222 LOC102723323 and rs641760 PITPNM2), increased age at natural menopause (rs641760 PITPNM2 and rs1986649 FOXO1), decreased age at natural menopause (rs11031731 THEM7P, WT1), increased BMI (rs11031731 THEM7P, WT1), and decreased BMI (rs66998222 LOC102723323 and rs1986649 FOXO1) (Table S10).
Discussion
UF pathogenesis involves a complex interplay between genetic predisposition and environmental factors.4 Examining how these interactions influence disease onset is critical.28,29 Among these, induced abortions and PID may act as modulators, influencing both the development and clinical course of UFs.14,19,30,31 While GWAS have identified numerous susceptibility loci,9,10,8,32 their interaction with environmental exposures remains understudied. Building on our prior work with this Central Russian cohort, where we demonstrated broad modification of UF-associated loci by reproductive factors,19 we now extend this paradigm through functional analysis of ten newly selected SNPs, specifically probing abortion- and inflammation-related modulation mechanisms.
We found that abortions and pelvic inflammation may contribute to UFs through three shared biological pathways: 1) disruption of hormonal balance following pregnancy termination, leading to unopposed estrogenic stimulation33,34; 2) altered myometrial healing, promoting fibroproliferative activity35,36; and 3) increased local inflammation, stimulating cytokine-driven growth pathways (e.g., IL-6, transforming growth factor-β).37 These factors can potentiate or mitigate genetic predispositions, depending on the locus involved (Fig. 5).
We identified a protective association of rs66998222 and rs641760 in women with a history of induced abortion. Prior GWAS linked rs66998222 (LOC102723323) with decreased UF risk,32 and our data confirm this in both the entire group and abortion-positive patients. Functional analysis revealed that the protective allele A at rs66998222 (LOC102723323) enhances binding of TFs involved in negative regulation of proliferation and cytokine responses, particularly relevant in an inflammatory microenvironment (GO:0008285, GO:0071345).13
The rs641760 (PITPNM2) locus, also protective in this subgroup, modulates the expression of genes such as CDK2AP1, ABCB9, and KMT5A, which are involved in progesterone-regulated uterine cell growth and tumor regulation, respectively.38–40 Other target genes include ARL6IP4 and MPHOSPH9, both of which influence estrogen signaling, metabolism, and mTOR activation—processes amplified during hormonal fluctuations post-abortion.41 Bioinformatic analysis via the RSKP portal confirmed the protective effect of the T allele of rs641760 (PITPNM2) and the A allele of rs66998222 (LOC102723323) against UFs and heavy menstrual bleeding, as well as their impact on UF-associated phenotypes such as BMI, age at menarche, and age at menopause.
The rs11031731 locus conferred elevated UF risk in the general cohort and the abortion-positive subgroup, but not in patients with PID. This disparity may reflect the dominant role of chronic inflammation in the PID group, overshadowing the genetic effect.12 Conversely, in patients with a history of induced abortions, the increased UF risk appears to result from the combined effects of inflammation and abrupt changes in steroid sex hormone levels following the procedure. Early termination of pregnancy fails to suppress the tumor-stimulating effects of estrogen during later stages of pregnancy,14 which can subsequently drive UF growth through various genetic and signaling pathways.33,34 Functionally, the risk allele A creates TF binding sites that upregulate pathways related to adipogenesis and proliferation (GO:1904177, GO:0008283). These mechanisms are amplified in patients with higher BMI or estrogen exposure, explaining the heightened risk in specific environmental contexts.42 According to the RSKP portal, this allele increases the risk of not only UF but also elevated BMI, an important UF risk factor.4
This locus demonstrated a protective effect in patients without a history of abortion. This protective effect may stem from the SNP’s role in regulating tissue healing and restoration.35,36 Functional annotation showed that the T allele influences CD44 expression, a gene implicated in tissue remodeling.43 The SNP also promotes TF binding linked to vascular growth, angiogenesis (GO:0010575, GO:0045766, GO:0001525), and myometrial homeostasis (GO:0042127, GO:0007515), potentially buffering against fibrotic changes in hormonally stable or low-inflammation environments. In patients without a history of induced abortions, this allele may facilitate balanced healing processes, contrasting with the more chaotic cellular responses triggered by injury or inflammation.36 This ability to regulate tissue repair and maintain homeostasis reduces UF risk by preventing hyperplastic changes in the myometrium.35
Although not modulated by abortion or PID, the TT genotype of rs1986649 (FOXO1) was associated with earlier onset and accelerated UF growth. FOXO1 has long been studied as a key gene associated with UF development.44–46FOXO1 regulates progesterone (P4) signaling and inflammatory pathways, such as COX-2 and IL-6, linking endocrine and immune mechanisms.11,47 eQTL data indicate that it influences the expression of MRPS31 (p53 interaction),48ELF1 (increased tumor risk),49KBTBD7 (MAPK/AP-1 inflammation),50SLC25A15 (tumor formation under hypoxic conditions),51 and other genes implicated in tumor metabolism and UF risk. The risk allele T upregulates transforming growth factor-β/SMAD (GO:0071636, GO:0060395, GO:0071560, GO:0017015), SREBP (GO:0032933), and inflammatory IL-4/IL-5 (GO:0032674, GO:0032753) pathways, driving fibrosis and proliferation.35,52–54 According to the RSKP portal, the SNP rs1986649 (FOXO1) is associated with higher UF risk, heavy menstrual bleeding, increased age at natural menopause, and earlier age at menarche, extending the duration of reproductive years and thereby elevating UF risk.55
This study has several limitations. First, the limited number of analyzed SNPs may have reduced detectable genetic associations. Second, TaqMan genotyping constraints excluded some SNPs due to probe design issues. Third, the lack of an independent Russian cohort prevented replication analysis. Fourth, the binary classification of induced abortions (present/absent) precluded analysis of procedure frequency effects on genetic risk modulation.