Yingxiao Li1
Haoyue Chen
Hao Zhang
Chuanliang Zhao- 1Department of Gynecology, The Affiliated Taian City Central Hospital of Qingdao University, Tai’an, Shandong, China
- 2Department of Rehabilitation Medical Center, The Affiliated Taian City Central Hospital of Qingdao University, Tai’an, Shandong, China
- 3Hydrogen Medical Research Center, The Affiliated Taian City Central Hospital of Qingdao University, Tai’an, Shandong, China
- 4Department of Orthopedics, The Affiliated Taian City Central Hospital of Qingdao University, Tai’an, Shandong, China
- 5Medical Integration and Practice Center, Shandong University School of Medicine, Shandong University, Jinan, Shandong, China
Background: Oxidative stress has been implicated in the pathogenesis of uterine leiomyoma (ULM) with an increasing incidence. This study aimed to identify potential oxidative stress-related biomarkers in ULM using transcriptome data integrated with Mendelian randomization (MR) analysis.
Methods: Data from GSE64763 and GSE31699 in the Gene Expression Omnibus (GEO) were included in the analysis. Oxidative stress-related genes (OSRGs) were identified, and the intersection of differentially expressed genes (DEGs), Weighted Gene Co-expression Network Analysis (WGCNA) genes, and OSRGs was used to derive differentially expressed oxidative stress-related genes (DE-OSRGs). Biomarkers were subsequently identified via MR analysis, followed by Gene Set Enrichment Analysis (GSEA) and immune infiltration analysis. Nomograms, regulatory networks, and gene-drug interaction networks were constructed based on the identified biomarkers.
Results: A total of 883 DEGs were identified between ULM and control samples, from which 42 DE-OSRGs were screened. MR analysis revealed four biomarkers: ANXA1, CD36, MICB, and PRDX6. Predictive nomograms were generated based on these biomarkers. ANXA1, CD36, and MICB were significantly enriched in chemokine signaling and other pathways. Notably, ANXA1 showed strong associations with follicular helper T cells, resting mast cells, and M0 macrophages. CD36 was positively correlated with resting mast cells, while MICB was negatively correlated with macrophages. Additionally, ANXA1 displayed strong binding energy with amcinonide, and MICB with ribavirin.
Conclusion: This study identified oxidative stress-related biomarkers (ANXA1, CD36, MICB, and PRDX6) in ULM through transcriptomic and MR analysis, providing valuable insights for ULM therapeutic research.
1 Introduction
Uterine leiomyoma (ULM), also referred to as uterine fibroids or myomas, represents the most prevalent benign tumor in the female reproductive system (1). Epidemiological evidence indicates that the incidence of ULM increases with age, with a particularly high prevalence among women of reproductive age. In women over 35 years, the incidence rate ranges from approximately 20% to 40% (2). Several factors, including environmental and psychological influences, contribute to ULM development. Notably, the elevation of estrogen receptor concentrations is widely accepted as a primary driver of ULM formation (3). Furthermore, familial clustering observed in numerous ULM cases suggests a substantial genetic component in its pathogenesis (4). Despite these insights, the precise molecular mechanisms underlying ULM remain incompletely understood, necessitating further investigation into potential pathogenic pathways and biomarkers to improve patient quality of life.
Oxidative stress refers to the disruption of the equilibrium between reactive oxygen species (ROS) production and the body’s endogenous antioxidant defenses (5). While basal levels of ROS are crucial for maintaining cellular homeostasis, excessive ROS can inflict structural and genetic damage, including protein and lipid peroxidation, which can impair cellular function and trigger apoptosis, contributing to disease progression (6). The uterus, within the female reproductive system, is particularly vulnerable to oxidative stress and subsequent DNA damage, which may be pivotal in ULM development (7, 8). Whole genome sequencing has elucidated the genetic foundation of ULM, revealing that over 70% of leiomyoma (LM) cases harbor mutations in the mediator complex subunit 12 (MED12), with mutation rates positively correlating with tumor quantity. ROS is implicated in both promoting MED12 mutations and facilitating tumor growth (9). Research suggests that uterine muscle contractions and vascular changes during the menstrual cycle induce hypoxic microenvironments within the uterine muscle, which subsequently trigger gene expression promoting leiomyoma proliferation under low oxygen conditions (10). Leiomyoma cells exhibit lower expression of antioxidant enzymes, such as superoxide dismutase (SOD) and catalase (CAT), compared to normal uterine muscle cells, particularly in hypoxic environments. ULM is characterized by elevated oxidative stress and insufficient antioxidant defense capacity (11). Therefore, investigating the role of oxidative stress in ULM pathogenesis is essential for advancing diagnostic and therapeutic strategies. Biomarkers of oxidative stress may emerge as valuable diagnostic tools for identifying patients with ULM. While medical intervention is often necessary for ULM treatment, antioxidants have shown potential in both prevention and therapeutic management, though further research is required. In summary, oxidative stress represents a critical and emerging aspect of ULM pathology, playing a central role in its progression. Further experimental and clinical research is needed to elucidate the cellular and molecular mechanisms linking oxidative stress to ULM.
Mendelian randomization (MR) represents a sophisticated approach in epidemiology and related fields for establishing causal relationships, utilizing genetic variation to evaluate the causal influence of potential risk factors on health outcomes. By leveraging the random assortment of alleles, MR effectively mitigates confounding factors that often obscure the exposure-outcome relationship in observational studies, enabling a more reliable and unbiased assessment of exposure effects on various outcomes compared to traditional methodologies (12, 13). MR analysis is contingent upon three fundamental assumptions: relevance (the genetic variant must be strongly associated with the exposure), independence (the genetic variant should be independent of confounders that affect the exposure-outcome relationship), and exclusion restriction (the genetic variant influences the outcome solely through its effect on the exposure, with no alternative direct causal pathways). Ensuring these assumptions are met is critical for the validity of causal inferences derived from MR studies (14). In parallel, the advent of microarray technology has revolutionized biological research by enabling high-throughput transcriptome analysis, transitioning investigations from single-gene studies to genome-wide expression profiling. Over the past two decades, microarray data have contributed substantially to advancing diverse areas of biological research, offering significant advantages for expression profiling and large-scale data analysis (15).
This study employed publicly available datasets, including GSE64763 and GSE31699, which are associated with ULM, alongside 436 oxidative stress-related genes. We chose GSE64763 because the dataset contains rich gene expression data covering multiple disease states and normal control samples, which is critical for identifying disease-associated genes and pathways. In addition, Dai Huang et al. (16) and Yumin Ke et al. (17) also used the GSE64763 dataset to conduct ULM-related research. The GSE31699 dataset, on the other hand, focuses on molecular marker studies of ULM, providing in-depth gene expression analysis, which is important for studying disease progression and treatment response mechanisms. Research by Lei Cai et al. (18) is based on this dataset. Genome-wide association study (GWAS) data for ULM and relevant candidate genes were sourced from the IEU OpenGWAS database. A comprehensive analytical framework was utilized, comprising differential expression analysis, weighted gene co-expression network analysis (WGCNA), MR, receiver operating characteristic (ROC) curves, and other statistical methodologies to identify oxidative stress-related biomarkers in ULM. Gene set enrichment analysis (GSEA) was performed to elucidate functional pathways associated with these biomarkers, while immune infiltration analysis was conducted to examine their role in the immune microenvironment. Additionally, potential therapeutic agents targeting the expression of these biomarkers were predicted (Supplementary Figure S1). Investigating the mechanisms underlying oxidative stress-related biomarkers in ULM is pivotal for advancing diagnostic and therapeutic strategies, as well as for guiding future research in this domain.
2 Materials and methods
2.1 Data source
The datasets GSE64763 and GSE31699 were extracted from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/). GSE64763, used as the training set, contains transcriptomic microarray sequencing data from 25 uterine leiomyoma (ULM) samples and 29 control samples, based on the GPL571 sequencing platform (19). GSE31699, serving as the validation set, consists of 16 ULM and 16 control samples on the GPL6947 platform (20). A total of 436 oxidative stress-related genes (OSRGs) (21) were sourced from the Molecular Signatures Database (MSigDB) (https://www.gsea-msigdb.org/gsea/msigdb) (22), using the ‘msigdbr’ (v7.5.1) R package.
2.2 Weighted gene co-expression network analysis
To identify the module most strongly correlated with ULM, the R package ‘wgcna’ (v1.71) (23) was employed to perform WGCNA analysis on the GSE64763 samples. Initially, all samples were clustered to detect outliers, which were subsequently removed. Soft thresholding and scale-free network coefficients were selected to ensure the construction of a scale-free network. The module with the highest correlation to ULM was then identified.
2.3 Differential expression analysis and functional enrichment analysis
The differentially expressed genes (DEGs) between ULM and normal samples were identified using the ‘limma’ (v3.50.1) R package (24), applying the cut-off criteria of adj.p.value < 0.05 and |log2(FoldChange)| > 0.5. The DEGs were visualized using a volcano plot with the ‘ggplot2’ package (v3.4.1) (25), and a heatmap was generated with the ‘ComplexHeatmap’ package (v2.14.0) (26). Next, upregulated DEGs were intersected with genes from the ULM module positively correlated with ULM, while downregulated DEGs were intersected with negatively correlated ULM module genes. And these two sets of intersecting genes were combined to obtain the ULM-associated genes, which were intersected with 436 OSRGs to yield the differentially expressed OSRGs (DE-OSRGs). Subsequently, gene ontology (GO) and kyoto encyclopedia of genes and genomes (KEGG) enrichment analyses based on GO and KEGG databases (https://www.geneontology.org, https://www.genome.jp/kegg/) were performed using the ‘clusterProfiler’ (v4.2.2) R package (27) to explore the biological functions and pathways associated with these candidate genes.
2.4 Mendelian randomization analysis
MR analysis was conducted to investigate the causal relationship between the exposure factors (DE-OSRGs) and the outcome (ULM). Using the extract_instruments function of the R package ‘TwoSampleMR’ (v0.5.6) (28), we read the data related to the exposure factors and screened the instrumental variables (SNPs) according to the following criteria: p = 5*10^-8 to ensure that the instrumental variables were significantly correlated with the exposure factors; Enable clump = TRUE to exclude the effect of linkage disequilibrium (LD); set r2 = 0.001 and kb = 200. Next, the extract_outcome_data function was used to read the outcome data and filter out instrumental variables that were associated with exposure factors but not with the outcome by setting proxies = TRUE and rsq = 0.8. In addition, we calculated the F-statistics of the instrumental variables, and when F < 10, it indicated that the instrumental variables were weak instrumental variables and needed to be excluded. The above correlation analysis and the removal of weak instrumental variables ensured that these instrumental variables satisfied the correlation assumptions. Subsequently, five univariate MR methods were employed: MR Egger, Weighted Median, Inverse Variance Weighted (IVW) (29), Simple Mode, and Weighted Mode (30). The IVW method is based on Mendelian laws of inheritance, which ensure the random distribution of locus alleles in the population. This property makes genetic variants (e.g., SNPs) ideal instrumental variables for probing causal associations between exposures and diseases. IVW is effective in reducing confounding factors and reverse causality bias in traditional observational studies, and has evolved over the years into a mature and reliable method for causal inference. Therefore, IVW has become our primary method for identifying biomarkers causally associated with ULM. A p-value less than 0.05 was considered evidence of a causal association. A scatter plot was used to assess the relationship between exposure factors and outcomes, while a forest plot identified the predictive exposure factors for each single nucleotide polymorphism (SNP) in relation to the outcome. Finally, a funnel plot was utilized to evaluate whether the MR analysis adhered to Mendel’s second law of randomization.
2.5 Sensitivity analysis
To assess the robustness of the MR results, a sensitivity analysis was performed. Firstly, the heterogeneity analysis was carried out by using the mr_heterogeneity function in the R package ‘TwoSampleMR’ (v0.5.6) (28), and if the Q p-value was lower than 0.05, it indicated the existence of heterogeneity. Next, the horizontal polytropy test was performed by the mr_pleiotropy_test function, and if the p-value was greater than 0.05, it indicated that there were no confounders. Finally, the mr_leaveoneout function was used to perform a leave-one-out (LOO) analysis (step-by-step culling of each SNP and calculation of the meta effect of the remaining SNPs) and the inverse variance weighted method was used to assess the changes in the results, which were visualized by forest plots. The exclusivity assumption was satisfied by sensitivity analyses that ensured the genetic variation influenced outcomes only through exposure factors, excluding the role of other potential pathways. At the same time, the independence assumption was satisfied by excluding instrumental variables associated with confounders. In addition, the proportion of variance explained may also affect the validity of the MR analysis. In general, a higher proportion of variance explained implies stronger instrumental variables, which improves the confidence of the results. If SNPs explain only a small amount of variance, this may lead to bias and invalid causal inferences.
2.6 Receiver operating characteristic curve
To evaluate the diagnostic potential of the identified biomarkers for ULM, an ROC analysis was performed for each individual gene using the R package ‘pROC’ (v1.18.0) (31). The diagnostic accuracy increased as the area under the curve (AUC) approached 1, while an AUC of 0.5 indicated no diagnostic value.
2.7 Establishment of alignment diagram
Based on biomarker expression, an alignment diagram was constructed using the ‘rms’ R package (32). The predictive accuracy of this diagram was assessed using a calibration curve, decision curve analysis (DCA), and clinical impact curve (CIC).
2.8 Functional similarity and correlation analysis of biomarkers
To investigate the functional similarities among biomarkers, functional annotations were carried out with the ‘GOSemSim’ R package (v2.24.0) (33), focusing on three aspects: biological processes (BP), cellular components (CC), and molecular functions (MF). The average similarity scores of genes at different levels were calculated, and biomarkers were ranked based on these average similarity values. Additionally, Spearman correlation analysis of biomarkers was conducted using the ‘corrplot’ R package (v0.92) (34).
2.9 Gene set enrichment analysis
First, the correlation coefficients between each biomarker and other genes were calculated, after which genes were ranked according to these coefficients, yielding a list of genes relevant to each biomarker. The C2: KEGG gene sets in the MSigDB database served as the reference gene sets for GSEA enrichment analysis.
2.10 Immunoinfiltration analysis
To examine the relationship between biomarkers and immune cells, the abundance of 22 types of immune-infiltrating cells in the training set samples was estimated. Wilcoxon tests were then applied to calculate differences in immune cell infiltration between the ULM and control groups, identifying differential immune-infiltrating cells. Spearman correlation analysis was conducted to evaluate the relationships between differential immune-infiltrating cells, and between the biomarkers and the differential immune-infiltrating cells.
2.11 Regulatory network
The Network Analyst platform (https://www.networkanalyst.ca/) was utilized to predict transcription factors (TFs) regulating the identified biomarkers. A TF-mRNA regulatory network was constructed using Cytoscape. Additionally, to explore the regulatory mechanisms of biomarkers in ULM, corresponding miRNAs were predicted using the TarBase (http://www.microrna.gr/tarbase) and miRTarBase (https://mirtarbase.cuhk.edu.cn) databases. The miRNA-mRNA relationship pairs obtained from both databases were integrated, and miRNet was then used to predict lncRNAs corresponding to the miRNAs. Ultimately, a lncRNA-miRNA-mRNA regulatory network was established.
2.12 Disease and drug analysis
To assess the role of biomarkers in other uterine diseases, their relationships were analyzed through the CTD database (https://ctdbase.org). The DGidb database (https://dgidb.org) was then employed to predict small-molecule drugs associated with the biomarkers, and a gene-drug interaction network was constructed. Protein structures of biomarkers were retrieved from the PDB database (https://www1.rcsb.org/), while 3D structures of the therapeutic compounds were sourced from the NCBI PubChem Compound database (https://www.ncbi.nlm.nih.gov/pccompound/). Molecular docking was performed using CB-Dock to identify the best binding conformations between protein targets and drugs. Pymol was used for the visualization of these interactions. Binding energies below -5 kcal/mol were considered indicative of strong binding affinity between the compound and its target.
2.13 Expression and validation of biomarkers
The expression levels of the biomarkers were assessed in both ULM and control samples. To further verify the accuracy of these results, the biomarkers were validated in the independent validation set.
2.14 Real-time quantitative polymerase chain reaction
A total of five paired normal and ULM samples were collected from The Affiliated Taian City Central Hospital of Qingdao University. All participants provided informed consent, and the study received approval from the hospital’s ethics committee.
Total RNA was extracted from the 10 samples using TRIzol reagent (Invitrogen, China), following the manufacturer’s protocol. RNA concentrations were measured using the NanoPhotometer N50. cDNA was synthesized through reverse transcription using the SureScript First-strand cDNA synthesis kit (Servicebio, China). RT-qPCR was performed on the CFX Connect Thermal Cycler (Bio-Rad, USA), and relative mRNA expression levels were quantified using the 2-ΔΔCT method. The sequences of all primers are available in Supplementary Table 1.
2.15 Statistical analysis
The data underwent processing and analysis utilizing R software (v4.2) (35). The Wilcoxon test was used for comparison between groups, and statistical significance was assessed based on a p-value below 0.05.
3 Results
3.1 Identification of key module genes associated with ULM
In the GSE64763 dataset, all samples were clustered without notable outliers (Supplementary Figure 2). With R2 set to 0.85, the soft threshold β was determined to be 7 (Figure 1A), enabling the construction of a scale-free network. Thirteen co-expression modules were identified through WGCNA (Figure 1B). Among these, the MEgreen module (cor = 0.68, 1305 genes) exhibited a strong positive correlation with ULM, while the MEbrown (cor = -0.8, 1724 genes) and MEpurple (cor = -0.67, 261 genes) modules showed strong negative correlations with ULM, marking them as key ULM-related modules (Figure 1C).
Figure 1. Identification of key module genes associated with ULM. (A) Soft threshold filtering. (B) Gene dendrogram and modules before merging. (C) Heatmap of module correlations with phenotypes.
3.2 Differential gene expression analysis and functional enrichment analysis of DE-OSRGs
In the GSE64763 dataset, 883 DEGs were identified between ULM and control samples, consisting of 372 upregulated and 511 downregulated genes (Figures 2A, B). Subsequently, 776 ULM-related genes were identified, including 301 upregulated genes positively correlated with ULM (Figure 2C) and 475 downregulated genes negatively correlated with ULM (Figure 2D). Through the intersection of 776 ULM-related genes and 436 OSRGs, 42 DE-OSRGs were identified (Figure 2E). GO analysis indicated that these DE-OSRGs were predominantly involved in cellular response to hydrogen peroxide, cellular response to toxic substance, cellular detoxification, response to oxidative stress, cellular response to oxidative stress, response to reactive oxygen species, cellular response to chemical stress, cellular response to reactive oxygen species, response to hydrogen peroxide and cellular oxidant detoxification pathways (Figure 2F). KEGG pathway analysis revealed enrichment in the IL-17 signaling pathway, endocrine resistance, hepatitis B, lipid and atherosclerosis, proteoglycans in cancer, relaxin signaling pathway, leishmaniasis, fluid shear stress and atherosclerosis, malaria, and TNF signaling pathway (Figure 2G).
Figure 2. Differential expression analysis and enrichment analysis. (A) Volcano plot of DEGs. (B) Heatmap of DEGs. (C) Venn diagram showing 301 intersection genes from up-regulated DEGs and ULM positively correlated module genes. (D) Venn diagram showing 475 intersection genes from down-regulated DEGs and ULM negatively correlated module genes. (E) Venn diagram showing 42 DE-OSRGs from the intersection of 776 ULM-related genes and 436 OSRGs. (F) GO analysis of DE-OSRGs. (G) KEGG pathway analysis of DE-OSRGs.
3.3 Identification of four biomarkers through MR analysis
MR analysis was conducted to assess the relationship between DE-OSRGs and ULM. F-statistics for SNP were shown in Supplementary Table S2. Based on the IVW method, four biomarkers—ANXA1 (p = 0.0474, odds ratio (OR) = 1.021), CD36 (p = 0.0013, OR = 1.074), MICB (p = 0.0025, OR = 1.032), and PRDX6 (p = 0.0215, OR = 1.042)—demonstrated significant causal relationships with ULM, all with p < 0.05. The OR values greater than 1 suggest these biomarkers are potential risk factors for ULM (Figure 3A). Scatter plots for the five algorithms showed positive slopes, confirming consistent and robust positive causal relationships with ULM (Figures 3B–E). Moreover, the MR effect sizes in the forest plots were all above 0, further supporting these biomarkers as ULM risk factors (Figures 3F–I). The funnel plot revealed uniform SNP distributions for ANXA1, CD36, MICB, and PRDX6, indicating consistency with Mendel’s second law of random assortment (Figures 4A–D).
Figure 3. MR analysis. (A) Forest plot showing MR analysis of ULM candidate genes. (B–E) Scatter plots showing positive causal relationships of ANXA1, CD36, MICB, and PRDX6 with ULM. (F–I) Forest plots predicting SNP loci exposure factors for outcome diagnosis.
Figure 4. Sensitivity analysis for MR results. (A–D) Funnel plots showing uniform SNP distributions for biomarkers ANXA1, CD36, MICB, and PRDX6. (E–H) Leave-one-out analysis for ANXA1, CD36, MICB, and PRDX6.
Sensitivity analysis revealed that the heterogeneity (Table 1) and pleiotropy (Table 2) tests for ANXA1, CD36, and MICB yielded p-values greater than 0.05, suggesting no significant heterogeneity or horizontal pleiotropy. Although the p-value for the heterogeneity test of PRDX6 was 0.0153, the IVW method remained unaffected, ensuring reliable MR results. The leave-one-out sensitivity tests further demonstrated that the SNPs for the exposure factors and outcomes were consistently aligned to the right of 0, with no substantial deviations, confirming stable results without excessive sensitivity to individual SNPs (Figures 4E–H). Therefore, ANXA1, CD36, MICB, and PRDX6 were confirmed as reliable biomarkers for ULM.
Table 1. MR heterogeneity test.
Table 2. MR pleiotropy test.
3.4 Diagnostic value of biomarkers
In the GSE64763 dataset, the AUC values for ANXA1, CD36, MICB, and PRDX6 were 0.927, 0.853, 0.768, and 0.849, respectively (Figure 5A). Similarly, in the GSE31699 dataset, all biomarkers exhibited AUC values exceeding 0.7 (Figure 5B), confirming their diagnostic potential.
Figure 5. ROC curve analysis. (A) ROC curve of biomarkers in the training set (GSE64764). (B) ROC curve of biomarkers in the validation set (GSE31699).
3.5 Establishment and validation of a predictive alignment chart
To enhance predictive efficiency, these biomarkers were used to construct a predictive alignment chart in the training set (Figure 6A). The p-value from the Hosmer-Lemeshow test was 0.54, indicating no significant difference between predicted and ideal values (Figure 6B). The DCA demonstrated that the alignment chart had a graceful benefit rate, underscoring its refined predictive capability (Figure 6C). CIC analysis further evaluated the clinical utility of the alignment chart, revealing a favorable net benefit across a wide and practical threshold probability range, suggesting the chart’s significant predictive value (Figure 6D).
Figure 6. Biomarker alignment diagram. (A) Alignment diagram model based on biomarkers. (B) Calibration curves. (C) DCA curves. (D) Clinical impact curves evaluating the alignment diagram’s predictive power.
3.6 Correlation and pathway analysis of biomarkers
To explore biomarker correlations, GO analysis indicated that MICB and ANXA1 shared higher average functional similarity (Figure 7A). CD36 and ANXA1 exhibited the strongest positive correlation with a coefficient of 0.666, while MICB and PRDX6 displayed the strongest negative correlation (cor = -0.636) (Figure 7B). Additionally, GSEA revealed that ANXA1, CD36, and MICB were significantly enriched in chemokine signaling pathways and neural active ligand-receptor interactions (Figures 7C–E), while PRDX6 was prominently enriched in oxidative phosphorylation and cytokine-cytokine receptor interactions (Figure 7F).
Figure 7. Correlation and pathway analysis of biomarkers. (A) Functional similarity analysis of biomarkers. (B) Heatmap of biomarker correlations. (C–F) GSEA for ANXA1, CD36, MICB, and PRDX6.
3.7 Correlation between immune infiltration and biomarkers
The relationship between biomarkers and immune infiltration was investigated by calculating the abundance of 22 types of immune-infiltrating cells (Figure 8A). Notably, follicular helper T cells and neutrophils were significantly more abundant in ULM samples compared to controls, while M0 macrophages, resting mast cells, and monocytes were more prevalent in the control group (Figure 8B). Spearman correlation analysis revealed a strong positive association between M0 macrophages and neutrophils, while M0 macrophages showed significant negative correlations with resting mast cells and monocytes (Figure 8C). Further analysis showed that ANXA1 was significantly associated with follicular helper T cells (p = 0.0097), resting mast cells (p = 0.0049), and M0 macrophages (p = 0.0038). Specifically, ANXA1 was positively correlated with resting mast cells (cor = 0.48) and negatively correlated with M0 macrophages (cor = -0.49) and follicular helper T cells (cor = -0.44). CD36 was significantly positively correlated with resting mast cells (p = 0.0121, cor = 0.43), whereas MICB exhibited a significant negative correlation with M0 macrophages (p = 0.0297, cor = -0.38) (Figure 8D).
Figure 8. Immune infiltration correlation with biomarkers. (A) Immune cell infiltration abundance in the training set (GSE64763). (B) Boxplots comparing differential immune-infiltrating cell abundance between ULM and normal. (C) Heatmap of correlations between differentially immune-infiltrating cells. (D) Lollipop plot showing correlations between immune-infiltrating cells and biomarkers. *p < 0.05; **p < 0.01.
3.8 Regulatory network of biomarkers
To investigate the potential mechanisms of the biomarkers in ULM, 36 TFs were predicted. Among these, GATA2 was associated with CD36, MICB, and PRDX6, while NR3C1 and YY1 were linked to CD36 and ANXA1 (Figure 9A). Additionally, 17 miRNA-mRNA interaction pairs were identified from two databases (Figure 9B), and 324 lncRNAs were predicted. Ultimately, a ceRNA network was constructed, involving 4 mRNAs, 15 miRNAs, and 324 lncRNAs. PRDX6 was regulated solely by hsa-mir-24-3p and various lncRNAs, including PVT1, OLMALINC, and VASH1-AS1. ANXA1 and CD36 were regulated by hsa-mir-335-5p and hsa-mir-26b-5p (Figure 9C).
Figure 9. Regulatory network of biomarkers. (A) Network of TFs and biomarkers (TFs in purple, biomarkers in orange). (B) miRNA-mRNA interaction pairs. (C) ceRNA regulatory network (biomarkers in orange, miRNAs in green, lncRNAs in pink).
3.9 Associated between biomarkers and uterine diseases
The study also revealed associations between the 4 biomarkers and various uterine diseases. CD36 was exclusively linked to uterine prolapse, uterine cervicitis, and uterine cervical diseases, while the other biomarkers were associated with conditions such as uterine hemorrhage, uterine anomalies, uterine cervical dysplasia, uterine neoplasms, and uterine cervical neoplasms (Figure 10A). Small-molecule drugs corresponding to the biomarkers were predicted as well. ANXA1 exhibited the strongest interaction with amcinonide (interaction score = 1.77), CD36 had a notable interaction with ABT-510 (interaction score = 15.46), and MICB interacted with ribavirin (interaction score = 1.93) (Table 3). A gene-drug interaction network was then established (Figure 10B). Subsequently, molecular docking was performed between the biomarkers and the drugs with the highest interaction scores. ANXA1 and amcinonide had the lowest binding energy of -8.6 kcal/mol (Figure 10C), while MICB and ribavirin exhibited a binding energy of -6.2 kcal/mol (Figure 10D), indicating strong binding affinities.
Figure 10. Gene-disease and gene-drug analysis. (A) Relationships between biomarkers and intrauterine diseases. (B) Relationships between biomarkers and small molecule drugs (biomarkers in orange, drugs in green). (C) Molecular docking between ANXA1 and amcinonide. (D) Molecular docking between MICB and ribavirin.
Table 3. Correspondences between hub genes and drugs.
3.10 Expression of biomarkers
In terms of expression, the biomarkers showed significant differential expression in the training set. Specifically, ANXA1, CD36, and MICB were markedly downregulated in ULM samples, whereas PRDX6 was upregulated (Figure 11A). These expression trends were consistent in both the validation set and RT-qPCR results, confirming the findings from the training set (Figures 11B, C).
Figure 11. Biomarker expression analysis. (A) Violin plots of biomarker expression in the training set (GSE64763). (B) Violin plots of biomarker expression in the validation set (GSE31699). (C) RT-qPCR results for biomarkers in clinical samples. *p < 0.05; **p < 0.01; ***p < 0.001; ****p < 0.0001.
4 Discussion
ULM is the most prevalent benign tumor in the female reproductive system, often asymptomatic and incidentally detected during routine medical examinations. Although the exact etiology of ULMs remains unclear and specific tumor markers for diagnosis are lacking (36), growing evidence implicates oxidative stress as a key factor in their development and progression (11, 37). Oxidative stress can induce abnormal proliferation of ULM cells, promote DNA synthesis, and alter cell cycle dynamics, thereby driving tumor growth (38–40). Additionally, it influences cell apoptosis, inflammatory responses, and angiogenesis, further exacerbating ULM pathogenesis (10, 41–43). However, the precise roles of oxidative stress-related genes in ULMs remain insufficiently understood. In this study, transcriptional datasets from the GEO database were analyzed, identifying four biomarkers (ANXA1, CD36, MICB, and PRDX6) associated with oxidative stress in ULMs through differential analysis, WGCNA, and MR. MR analysis revealed that these biomarkers exert significant effects on ULM. ROC curve analysis further underscored their diagnostic potential for ULMs. Additionally, a predictive alignment chart incorporating these biomarkers was developed, effectively forecasting ULM risk.
In the NCBI (National Center for Biotechnology Information), there are more than 10 datasets on ULM. These data sets are used in different studies for different purposes and sample characteristics. For example, the study by José A Castro-Martínez et al. (44) used two datasets, GSE12814 and GSE23112. GSE23112 focused on microRNA expression, while GSE12814 was associated with molecular signatures of the del(7q) UL subgroup. Through comprehensive transcriptome analysis, this study provides an in-depth understanding of the molecular mechanisms of uterine diseases and provides scientific support for future diagnosis and treatment strategies. Another study (45) combined four datasets, GSE64763, GSE45188, GSE30673, and GSE593, to deepen our understanding of the pathogenesis of ULM and to identify genes involved in the development of ULM. Among them, GSE64763 is consistent with the one selected in this study; GSE45188 showed the DNA methylation profile of myometria in uterine fibroids and non-fibroids. GSE30673 focuses on MED12 mutation; GSE593 contained genome-wide expression profiles of 5 uterine fibroids and 5 normal uterine tissues. Compared with other datasets, we finally decided to use the two datasets GSE64763 and GSE31699 after comprehensively considering factors such as research design, requirement matching degree and sample size.
ANXA1, a 37-kDa protein belonging to the annexin A family, plays multifunctional roles in both innate and adaptive immunity (46). In recent years, the role of ANXA1 in tumors has garnered increased attention (47), with studies showing its expression varies by tissue type (48). For instance, ANXA1 is overexpressed in melanoma, hepatocellular carcinoma, gastric cancer, and lung cancer (49), while its expression is diminished in prostate and esophageal cancers (50). Jamaluddin et al. conducted genetic and proteomic analyses on ULM and adjacent normal muscle tissues to characterize the expression patterns of extracellular matrix (ECM) proteins in MED12 mutation-positive and -negative ULMs. Among ECM-related proteins, ANXA1 was found to be downregulated in ULMs of various sizes (51), aligning with the findings of this study and suggesting potential therapeutic targets.
MHC Class I polypeptide-related sequence B (MICB), a non-classical HLA-I gene within the human major histocompatibility complex, is highly polymorphic, with 40 polymorphic sites and 26 protein variants identified to date (52–54). MICB is a stress-inducible protein primarily expressed on cell membranes, acting as a ligand for NK cell receptors (53). Its surface expression is minimal in normal cells but markedly upregulated in tumor cells or cells infected with viruses (55). This study is the first to implicate MICB in ULM, laying the groundwork for future investigations. However, further research is needed to elucidate the underlying mechanisms of this interaction.
CD36, a transmembrane protein and member of the class B scavenger receptor family, interacts with other transmembrane proteins on the cell surface to mediate ligand binding and signal transduction (56). It is widely expressed on the surface of microvascular endothelial cells and plays a pivotal role in regulating endothelial cell function (57). In tumor microvasculature, CD36 binds to thrombospondin-1 (TSP-1), mediating endothelial cell apoptosis within tumor vasculature (58). Knapp et al., in their study on energy substrate transport proteins in ULM, discovered that CD36 expression was reduced in patients with ULM compared to matched healthy muscle layers (59). These findings align with our study results.
PRDX6, a peroxidase enzyme composed of 224 amino acids, is expressed in various tissues and organs. Recent research has highlighted its significant role in cancer. PRDX6 is notably overexpressed in cervical cancer tissues, where its overexpression promotes proliferation, migration, and invasion of cancer cells while inhibiting apoptosis (60). Our findings similarly indicate elevated PRDX6 expression in ULM. Research has demonstrated that PRDX6 knockout in HepG2 cells results in slowed cell division, mitochondrial dysfunction, and cell cycle arrest at the G2/M phase (61).
GSEA of ANXA1, MICB, CD36, and PRDX6 revealed common enrichment in signaling pathways, including chemokine signaling and cytokine-cytokine receptor interactions. These pathways are hypothesized to play a role in ULM development. Xia et al., in their comprehensive analysis of four ULM-related mRNA datasets, suggested that extracellular matrix-receptor interactions could lead to abnormal extracellular matrix formation in ULM. Additionally, focal adhesion and cell adhesion molecules were implicated in ULM pathogenesis. Xia et al. further proposed that chemokine signaling and cytokine-cytokine receptor interactions are likely involved in ULM progression (45). Prates et al. demonstrated that ANXA1 influences cervical cancer development by activating the transcriptional expression of formyl peptide receptors (FPRs) and the inhibitor of DNA binding 1 (ID1) (62). FPR plays a central role in the chemokine signaling pathway, while ID1, a DNA transcription regulator, is closely linked to DNA replication. Our findings are consistent with their results, reinforcing the involvement of these pathways in ULM development.
Through in-depth analysis of the functional and pathway information enriched by GO and KEGG, we found several clues closely related to the pathogenesis or related hypotheses of ULM. Studies have shown that ROS, as a pro-inflammatory mediator, can regulate cell proliferation and is known to activate the MAPK/ERK pathway in endometriosis (63), suggesting that ROS may also be involved in the pathological process of ULM, although this hypothesis needs to be further verified. Notably, ganirelix as a potential therapeutic agent has demonstrated significant toxicity to ULM cells, suggesting that it may induce tumor reduction in a broad patient population (64). In addition, chronic inflammation caused by visceral fat plays a key role in cell differentiation and proliferation, which is a necessary condition for the onset of ULM. Increased visceral fat deposits may alter cholesterol composition, thereby increasing the risk of subclinical atherosclerosis in patients with ULM (65). Fibromodulin (FMOD), as an important component of the extracellular matrix, not only plays a role in structure, but also is regarded as a novel tumor-associated antigen for a variety of malignant tumors, including ULM, and has the potential as a biomarker for cancer diagnosis and treatment (66). Relaxin has also been suggested to be involved in the pathogenesis of ULM, providing a new perspective for disease research (67). At the same time, TNF-α up-regulates MMP-2 expression and promotes cell migration by activating extracellular signal-regulated kinase (ERK) signaling pathway in uterine fibroid smooth muscle cells, which provides a theoretical basis for the development of ULM treatment strategies targeting TNF-α and MMPs inhibitors (68).
Immunoreactive cells, including lymphocytes, macrophages, granulocytes, and dendritic cells, are crucial participants in immune responses. Recent studies suggest that the composition and abundance of these cells may influence the development and progression of uterine ULMs (69). Using the CIBERSORT algorithm, we analyzed the immunoreactive cell infiltration characteristics in ULM versus normal tissues, revealing significant differences in five types of immunoreactive cells between the two groups. Supporting our findings, Zannotti et al. reported a greater abundance of macrophages within and around ULMs compared to distant myometrial tissue (69). The dysregulation of macrophage proliferation, infiltration, and accumulation contributes to pathological fibrosis and uncontrolled tissue repair, with key roles played by factors such as MCP-1, GM-CSF, TGF-β, and TNF-α. Similarly, Protic et al. observed a higher density of CD68-positive macrophages in ULMs and adjacent myometrium compared to distant myometrium (70). In this context, TNF-α secreted by macrophages may enhance activin A mRNA expression in ULM cells, potentially leading to excessive extracellular matrix production, tissue remodeling, and tumor growth. However, a clinical study by Liu et al. found no significant differences in MCT-positive mast cells or CD45-positive leukocytes between ULM and normal tissues in premenopausal women (71).
In another study, Wang et al. identified significant overexpression of miR-30a, miR-34, and miR-24, in addition to Let-7, in ULM. These miRNAs may inhibit ULM proliferation and differentiation by suppressing HMGA2a, acting as protective factors against the overexpression of oncogenes such as RAS and MYC (72). Extensive research has also explored miR-196a’s role in various cancers, showing its involvement in key biological processes related to tumorigenesis. Hu et al. demonstrated that miR-196a regulates the proliferation, invasion, and migration of esophageal squamous cell carcinoma by targeting ANXA1 (73). In breast cancer, miR-196a expression correlates with certain HOXC genes involved in tumor progression (74), while in cervical cancer, it targets netrin 4, influencing cell proliferation and migration (75). Although ULMs are benign, their abnormal cellular proliferation mirrors that of malignant tumors, suggesting that miR-196a may play a similar role in ULM pathogenesis. Furthermore, miR-155-5p exhibits dual roles, acting as an oncogene or tumor suppressor depending on the cellular environment and cancer type. Research by Navarro et al. indicated that miR-155-5p is implicated in the inflammatory processes often associated with ULMs (76).
Our research identified a strong binding affinity between ANXA1 and amcinonide, as well as between MICB and ribavirin. Previous studies have demonstrated that amcinonide, an anti-inflammatory agent with affinity for glucocorticoid receptors, can induce depigmentation in both black and albino mice. It significantly reduces the number of DOPA-positive epidermal melanocytes in these animals (77). Ribavirin, a broad-spectrum antiviral agent effective against HCV, HIV, and RSV, has been shown to decrease activin A levels while increasing follicle-stimulating hormone concentrations in the serum and liver of Wistar rats (78). The drug-binding analysis in this study suggests that amcinonide and ribavirin may serve as potential therapeutic agents for ULM, interacting with ANXA1 and MICB, respectively. While these findings hold promise for therapeutic interventions in ULM, further investigation into the precise mechanisms of action is necessary.
In conclusion, this study identified oxidative stress-related biomarkers in ULM and provided a comprehensive analysis of drug networks and immune cell infiltration, shedding light on the molecular mechanisms of these biomarkers in ULM. The findings offer potential new avenues for ULM treatment.
However, several limitations should be noted. The relatively small sample size used in this study may have introduced bias and affected both the statistical power and the biological relevance of the results. Future studies will aim to include larger sample sizes to further validate the accuracy of these biomarkers and their association with disease. At the same time, we will also consider implementing cross-validation strategies on a broader dataset to further consolidate and extend our research results. Additionally, the study focused exclusively on gene expression levels, which do not necessarily reflect direct biological effects. To validate these findings, further animal and clinical studies employing different experimental approaches will be required.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author/s.
Ethics statement
The studies involving humans were approved by Medical Ethics Committee of The Affiliated Taian City Central Hospital of Qingdao University. The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study. Written informed consent was obtained from the individual(s) for the publication of any potentially identifiable images or data included in this article.
Author contributions
YL: Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Software, Supervision, Writing – original draft. HC: Data curation, Formal analysis, Investigation, Methodology, Writing – original draft. HZ: Data curation, Formal analysis, Software, Validation, Visualization, Writing – original draft. ZL: Data curation, Formal analysis, Investigation, Methodology, Software, Validation, Visualization, Writing – original draft. LS: Investigation, Software, Writing – original draft. CZ: Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Software, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing.
Funding
The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This study received funding from Tai ‘an City science and technology innovation development project (N0. 2021NS351) and TCM science and technology project of Shandong Province (NO. M-2022079).
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fendo.2024.1373011/full#supplementary-material
Supplementary Figure 1 | Experimental technical roadmap illustrating study workflow.
Supplementary Figure 2 | Sample clustering and ULM gene feature heatmap using WGCNA.
Supplementary Table 1 | Primer sequences used in the study.
Supplementary Table 2 | F-statistic information for SNP.
References
1. Lin LC, Chang HY, Kuo TT, Chen HY, Liu WS, Lo YJ, et al. Oxidative stress mediates the inhibitory effects of Manzamine A on uterine leiomyoma cell proliferation and extracellular matrix deposition via SOAT inhibition. Redox Biol. (2023) 66:102861. doi: 10.1016/j.redox.2023.102861
2. Katon JG, Plowden TC, Marsh EE. Racial disparities in uterine fibroids and endometriosis: a systematic review and application of social, structural, and political context. Fertil Steril. (2023) 119:355–63. doi: 10.1016/j.fertnstert.2023.01.022
3. Ali M, Ciebiera M, Vafaei S, Alkhrait S, Chen HY, Chiang YF, et al. Progesterone signaling and uterine fibroid pathogenesis; molecular mechanisms and potential therapeutics. Cells. (2023) 12:1117. doi: 10.3390/cells12081117
4. Qu Y, Chen L, Guo S, Liu Y, Wu H. Genetic liability to multiple factors and uterine leiomyoma risk: a Mendelian randomization study. Front Endocrinol (Lausanne). (2023) 14:1133260. doi: 10.3389/fendo.2023.1133260
5. Zhou Y, You CG. Lipoxin alleviates oxidative stress: a state-of-the-art review. Inflammation Res. (2022) 71:1169–79. doi: 10.1007/s00011-022-01621-y
6. Szechynska-Hebda M, Ghalami RZ, Kamran M, Van Breusegem F, Karpinski S. To be or not to be? Are reactive oxygen species, antioxidants, and stress signalling universal determinants of life or death? Cells. (2022) 11:4105. doi: 10.3390/cells11244105
7. Tubbs A, Nussenzweig A. Endogenous DNA damage as a source of genomic instability in cancer. Cell. (2017) 168:644–56. doi: 10.1016/j.cell.2017.01.002
8. Xu X, Kim JJ, Li Y, Xie J, Shao C, Wei JJ. Oxidative stress-induced miRNAs modulate AKT signaling and promote cellular senescence in uterine leiomyoma. J Mol Med (Berl). (2018) 96:1095–106. doi: 10.1007/s00109-018-1682-1
9. Li Y, Xu X, Asif H, Feng Y, Kohrn BF, Kennedy SR, et al. Myometrial oxidative stress drives MED12 mutations in leiomyoma. Cell Biosci. (2022) 12:111. doi: 10.1186/s13578-022-00852-0
10. Miyashita-Ishiwata M, El Sabeh M, Reschke LD, Afrin S, Borahay MA. Differential response to hypoxia in leiomyoma and myometrial cells. Life Sci. (2022) 290:120238. doi: 10.1016/j.lfs.2021.120238
11. AlAshqar A, Lulseged B, Mason-Otey A, Liang J, Begum UAM, Afrin S, et al. Oxidative stress and antioxidants in uterine fibroids: pathophysiology and clinical implications. Antioxid (Basel). (2023) 12:807. doi: 10.3390/antiox12040807
12. Richmond RC, Davey Smith G. Mendelian randomization: concepts and scope. Cold Spring Harb Perspect Med. (2022) 12:a040501. doi: 10.1101/cshperspect.a040501
13. Sanderson E, Glymour MM, Holmes MV, Kang H, Morrison J, Munafo MR, et al. Mendelian randomization. Nat Rev Methods Primers. (2022) 2:6. doi: 10.1038/s43586-021-00092-5
14. Larsson SC, Butterworth AS, Burgess S. Mendelian randomization for cardiovascular diseases: principles and applications. Eur Heart J. (2023) 44:4913–24. doi: 10.1093/eurheartj/ehad736
15. Sun MA, Shao X, Wang Y. Microarray data analysis for transcriptome profiling. Methods Mol Biol. (2018) 1751:17–33. doi: 10.1007/978-1-4939-7710-9_2
16. Huang D, Xue H, Shao W, Wang X, Liao H, Ye Y. Inhibiting effect of miR-29 on proliferation and migration of uterine leiomyoma via the STAT3 signaling pathway. Aging (Albany NY). (2022) 14:1307–20. doi: 10.18632/aging.203873
17. Ke Y, You L, Xu Y, Wu D, Lin Q, Wu Z. DPP6 and MFAP5 are associated with immune infiltration as diagnostic biomarkers in distinguishing uterine leiomyosarcoma from leiomyoma. Front Oncol. (2022) 12:1084192. doi: 10.3389/fonc.2022.1084192
18. Cai L, Liao Z, Li S, Wu R, Li J, Ren F, et al. PLP1 may serve as a potential diagnostic biomarker of uterine fibroids. Front Genet. (2022) 13:1045395. doi: 10.3389/fgene.2022.1045395
19. Barlin JN, Zhou QC, Leitao MM, Bisogna M, Olvera N, Shih KK, et al. Molecular subtypes of uterine leiomyosarcoma and correlation with clinical outcome. Neoplasia. (2015) 17:183–9. doi: 10.1016/j.neo.2014.12.007
20. Navarro A, Yin P, Monsivais D, Lin SM, Du P, Wei JJ, et al. Genome-wide DNA methylation indicates silencing of tumor suppressor genes in uterine leiomyoma. PloS One. (2012) 7:e33284. doi: 10.1371/journal.pone.0033284
21. Liu Q, Yu M, Zhang T. Construction of oxidative stress-related genes risk model predicts the prognosis of uterine corpus endometrial cancer patients. Cancers (Basel). (2022) 14:5572. doi: 10.3390/cancers14225572
22. Gu X, Lai D, Liu S, Chen K, Zhang P, Chen B, et al. Hub genes, diagnostic model, and predicted drugs related to iron metabolism in Alzheimer’s disease. Front Aging Neurosci. (2022) 14:949083. doi: 10.3389/fnagi.2022.949083
23. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinf. (2008) 9:559. doi: 10.1186/1471-2105-9-559
24. Phipson B, Lee S, Majewski IJ, Alexander WS, Smyth GK. Robust hyperparameter estimation protects against hypervariable genes and improves power to detect differential expression. Ann Appl Stat. (2016) 10:946–63. doi: 10.1214/16-AOAS920
25. Liu S, Xie X, Lei H, Zou B, Xie L. Identification of key circRNAs/lncRNAs/miRNAs/mRNAs and pathways in preeclampsia using bioinformatics analysis. Med Sci Monit. (2019) 25:1679–93. doi: 10.12659/MSM.912801
26. Gu Z, Hubschmann D. Make interactive complex heatmaps in R. Bioinformatics. (2022) 38:1460–2. doi: 10.1093/bioinformatics/btab806
27. Wu T, Hu E, Xu S, Chen M, Guo P, Dai Z, et al. clusterProfiler 4.0: A universal enrichment tool for interpreting omics data. Innovation (Camb). (2021) 2:100141. doi: 10.1016/j.xinn.2021.100141
28. Hemani G, Zheng J, Elsworth B, Wade KH, Haberland V, Baird D, et al. The MR-Base platform supports systematic causal inference across the human phenome. Elife. (2018) 7:e34408. doi: 10.7554/eLife.34408
29. Burgess S, Butterworth A, Thompson SG. Mendelian randomization analysis with multiple genetic variants using summarized data. Genet Epidemiol. (2013) 37:658–65. doi: 10.1002/gepi.2013.37.issue-7
30. Xu J, Zhang S, Tian Y, Si H, Zeng Y, Wu Y, et al. Genetic causal association between iron status and osteoarthritis: A two-sample Mendelian randomization. Nutrients. (2022) 14:3683. doi: 10.3390/nu14183683
31. Xu C, Sun D, Wei C, Chang H. Bioinformatic analysis and experimental validation identified DNA methylation-Related biomarkers and immune-cell infiltration of atherosclerosis. Front Genet. (2022) 13:989459. doi: 10.3389/fgene.2022.989459
32. Liu TT, Li R, Huo C, Li JP, Yao J, Ji XL, et al. Identification of CDK2-related immune forecast model and ceRNA in lung adenocarcinoma, a pan-cancer analysis. Front Cell Dev Biol. (2021) 9:682002. doi: 10.3389/fcell.2021.682002
33. Yu G, Li F, Qin Y, Bo X, Wu Y, Wang S. GOSemSim: an R package for measuring semantic similarity among GO terms and gene products. Bioinformatics. (2010) 26:976–8. doi: 10.1093/bioinformatics/btq064
34. Liu Z, Wang L, Xing Q, Liu X, Hu Y, Li W, et al. Identification of GLS as a cuproptosis-related diagnosis gene in acute myocardial infarction. Front Cardiovasc Med. (2022) 9:1016081. doi: 10.3389/fcvm.2022.1016081
35. Chan BKC. Data analysis using R programming. Adv Exp Med Biol. (2018) 1082:47–122. doi: 10.1007/978-3-319-93791-5
36. Koltsova AS, Efimova OA, Pendina AA. A view on uterine leiomyoma genesis through the prism of genetic, epigenetic and cellular heterogeneity. Int J Mol Sci. (2023) 24:5752. doi: 10.3390/ijms24065752
37. Fletcher NM, Abusamaan MS, Memaj I, Saed MG, Al-Hendy A, Diamond MP, et al. Oxidative stress: a key regulator of leiomyoma cell survival. Fertil Steril. (2017) 107:1387–1394 e1381. doi: 10.1016/j.fertnstert.2017.04.015
38. Sugino N, Karube-Harada A, Taketani T, Sakata A, Nakamura Y. Withdrawal of ovarian steroids stimulates prostaglandin F2alpha production through nuclear factor-kappaB activation via oxygen radicals in human endometrial stromal cells: potential relevance to menstruation. J Reprod Dev. (2004) 50:215–25. doi: 10.1262/jrd.50.215
39. Agarwal A, Allamaneni SS. Role of free radicals in female reproductive diseases and assisted reproduction. Reprod BioMed Online. (2004) 9:338–47. doi: 10.1016/S1472-6483(10)62151-7
40. Prusinski Fernung LE, Al-Hendy A, Yang Q. A preliminary study: human fibroid stro-1(+)/CD44(+) stem cells isolated from uterine fibroids demonstrate decreased DNA repair and genomic integrity compared to adjacent myometrial stro-1(+)/CD44(+) cells. Reprod Sci. (2019) 26:619–38. doi: 10.1177/1933719118783252
41. Kirschen GW, AlAshqar A, Miyashita-Ishiwata M, Reschke L, El Sabeh M, Borahay MA. Vascular biology of uterine fibroids: connecting fibroids and vascular disorders. Reproductio. (2021) 162:R1–R18. doi: 10.1530/REP-21-0087
42. AlAshqar A, Reschke L, Kirschen GW, Borahay MA. Role of inflammation in benign gynecologic disorders: from pathogenesis to novel therapiesdagger. Biol Reprod. (2021) 105:7–31. doi: 10.1093/biolre/ioab054
43. Maghraby N, El Noweihi AM, El-Melegy NT, Mostafa NAM, Abbas AM, El-Deek HEM, et al. Increased expression of fibroblast activation protein is associated with autophagy dysregulation and oxidative stress in obese women with uterine fibroids. Reprod Sci. (2022) 29:448–59. doi: 10.1007/s43032-021-00810-0
44. Castro-Martínez JA, Vargas E, Díaz-Beltrán L, Esteban FJ. Comparative analysis of shapley values enhances transcriptomics insights across some common uterine pathologies. Genes (Basel). (2024) 15:723. doi: 10.3390/genes15060723
45. Xia L, Liu Y, Fu Y, Dongye S, Wang D. Integrated analysis reveals candidate mRNA and their potential roles in uterine leiomyomas. J Obstet Gynaecol Res. (2017) 43:149–56. doi: 10.1111/jog.2017.43.issue-1
46. Shen X, Zhang S, Guo Z, Xing D, Chen W. The crosstalk of ABCA1 and ANXA1: a potential mechanism for protection against atherosclerosis. Mol Med. (2020) 26:84. doi: 10.1186/s10020-020-00213-y
47. Foo SL, Yap G, Cui J, Lim LHK. Annexin-A1 - A blessing or a curse in cancer? Trends Mol Med. (2019) 25:315–27. doi: 10.1016/j.molmed.2019.02.004
48. Guo C, Liu S, Sun MZ. Potential role of Anxa1 in cancer. Future Oncol. (2013) 9:1773–93. doi: 10.2217/fon.13.114
49. Feng J, Xiao T, Lu SS, Hung XP, Yi H, He QY, et al. ANXA1−derived peptides suppress gastric and colon cancer cell growth by targeting EphA2 degradation. Int J Oncol. (2020) 57:1203–13. doi: 10.3892/ijo.2020.5119
50. Zhao X, Ma W, Li X, Li H, Li J, Li H, et al. ANXA1 enhances tumor proliferation and migration by regulating epithelial-mesenchymal transition and IL-6/JAK2/STAT3 pathway in papillary thyroid carcinoma. J Cancer. (2021) 12:1295–306. doi: 10.7150/jca.52171
51. Jamaluddin MFB, Ko YA, Kumar M, Brown Y, Bajwa P, Nagendra PB, et al. Proteomic profiling of human uterine fibroids reveals upregulation of the extracellular matrix protein periostin. Endocrinology. (2018) 159:1106–18. doi: 10.1210/en.2017-03018
52. Visser CJ, Tilanus MG, Schaeffer V, Tatari Z, Tamouza R, Janin A, et al. Sequencing-based typing reveals six novel MHC class I chain-related gene B (MICB) alleles. Tissue Antigens. (1998) 51:649–52. doi: 10.1111/j.1399-0039.1998.tb03008.x
53. Ferrari de Andrade L, Tay RE, Pan D, Luoma AM, Ito Y, Badrinath S, et al. Antibody-mediated inhibition of MICA and MICB shedding promotes NK cell-driven tumor immunity. Science. (2018) 359:1537–42. doi: 10.1126/science.aao0505
54. Groh V, Rhinehart R, Secrist H, Bauer S, Grabstein KH, Spies T. Broad tumor-associated expression and recognition by tumor-derived gamma delta T cells of MICA and MICB. Proc Natl Acad Sci U.S.A. (1999) 96:6879–84. doi: 10.1073/pnas.96.12.6879
55. Lanier LL. NKG2D receptor and its ligands in host defense. Cancer Immunol Res. (2015) 3:575–82. doi: 10.1158/2326-6066.CIR-15-0098
56. Silverstein RL, Febbraio M. CD36, a scavenger receptor involved in immunity, metabolism, angiogenesis, and behavior. Sci Signal. (2009) 2:re3. doi: 10.1126/scisignal.272re3
57. Yang M, Silverstein RL. CD36 signaling in vascular redox stress. Free Radic Biol Med. (2019) 136:159–71. doi: 10.1016/j.freeradbiomed.2019.02.021
58. Wang S, Blois A, El Rayes T, Liu JF, Hirsch MS, Gravdal K, et al. Development of a prosaposin-derived therapeutic cyclic peptide that targets ovarian cancer via the tumor microenvironment. Sci Transl Med. (2016) 8:329ra334. doi: 10.1126/scitranslmed.aad5653
59. Knapp P, Chabowski A, Posmyk R, Gorski J. Expression of the energy substrate transporters in uterine fibroids. Prostaglandins Other Lipid Mediat. (2016) 123:9–15. doi: 10.1016/j.prostaglandins.2016.02.002
60. Hu X, Lu E, Pan C, Xu Y, Zhu X. Overexpression and biological function of PRDX6 in human cervical cancer. J Cancer. (2020) 11:2390–400. doi: 10.7150/jca.39892
61. Lopez-Grueso MJ, Lagal DJ, Garcia-Jimenez AF, Tarradas RM, Carmona-Hidalgo B, Peinado J, et al. Knockout of PRDX6 induces mitochondrial dysfunction and cell cycle arrest at G2/M in HepG2 hepatocarcinoma cells. Redox Biol. (2020) 37:101737. doi: 10.1016/j.redox.2020.101737
62. Prates J, Moreli JB, Gimenes AD, Biselli JM, Pires D’Avila SCG, Sandri S, et al. Cisplatin treatment modulates Annexin A1 and inhibitor of differentiation to DNA 1 expression in cervical cancer cells. BioMed Pharmacother. (2020) 129:110331. doi: 10.1016/j.biopha.2020.110331
63. Streuli I, Santulli P, Chouzenoux S, Chapron C, Batteux F. Activation of the MAPK/ERK cell-signaling pathway in uterine smooth muscle cells of women with adenomyosis. Reprod Sci. (2015) 22:1549–60. doi: 10.1177/1933719115589410
64. Salas A, García-García P, Díaz-Rodríguez P, Évora C, Almeida TA, Delgado A. New local ganirelix sustained release therapy for uterine leiomyoma. Evaluation in a preclinical organ model. BioMed Pharmacother. (2022) 156:113909. doi: 10.1016/j.biopha.2022.113909
65. Vignini A, Sabbatinelli J, Clemente N, Delli Carpini G, Tassetti M, Zagaglia G, et al. Preperitoneal fat thicknesses, lipid profile, and oxidative status in women with uterine fibroids. Reprod Sci. (2017) 24:1419–25. doi: 10.1177/1933719116689598
66. Pourhanifeh MH, Mohammadi R, Noruzi S, Hosseini SA, Fanoudi S, Mohamadi Y, et al. The role of fibromodulin in cancer pathogenesis: implications for diagnosis and therapy. Cancer Cell Int. (2019) 19:157. doi: 10.1186/s12935-019-0870-6
67. Li Z, Burzawa JK, Troung A, Feng S, Agoulnik IU, Tong X, et al. Relaxin signaling in uterine fibroids. Ann N Y Acad Sci. (2009) 1160:374–8. doi: 10.1111/j.1749-6632.2008.03803.x
68. Wang Y, Feng G, Wang J, Zhou Y, Liu Y, Shi Y, et al. Differential effects of tumor necrosis factor-α on matrix metalloproteinase-2 expression in human myometrial and uterine leiomyoma smooth muscle cells. Hum Reprod. (2015) 30:61–70. doi: 10.1093/humrep/deu300
69. Zannotti A, Greco S, Pellegrino P, Giantomassi F, Delli Carpini G, Goteri G, et al. Macrophages and immune responses in uterine fibroids. Cells. (2021) 10:982. doi: 10.3390/cells10050982
70. Protic O, Toti P, Islam MS, Occhini R, Giannubilo SR, Catherino WH, et al. Possible involvement of inflammatory/reparative processes in the development of uterine fibroids. Cell Tissue Res. (2016) 364:415–27. doi: 10.1007/s00441-015-2324-3
71. Liu ZQ, Lu MY, Sun RL, Yin ZN, Liu B, Wu YZ. Characteristics of peripheral immune function in reproductive females with uterine leiomyoma. J Oncol. (2019) 2019:5935640. doi: 10.1155/2019/5935640
72. Wang T, Zhang X, Obijuru L, Laser J, Aris V, Lee P, et al. A micro-RNA signature associated with race, tumor size, and target gene activity in human uterine leiomyomas. Genes Chromosomes Cancer. (2007) 46:336–47. doi: 10.1002/gcc.20415
73. Hu C, Peng J, Lv L, Wang X, Zhou Y, Huo J, et al. miR-196a regulates the proliferation, invasion and migration of esophageal squamous carcinoma cells by targeting ANXA1. Oncol Lett. (2019) 17:5201–9. doi: 10.3892/ol.2019.10186
74. Milevskiy MJG, Gujral U, Del Lama Marques C, Stone A, Northwood K, Burke LJ, et al. MicroRNA-196a is regulated by ER and is a prognostic biomarker in ER+ breast cancer. Br J Cancer. (2019) 120:621–32. doi: 10.1038/s41416-019-0395-8
75. Zhang J, Zheng F, Yu G, Yin Y, Lu Q. miR-196a targets netrin 4 and regulates cell proliferation and migration of cervical cancer cells. Biochem Biophys Res Commun. (2013) 440:582–8. doi: 10.1016/j.bbrc.2013.09.142
76. Navarro A, Bariani MV, Park HS, Zota AR, Al-Hendy A. Report of exosomes isolated from a human uterine leiomyoma cell line and their impact on endometrial vascular endothelial cells. Pharm (Basel). (2022) 15:577. doi: 10.3390/ph15050577
77. Quevedo WC Jr., Holstein TJ, Dyckman J, Nordlund JJ. Influence of depigmenting chemical agents on hair and skin color in yellow (pheomelanic) and black (eumelanic) mice. Pigment Cell Res. (1990) 3:71–9. doi: 10.1111/j.1600-0749.1990.tb00325.x
Keywords: uterine leiomyoma, oxidative stress, biomarkers, Mendelian randomization, transcriptome
Citation: Li Y, Chen H, Zhang H, Lin Z, Song L and Zhao C (2024) Identification of oxidative stress-related biomarkers in uterine leiomyoma: a transcriptome-combined Mendelian randomization analysis. Front. Endocrinol. 15:1373011. doi: 10.3389/fendo.2024.1373011
Received: 03 February 2024; Accepted: 05 November 2024;
Published: 21 November 2024.
Edited by:
Claire Perks, University of Bristol, United KingdomReviewed by:
Ping Yin, Northwestern University, United StatesMunichandra Babu Tirumalasetty, New York University, United States
Praveen Koganti, Sanford Burnham Prebys Medical Discovery Institute, United States
Copyright © 2024 Li, Chen, Zhang, Lin, Song and Zhao. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Chuanliang Zhao, zhaochuanliang2003@163.com