Autoimmune Disease-Related Hub Genes are Potential Biomarkers and Associated with Immune Microenvironment in Endometriosis

Introduction

Endometriosis (EMT) is a disease in which endometrium-like tissue grows outside the endometrium. EMT is an inflammatory endometrial disease that affects women of reproductive age and is progesterone and oestrogen dependent.¹ Approximately 200 million women and teens globally are affected by it, or between 7 and 10% of women of reproductive age.^2,3 The staging of endometriosis is based on the location, number, size and degree of adhesions of the ectopic endometrium, but the severity of symptoms or recurrence of endometriosis is not related to stage. For example, women with revised American Society of Reproductive Medicine (rASRM) stage I disease may have severe pain, infertility, or both, whereas stage IV endometriosis may be asymptomatic. When symptoms of EMT are severe, it can cause recurrent dysmenorrhoea, persistent pelvic pain, abnormal menstruation and non-uterine amenorrhoea, leading to infertility in women.⁴ According to some reports, EMT affects 71% to 87% of women who experience chronic pelvic pain.²

The current ideas for understanding the origins of EMT are as following, Sampson’s “menstrual blood reverse flow” theory was initially put forth in 1927. Sampson came to believe that the majority of cases of EMT were brought on by endometrial fragments that were absorbed through the fallopian duct during menstruation, either by change from endometrial to ectopic endometrium, and by peritoneal and ovarian growth.⁵ The word “incumbent intimal determinism”, which refers to the patient’s incumbent endometrial traits that influence whether the endometrial fragments can develop in ectopic plants by blood reflux, was first proposed by Lang et al.⁶ It’s possible that there’s a correlation between EMT and the body’s immune mechanism dysfunction, according to pathophysiology ideas like autoantibodies⁷ and autoimmune regulation.⁸ Normal human immune system modulation can quickly identify and eliminate foreign objects and lessen inflammatory reactions. Therefore, it is believed that immune disorder is the root cause of EMT. Other studies have shown that genetic mutations are associated with the development of EMT, such as the dramatic increase in mutant allele frequencies (MAF) for Kirsten rats arcomaviral oncogene homolog (KRAS) mutations in endometriotic epithelium suggests that endometrial tissue carrying KRAS mutations is transported in a retrograde direction to the ovarian surface, where specific KRAS mutations give them a selective advantage at that tissues and other ectopic tissues, leading to the development and widespread distribution of EMT.⁹ Other theories of pathogenesis are the stem cell theory;¹⁰ malignant tendency;¹¹ and so forth.

Studies illustrated that there are correlations between autoimmune diseases (AID) and the occurrence of EMT.¹² A number of studies have found a significantly higher prevalence of autoimmune thyroid disease in women with EMT than in controls.¹³ One of the potential associations is autoimmunity. EMT shares many features with autoimmune diseases, including polyclonal B-cell activation, abnormal T- and B-cell function, and inflammatory tissue damage.¹⁴ Humoral immunity also explains its association with autoimmune thyroid disease. A high reactivity of some autoantibodies, such as thyroid peroxidase antibodies (TPOAb), has been found in patients with EMT. In addition, the incidence of positive antinuclear antibodies (ANAs) appears to be higher in patients with EMT than in patients with Graves’ disease. Another possible link between Graves’ disease and EMT is the estrogenic effect:¹⁵ Graves’ disease is five times more common in women than in men. Estrogen is an important regulator of the immune system, acting as a regulator of cytokine expression, antigen presentation and B-cell lymphocyte production. In addition, a large number of clinical studies have confirmed that the incidence of autoimmune thyroid diseases, systemic lupus erythematosus,¹⁶ inflammatory bowel disease and rheumatoid arthritis in EMT patients is significantly higher than that of normal people,¹⁷ and there are correlations in pathogenesis, clinical symptoms, disease progression and prognosis.¹⁸

This study used Weighted gene co-expression network analysis (WGCNA) for feature selection, screening for autoimmune disease related genes (AIDGs) involved in EMT, and verifying the diagnostic value of hub genes in EMT through machine learning algorithms. At the same time, this study also revealed the specific impact of AIDGs hub genes on the changes in the immune microenvironment of EMT patients.

Methods and Materials

Gene Expression Data and Autoimmune Disease Related Gene

A public functional genomics database called Gene Expression Omnibus (GEO, http://www.ncbi.nlm.nih.gov/geo) offers high-throughput data for gene expression and microarray analysis. The EMT-related datasets (GSE141549 and GSE7305) were downloaded from GEO database.^19,20 All datasets are belonged to Homo sapiens; the GSE141549 data platform was GPL13376, and the GSE7305 data platform was GPL570. GSE141549 contained 338 samples, 61 eutopic (EU) and 175 ectopic (EC) endometrium samples are retained by screening the sample, and GSE7305 contained 20 samples, including 10 EC and 10 EU samples. The current investigation comprised all samples from these two datasets; GSE141549 was used to screen for the hub gene, and GSE7305 was utilized to confirm the diagnostic value. Standardized data processing was carried out using the R package “limma”.²¹ A total of 4186 AIDGs were downloaded from the AIDGs Database (GAAD, http://gaad.medgenius.info/genes/).

Identification of DEGs

The samples were divided into an EC group and an EU group based on the expression matrix generated from the dataset GSE141549. Using the R package “limma” with absolute values of log fold-change (log FC) > 1 and adj.P.Val 0.05 as thresholds, differentially expressed genes (DEGs) between the several groups were found. In the EC group, DEGs that had elevated and downregulated logFC values, respectively, were found. The DEG analysis results are presented as volcano maps created with the “ggplot2” R package.²²

Weighted Gene Co-Expression Network Analysis

The GSE141549 dataset’s expression data were analyzed using WGCNA R package, and a co-expression network using EC and EU samples as clinical characteristics was created. In order to identify and eliminate outliers, sample clustering was carried out using the “goodSamplesGenes” function. The WGCNA package’s “pickSoftThreshold” function was used to determine the soft threshold’s optimal filtering power value, and the “blockwiseModules” function was used to create the one-step co-expression matrix. The adjacency matrix was obtained by first calculating the adjacency of the expression levels across all genes. Next, the adjacency matrix can be converted into a topological overlap matrix (topological overlap matrix, TOM), and calculate the corresponding difference degree (1-TOM), and divide the module according to the dynamic shear tree standard. When 1-TOM is used as the measurement value, the minimum capacity of the module is set to 30 and the shear height of the module is 0.20. With a minimum of 100 genes in each module, the dynamic tree cutting method was employed to identify various modules. Then, to merge related modules, a threshold of 0.2 was chosen. The implementation of a correlation analysis between module membership and gene significance followed. The key module genes, differential genes (DEGs) and AIDGs set were screened and crossed to obtain the autoimmune disease gene-related differential genes (AID-DEGs). In addition, Venn diagrams were used to screen to hub genes.

Functional Enrichment Analysis

The AID-DEGs were assessed using the GO functional enrichment research, which includes biological process (BP), molecular function (MF), cellular component (CC), and KEGG pathway analysis. KEGG is a database that provides information on the molecular and higher-level functions of genes, including biochemical pathways. The “clusterProfiler” R package handled the annotation and visualization.²³ A hypergeometric test was used to do the enrichment analysis. The cut-off number for a statistically significant difference was set at P<0.05. “GOplot” packages were used to screen for enrichment of genes with more than 40 replicates in the enrichment pathway in the first 6 biological pathways.

Construction of the Protein–Protein Interaction Network

The STRING database (https://string-db.org) was utilized to generate the protein-protein interaction network (PPI) of the altered genes.²⁴ “Medium confidence (0.400)” interaction scores were established as the minimum requirements. The PPI network was then visualized using Cytoscape software.²⁵ Important PPI network modules were also found using the molecular complex detection (MCODE)²⁶ plug-in with degree cutoff = 2, maximum depth = 100, and k-score = 2. The Hub Genes were discovered via Cytohubba,²⁷ a plug-in for the Cytoscape software.

Analyze the Differences in Expression of Hub Genes and Their Diagnostic Value

Analyze the difference of hub gene expression in EC group and EU group by the “ggpubr” package, and “pROC” package²⁸ was used to analyze the diagnostic value of hub gene in EC group and EU group. GSE7305 as the verification set for verification the differences in expression and their diagnostic value.

Random Forest (RF) and Support Vector Machine (SVM) Algorithm

The Random Forest structure is utilized to determine the mean decrease impurity, commonly referred to as the Gini significance. Decision Trees are part of Random Forest. Every decision tree has internal nodes and leaves. The internal node decides how to divide the data set into two subgroups with related responses using the selected characteristic. For classification or information gain tasks, gini impurity or information gain criteria may be utilized, while variance reduction may be used for regression activities. Each feature’s impact on the impurity of the split can be measured (the feature with highest decrease is selected for internal node). We can monitor the average impurity decrease for each feature. The relative importance of a feature is determined by taking the average of all the trees in the forest. The Scikit-Learn Random Forest implementation provides this method (for both classifier and regressor). SVM was employed to evaluate the diagnostic value of specific genes taken together.

Detection of Immune Infiltrating Cells in EC and EU Endometriosis

The relative proportions of 22 immunological components were calculated using the CIBERSORT algorithm. A gene expression matrix and particular gene sets from 22 different types of immune cells were used in the simulation calculation, which was carried out 1000 times. These immune cells’ respective composition ratios in each tissue were calculated. Using the ESTIMATE algorithm, immune and stromal scores for each tissue were calculated.

Human Tissue Collection

The tissue samples in this study were obtained from 8 patients with EC and EU tissues paraffin block from post-operation patients with EMT in Anhui Maternal and Child Health Hospital from January 2020 to October 2022.

Immunohistochemical Examination (IHC)

EC and EU endometrium were collected, and paraffin sections were made by dewaxing, antigen repair and sealing, and then incubated overnight at 4°C with HRP-labeled polyclonal antibodies against CXCL12, PECAM1, NGF, CTGF and WNT5A. Goat anti-rabbit secondary antibody was incubated at room temperature for 50 min, and freshly prepared DAB color development solution was added dropwise, and the color development time was controlled under the microscope. Then, the results of CXCL12, PECAM1, NGF, CTGF and WNT5A staining in each tissue sample were quantified using Image Pro-Plus 6.0 software, and the average optical density (AOD) values were calculated.

Statistical Analysis

Cytoscape (version 3.9.1) and R (version 4.2.1) were utilized for the statistical analysis in this work. The bilateral test also had a P<0.05 statistical significance level. For determining the differences in variables among groups, the Student’s t-test, Kruskal–Wallis, Wilcoxon, or one-way ANOVA test was used. The Area Under Curve (AUC) values were used to determine each remarkable gene’s predictive potential. The correlation test Spearman and Pearson were applied to determine the associations between variables.

Results

Identify Differential Genes (DEGs)

We used high-throughput sequencing technology to perform whole-transcriptome sequencing on EC and EU samples from GSE141549 in order to investigate the pathophysiology of EMT. Firstly, principal component analysis (Figure 1A) was utilized to examine whether there was a difference between the EU and EC using “factoextra” tools based on sequencing data. Using “limma” packages, we also found 926 DEGs in EC group compared to EU group. 375 genes were down-regulated, and 551 up-regulated (Figure 1B).

Figure 1 DEGs were obtained by differential analysis of the eutopic (EU) and ectopic (EC) endometrium samples from the GSE141549 dataset. (A) Principal component analysis: red represents the EC group, and blue represents the EU group. (B) Volcano plot: gray denotes non-differentially expressed genes, while red and blue stand for up- and down-regulated DEGs, respectively.

WGCNA Analysis the Most Relevant to Endometriosis Module

To further explore genes related to EMT, we utilized the GSE141549 dataset and performed WGCNA. Genes with the top 10,000 median absolute values were filtered to construct the WGCNA network. Cluster analysis excludes samples with obvious outliers and erroneous clustering (Figure 2A). The remaining 218 samples were analyzed for the next step (Figure 2B). Choose 11 as the optimal soft threshold, R² = 0.9 (Figure 2C). In accordance with the dynamic tree species cutting method, the genes were separated into several modules based on the appropriate soft threshold. After merging, a total of nine modules (Figure 2D), including the dark-red module (3157 genes), the light-cyan module (215 genes), the light-yellow module (493 genes), the salmon module (1580 genes), the light-green module (143 genes), the dark-green module (497 genes), the black module (423 genes), the cyan module (1721 genes), and the grey module (1771 genes) were produced. And then, we calculated the module’s association with the clinical features. As shown in Figure 2E and F, the largest positive link with the EC was found in the cyan module.(cor = 0.86, P = 3 e–66). To determine the relationship between the color module and gene significance, a thorough calculation was done. The cyan module was closely related to EC (cor=0.9, P<1e-200; Figure 2G). Next, we will analyze the cyan module’s subsequent step.

Figure 2 Continued.

Figure 2 Continued.

Figure 2 Continued.

Figure 2 Identifying the essential modules linked to endometriosis tissues. (A and B) Sample cluster tree: the outlying samples are removed through the clustering book, blue represents EU, red represents EC; the upper half is a cluster, and the lower half is phenotype. (C) Screening for appropriate soft thresholds, the scale-free topology fit index value was 0.90, as indicated by the red line. (D) Dendrograms of genes with different similarity based on module colors of topological overlap and assignments. (E) Heat map of the correlation between modules and clinical characteristics. (F) Correlation of modules and EC endometrium. (G) Scatterplots of gene significance (GS) vs module membership (MM) in the cyan module.

Screening for AID-DEGs Associated with Endometriosis

Through the screening and intersection of genes in the cyan module (n=1721), gene sets downloaded from the GAAD (n=4186) and DEGs (n=926) according to gene names, a total of 114 genes related to EMT associated with Autoimmune diseases (AID-DEGs) were obtained. The results shown in the Venn diagram above (Figure 3).

Figure 3 Screening of autoimmune gene-related genes for key module genes involved in the development of endometriosis.

Function and Pathway Analysis of AID-DEGs Associated with Endometriosis

Positive regulation of “response to external stimulus”, “leukocyte migration”, “regulation of actin cytoskeleton organization”, “regulation of leukocyte chemotaxis”, “transmembrane receptor protein serine/threonine kinase signaling pathway”, “positive regulation of leukocyte migration”, “cellular extravasation”, and “regulation of angiogenesis” were all enriched in the 114 AID-DEGs. GO biological processes (Figure 4B). In molecular function (MF), “receptor ligand activity” and “signaling receptor activator activity” were the most significantly enriched pathways. “Collagen containing extracellular matrix” pathways were prominent in the cellular component (CC) pathways (Figure 4A). The AID-DEGs were shown to participate in the pathways “Cell adhesion molecules”, “Cytokinecytokine receptor interaction”, and “Leukocyte transendothelial migration”, according to KEGG pathway analysis (Figure 4D–E). Genes assigned to the top 6 biological process GO terms are displayed in (Figure 4C).

Figure 4 Continued.

Figure 4 Continued.

Figure 4 Functional enrichment analysis was performed for the AID-DEGs. (A) GO results of the analysis of AID-DEGs to show the enrichment results of BP, CC, and MF. (B) Circle diagram: showing the top 10 AID-DEGs. (C) The top 6 biological process GO terms and their assigned genes. (D–E) Results of the KEGG enrichment analysis of the AID-DEGs.

PPI Network and Identifying Hub Genes

The STRING online database (version 11.5) was used to import 114 AID-DEGs in total to create the PPI network. Cytoscape’s (version 3.9.1) processing of the protein network is displayed in (Figure 5A). The Cytohubba (version 0.1) plug-ins of Cytoscape were used to assess the 114 genes with interactions using five different types of algorithms (MCC, MNC, Degree, EPC, and Betweenness) (Table 1). A Venn diagram displaying the top 10 genes identified by each algorithm is depicted in Figure 5B. Five hub genes were identified from the overlapped region are Connective Tissue Growth Factor (CTGF), C-X-C Motif Chemokine Ligand 12 (CXCL12), Nerve Growth Factor (NGF), Platelet And Endothelial Cell Adhesion Molecule 1(PECAM1), Wnt family member 5A (WNT5A).

Table 1 Top 10 Essential Genes for 5 Algorithms (MCC, MNC, Degree, EPC, and Betweenness)

Figure 5 Hub genes of AID-DEGs were selected. (A) Protein interaction analysis of 114 genes using the STRING online website. (B) Venn diagram intersection of the top 10 genes.

Analysis and Exploration Differences in the Expression and the Diagnostic Value of Hub Genes

The differences between the expression of hub genes in the EU and EC were examined using the GSE141549 dataset. As illustrated in Figure 6A, the expression levels of CTGF (P<2.2e-16), CXCL12 (P<2.2e-16), NGF (P<2.2e-16) and PECAM1 (P<2.2e-16) were higher in EC than in EU. We analysed the expression of hub genes in EC and EU tissues in different types of endometriosis and found that the expression of hub genes in EC tissues in different types of endometriosis was consistent relative to the endometrium, as shown in Supplementary Figure 1. In particular, EU exhibits higher levels of WNT5A (P<1.9e-14) expression than EC. As shown in Figure 6D, the expression of hub genes in GSE7305 was the same as GSE141549. Finally, using the “pROC” R package, this work evaluated the diagnostic value of the hub genes in EMT with GSE141549, with the following results (Figure 6B): WNT5A:AUC=0.8215, CTGF:AUC=0.9048, NGF:AUC=0.9318. PECAM1:AUC=0.9295, CXCL12:AUC=0.9327. By displaying ROC curves, the diagnostic value of five hub genes was confirmed (Figure 6C), the result as follows: PECAM1: AUC=0.97, CXCL12:AUC=0.99, WNT5A:AUC=0.82, CTGF:AUC=0.79, NGF:AUC=1. The validation results were consistent with the analysis results, and the hub genes was expressed differently in EMT and had diagnostic value.

Figure 6 Continued.

Figure 6 Analysis of differential expression and diagnostic value of hub genes in EU and EC. (A): Analysis of CTGF, CXCL12, NGF, PECAM1, and WNT5A in EU and EC. (B): At GSE141549, the diagnostic value of CTGF, CXCL12, NGF, PECAM1, WNT5A in EU and EC were analyzed. (C):The diagnostic value of CTGF, CXCL12, NGF, PECAM1, and and WNT5A was verified using the GSE7305 dataset. (D): To verify the expression difference of CTGF, CXCL12, NGF, PECAM1, and WNT5A, using the GSE7305 dataset. (E): Explore the difference in the importance of hub genes for the diagnosis of EU and EC by using the random forest algorithm. (F): Value of joint diagnosis of hub genes using SVM algorithm.

Results of Diagnostic Value Testing of Hub Genes via RF and SVM Algorithms

The gini importance of the random forest was used to compare the importance of CXCL12, PECAM1, NGF, CTGF, and WNT5A for the diagnosis of the EMT, as shown in the figure (Figure 6E), CXCL12: MeanDecreaseGini=26.127; PECAM1: MeanDecreaseGini=20.291; NGF: MeanDecreaseGini=24.110;CTGF:MeanDecreaseGini=11.792; WNT5A:MeanDecreaseGini = 7.766. Using the SVM algorithm, the 236 samples were divided into the training set and the test set according to the ratio of the train: test 7:3, and the joint diagnosis model of hub genes for EMT was constructed, and the ROC curve was drawn (Figure 6F) AUC=0.945.

Relationship Between AID-DEGs Hub Genes and Immune Infiltrating Cells in EC and EU

Immune traits were evaluated in relation to immune-checkpoint, HLA, receptor, chemokine, and abundance levels of immune cells. When compared to the EU group, the levels of the majority of immunological checkpoints, HLAs, receptors, and chemokines are higher in the EC group (Figure 7A). We determined the proportions of 22 immunological components between EC and EU endometrial tissues using the CIBERSORT method (Figure 7B). Figure 7C showed the intricate interactions between these immune compositions, with resting macrophages and mast cells particularly associated with the majority of immune compositions. Additionally, compared to EU, most immune cells were more abundant in EC endometrium tissues. The five hub genes CXCL12, PECAM1, NGF, and CTGF were positively correlated with monocytes, macrophages, B cell memory, mast cells, and negatively correlated with regulatory T cells and NK cells, according to further research. WNT5A was substantially connected with regulatory T cells and NK cells, but was not adversely correlated with monocytes, macrophages, or mast cells (Figure 7D).

Figure 7 Continued.

Figure 7 The association of hub genes with infiltrating immune cells. (A) Display the 22 immune cells parts for each sample. (B and C) Heatmap, box plots showing the differential infiltration of immune cells between EU and EC groups. (“*” represents P<0.05, “***” represents P<0.001, “****” represents P<0.0001) (D) Lollipop chart displaying the relationship between immune cells and hub genes.

Hub Genes Expression in EC and EU Tissues

The results of immunohistochemical staining revealed that CXCL12 and CTGF were mainly expressed in the cytoplasm and mainly localized in endometrial glands and mesenchymal cells; WNT5A was mainly located in the cytoplasm of glands; PECAM1 was mainly located in the cytoplasm of vascular endothelial cells and mesenchymal cells; NGF was mainly located in the cytoplasm of endometrial glands (Figure 8A). The staining result of hub genes in EC and EU endometrium was analyzed by image pro plus software, and the AOD values were calculated. CXCL12, NGF, CTGF, and PECAM1 showed higher AOD in EC endometrium than in orthotopic endometrium (p<0.05) (Figure 8B), while WNT5A showed low expression in EC endometrium.

Figure 8 Continued.

Figure 8 Continued.

Figure 8 Immunohistochemical staining results of the EU and EC endometrium. (A) Immunohistochemical staining results of CXCL12, NGF, WNT 5 A, CTGF and PECAM1 at 400 × magnification in the EU and EC endometrium. (B) Difference in the average optical density (AOD) values of CXCL12, NGF, WNT5A, CTGF, and PECAM1 in the EU and EC endometrium.

Hub Genes Expression in Menstrual Cycle

We analyzed the relationship between the screened hub genes and menstrual cycle and found no significant difference in the expression of CXCL12, CTGF, NGF and WNT5A during the menstrual, proliferative and secretory phases, indicating that the menstrual cycle had no effect on the expression of CXCL12, CTGF, NGF and WNT5A; while the expression of PECAM1 was higher in the menstrual phase than in the secretory and proliferative phases, as shown in the Supplementary Figure 2.

Discussion

EMT is a non-malignant gynecological illness with an uncertain mechanism. Laparoscopy is still the gold standard for diagnosis EMT, and no particular treatment targets or biomarkers have been found.^29,30 Short menstrual cycles, low body mass index, low waist-to-hip ratio, and low gestational age are risk factors for endometriosis, but studies have also found an increased risk of melanoma, non-Hodgkin’s lymphoma, thyroid cancer, and endometrial cancer in endometriosis.³¹ The etiology of EMT is not known, but theories about it include the “implantation theory” and “existing endometrial determinism”. The pathophysiology of EMT is not known, but it is widely believed to be related to immunity, genetics, retrograde menstruation, and lymphatic blood. A meta-analysis showed that women with endometriosis are at increased risk of several autoimmune diseases and that EMT appears to be associated with autoimmune diseases of the body.³²

Autoimmune diseases³³ refers to the body to antigen immune response to their organs or tissue damage caused by a kind of disease, can be divided into systemic autoimmune diseases (such as rheumatoid arthritis, systemic lupus erythematosus) and organ-specific autoimmune diseases (such as chronic lymphocytic thyroiditis, hyperthyroidism, ulcerative colitis, mainly according to the typical clinical manifestations and related examination diagnosis, and the high incidence of women, the course of delay. These are primarily determined by the typical clinical manifestations and associated examination diagnosis It has been demonstrated in animal models that immunosuppressants such as thalidomide and triptolide, as well as biological agents such as etanercept and infliximab, are beneficial in reducing the size of ectopic lesions and EMT in the treatment of autoimmune disorders.^34,35

Through this study, we screened out the AIDGs hub genes related to EMT, CXCL12, PECAM1, NGF, CTGF, and WNT5A, and analyzed the diagnostic value of hub genes and the correlation with immune microenvironment. First, by analyzing the immune microenvironment of the EU and EC endometrium, it was found that macrophages, mast cells and T cells CD4 memory quiescent in the EC endometrium had a high percentage, while NK cells and plasma cells accounted for a relatively small percentage. By investigating the relationship between the hub genes and the immune environment, we found that CXCL12, PECAM1 and NGF were all strongly associated with mast cells, macrophages and activated NK cells, respectively. CXCL12, PECAM1 and NGF were highly expressed in EC endometrium. Therefore, we speculate that CXCL12, PECAM1 and NGF may promote EMT by activating macrophages, mast cells and inhibiting NK cell activation. In contrast, WNT5A was lowly expressed in EC endometrium and was positively correlated with NK cells according to immune infiltration analysis. It is hypothesized that NK cells are absent in EC endometrium due to the low expression of WNT5A.

CXCL12, an essential homeostatic chemokine, is also referred as stromal cell-derived factor-1 (SDF-1).³⁶ In physiological and pathological processes such as embryogenesis, hemocyogenesis, angiogenesis, and inflammation, CXCL12 plays a critical role.³⁷ Hematopoietic progenitor cells, stem cells, endothelial cells, and the majority of leukocytes can all be activated and induced to migrate.³⁸ The CXCL12/CXCR 4 axis is extensively expressed in EMT, and activation of the AKT pathway through the CXCL12/CXCR 4 axis promotes the invasiveness of uterine EC endometrial cells and the development of ectopic focal blood vessels.³⁹ PECAM-1 (CD31)⁴⁰ can be identified on a range of cells in the vascular compartment.⁴¹ Due to its involvement in numerous mechanisms, including leukocyte movement through intracellular junctions, cell death, angiogenesis, etc, PECAM-1 plays a variety of diverse and broad roles in inflammation.⁴² Inhibiting PECAM1 in animal models of arthritis can stop the advancement of inflammation. Additionally, it has been discovered that through stimulating the NF-κB signaling pathway, PECAM1 and CXCR 4 can interact to increase the inflammatory response in pulpitis samples.⁴³ A literature review found that NGF is associated with severe pain in EMT.⁴⁴ Due to the recruitment and polarization of macrophages in EMT tissue. Macrophages express different nerve growth factors (NGF) and NGF causes neurogenic imbalance in an estrogen-dependent manner.⁴⁵ Patients with EMT have neuropathic pain, which is a result of estrogen-controlled neuroimmune interactions, which sensitize peripheral neurons. CTGF is a diagnostic and therapeutic target for rheumatoid arthritis (RA) due to its role in cell adhesion, fibrosis, and angiogenesis in the early stages of the disease.^46,47 Accordingly, it is presumed that CTGF may influence EMT to transform from a stromal fibroblast phenotype to a myofibroblast-like phenotype, which is responsible for the characteristics of EMT.⁴⁸ WNT5A is an example of a ligand that can control cell motility and polarity by signal transduction and can also function in macrophages as a pro-inflammatory agent to cause inflammation.⁴⁹ It was discovered to be expressed in a number of autoimmune disorders and was linked to disease development, such as systemic lupus erythematosus studies that linked the expression of WNT5A protein to disease progression.⁵⁰ It had discovered that WNT5A can transmit RYK / RhoA / ROCK signaling and WNT / Ca ^{2 +} expression through inflammatory mediators, which is favorable to RA migration and invasion.⁵¹ Recent research demonstrates that WNT5A reduces the development of monocytes into inflammatory M1 macrophages and generates an immunosuppressive macrophage phenotype. This is consistent with what we discovered.⁵²

Admittedly, there are some limitations in this study. First, we did not validate the immune microenvironment components of EC and EU endometrium in vivo. Second, the correlation between hub genes and immune microenvironment was also not validated in vivo. Therefore, in the follow-up study, we will continue to explore the association of hub genes with endometriosis severity and pharmacological treatment, which will also allow us to explore new ideas about the pathogenesis and immunotherapy of EMT.

Conclusion

In conclusion, in our study, bioinformatics and experimental method were used to identify hub genes (CXCL12, PECAM1, NGF, CTGF, and WNT5A) for autoimmune diseases associated with the immune microenvironment of EMT. Data analysis and in vitro experiments confirmed that hub genes are differentially expressed in EU and EC. Furthermore, hub genes are important for the diagnosis of EMT and the immune microenvironment, although our study also has some limitations. However, we must pay more attention to hub genes and further investigate the molecular mechanisms underlying the role of hub genes in EMT and autoimmune diseases, as well as their immunotherapeutic effects on EMT, to facilitate further development of immunotherapy in EMT.