Gene expression profiles from keloid patients were obtained from the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo) using datasets GSE145725 and GSE44270, both derived from Homo sapiens. GSE145725 includes microarray data from skin and scar tissues, comprising 19 samples: 9 keloid fibroblast samples (Keloid group) and 10 normal fibroblast samples (Normal group). The analysis was performed using the GPL16043 GeneChip^® PrimeView™ Human Gene Expression Array platform. GSE44270 consists of microarray gene expression profiles from skin and scar tissues, including a total of 12 samples: 9 keloid fibroblasts (Keloid group) and 3 normal fibroblasts (Normal group). The data were generated using the GPL6244 (HuGene-1_0-st) Affymetrix Human Gene 1.0 ST Array platform.

To integrate the GSE145725 and GSE44270 datasets, we first identified common genes present in both datasets. After aligning the datasets based on these common genes, they were merged into a single expression matrix. The batch effects between the two datasets were corrected using the ‘ComBat’ function from the ‘sva’ package in R (version 3.48.0). This step was crucial to ensure consistency in the Combined-dataset, enabling accurate downstream analyses.

Boxplots were generated to visualize the distribution of expression values across samples before and after batch correction, demonstrating that the batch effects were effectively eliminated in the merged dataset.

After batch correction, differential expression analysis was conducted on the combined dataset of 18 keloid and 13 normal samples using the ‘limma’ package in R (version 3.56.2). A linear model was used to compare the keloid and normal groups, with contrast matrices generated to evaluate expression differences. Differentially expressed genes (DEGs) were identified based on the criteria: |logFC| > 1 and adjusted P-value <0.05. Genes with logFC >1 were classified as upregulated, while those with logFC < −1 were classified as downregulated.

A predefined list of 1,596 oxidative stress-related genes (OSRGs) was obtained from the GeneCards website (https://www.genecards.org/, Version 5.21, updated on August 5, 2024). The DEGs were then intersected with these OSRGs to identify OSRDEGs. The OSRDEGs were then visualized using Venn diagrams, volcano plots, differential sorting plots, heatmaps, and a chromosome location map, all of which were generated using the Xiantaoxueshu online tool (www.xiantao.love, version 2022.1).

The STRING database (https://string-db.org/, version 11.0) was utilized to construct a PPI network, identifying interactions among the OSRDEGs. Key hub genes were identified based on their connectivity within the network. Furthermore, the GeneMANIA tool (http://genemania.org/, version 3.6) was employed to predict networks of genes functionally related to these hub genes, providing additional insights into their interactions.

ROC analysis was carried out to evaluate the diagnostic accuracy of these hub genes. The Area Under the Curve (AUC) values were calculated to determine the potential of these genes in distinguishing keloid from normal tissues, further refining the selection of hub genes with high diagnostic relevance.

To investigate the interactions between hub genes and miRNAs, we used the miRDB database (http://mirdb.org, version 6.0). For TFs interacting with hub genes, the CHIPBase database (https://rna.sysu.wsu.cn/chipbase/, version 3.0) and the HTFtarget database (http://bioinfo.life.hust.edu.cn/hTFtarget/) were consulted. Additionally, drug-gene interaction analysis was performed using the Comparative Toxicogenomics Database (CTD, http://ctdbase.org, version 5.0) to identify potential therapeutic targets. All interaction networks, including mRNA-TF and drug-gene interactions, were visualized using Cytoscape software (version 3.10.1).

Quantitative polymerase chain reaction (qPCR) was conducted to measure the mRNA expression levels of genes in keloid and normal skin tissues. Total RNA was extracted using TRizol reagent (Invitrogen). cDNA was synthesized using the PrimeScript™ RT Reagent Kit with gDNA Eraser (Takara) according to the manufacturer’s protocol. qPCR was performed using the 2× Universal SYBR Green Fast qPCR Mix (Abclonal) with the following reaction conditions: 1 cycle of 95°C for 30 min, 40 cycles of 95°C for 5 s, and 60°C for 30 s. All PCR primers were designed using NCBI Primer (Supplementary Table S1) and were synthesized by Beijing Tsingke Biotech Co., Ltd. The relative expression levels of various genes were analyzed using qPCR, comparing the expression between keloid and normal skin tissues.

All statistical analyses were conducted using R software (version 4.2.1) to assess gene expression differences and their diagnostic relevance. Gene expression differences between groups were evaluated using t-tests for two-group comparisons, while one-way ANOVA was applied for comparisons across multiple groups. Correlations between gene expression and immune cell infiltration were assessed using Pearson or Spearman correlation coefficients, depending on the data distribution. A p < 0.05 was considered statistically significant, and 95% confidence intervals (CI) were provided for key metrics to quantify the uncertainty of estimates. The research flowchart for this study is shown in Figure 1.

Boxplots were generated to visualize the distribution of expression values across samples before and after batch correction, demonstrating that the batch effects were effectively eliminated in the merged dataset (Supplementary Figure S1). Additionally, a predefined list of 1,596 oxidative stress-related genes (OSRGs) was obtained from the GeneCards website (Supplementary Table S1).

Through the intersection analysis of DEGs and OSRGs, 13 OSRDEGs were identified (Figure 2A). These genes are prominently highlighted in the volcano plot (Figure 2B) and are mainly located on chromosomes 1 and 5 (Figure 2C). The heatmap illustrates the differential expression of these 13 OSRDEGs between keloid and normal skin tissues (Figure 2D), where F3, FOXL2, EDN1, and NTF3 show higher expression in keloid tissues, while BDNF, FLT1, VCAM1 and ADRB2 are more highly expressed in normal skin. Additionally, the ranked plot of differential expression shows that F3, FOXL2, NTF3 and EDN1 exhibit larger fold changes and higher statistical significance, whereas BDNF and FLT1 show smaller expression changes (Figure 2E).

GeneRatio analysis revealed significant enrichment in biological processes, including the “positive regulation of the MAPK cascade” (GeneRatio 0.4, p < 0.01) and the “regulation of mononuclear cell proliferation” (GeneRatio 0.35, p < 0.03). Cellular components like the “external side of the plasma membrane” and molecular functions such as “neurotrophin receptor binding” and “growth factor activity” were also enriched. KEGG pathway analysis highlighted key pathways, such as “MAPK signaling” (GeneRatio 0.3, p < 0.02) and “fluid shear stress and atherosclerosis” (Figure 3A). The relationship between OSRDEGs and enriched GO and KEGG terms revealed that genes such as F3, FOXL2, NTF3, and EDN1 were particularly linked to these pathways, notably the “MAPK signaling pathway” (Figure 3B). The analysis combining logFC and Z-score underscored the enrichment of key biological processes, further emphasizing pathways such as “neurotrophin receptor binding” and “positive regulation of the MAPK cascade” (Figures 3C, D).

PPI network analysis (minimum required interaction score = 0.700) revealed complex interactions among the 13 OSRDEGs, with EDN1 (interaction score = 0.85), ADRB2 (interaction score = 0.82), NTF3 (interaction score = 0.80), and BDNF (interaction score = 0.78) emerging as central nodes, suggesting their critical roles in regulating keloid-related pathways (Figure 11A). These interactions outline a complex signaling network essential for cellular communication and keloid pathology. GeneMANIA further refined this network (Figure 11B), where thicker lines represent stronger interactions between nodes. Key genes such as EDN1, ADRB2, NTF3, and BDNF once again emerged as central hubs, with significant interaction strengths. Ultimately, 7 hub genes—EDN1, NTF3, VCAM1, ADRB2, BDNF, F3 and FLT1—were identified, all of which had interaction scores above 0.70, indicating their potential roles in keloid formation and progression. These genes will be the focus of further investigation to determine their specific contributions to the disease.

In the violin plot (Figure 12A), EDN1, NTF3, ADRB2, and VCAM1 were found to be significantly upregulated in keloid tissue, with very high statistical significance (***, p < 0.001). Compared to the normal skin group, these genes showed markedly elevated expression in keloid tissue, suggesting their potential importance in the pathogenesis of keloids. Next, ROC curve analysis was used to rank the AUC values of the 7 hub genes. The highest AUC was observed for ADRB2 (AUC = 0.979), followed by NTF3 (AUC = 0.940), VCAM1 (AUC = 0.910), and EDN1 (AUC = 0.850). According to AUC classification criteria, ADRB2, NTF3, and VCAM1—with AUC values above 0.9—exhibited high diagnostic accuracy, while EDN1, with an AUC between 0.7 and 0.9, showed moderate diagnostic accuracy. Through these analyses, we further refined our focus to 4 key genes: ADRB2, NTF3, VCAM1, and EDN1.

Our analysis revealed the interaction networks associated with the hub genes EDN1, NTF3, ADRB2, and VCAM1, including miRNA, transcription factor, and drug-gene interactions. For example, miRNAs such as hsa-miR-335-5p and hsa-miR-30b-5p were found to be associated with multiple genes (Figure 13A). Key transcription factors, such as RAD21 and SP1, may regulate the expression of these hub genes (Figure 13B). In the drug-gene interaction networks, several compounds were predicted to interact with these hub genes. For instance, valproic acid is predicted to interact with NTF3, while dexamethasone is associated with VCAM1, suggesting that these genes could serve as potential therapeutic targets for modulating gene activity in the treatment of keloids (Figures 13C–F).

The immune cell infiltration and correlation analysis in keloid tissues revealed significant associations between the expression of the hub genes ADRB2, NTF3, VCAM1, EDN1 and various immune cell types. The heatmap (Figure 14A) showed that among all immune cells, activated and resting mast cells exhibited the most significant correlations with these hub genes. Specifically, EDN1 had the strongest positive correlation with activated mast cells (Figure 14E, r = 0.460, p = 0.022), while NTF3 showed the strongest negative correlation with resting mast cells (Figure 14F, r = −0.425, p = 0.017). As a result, EDN1 and NTF3 were identified as the most important hub genes, warranting further investigation into their roles in the formation and progression of keloids.

In keloid tissues, EDN1 and NTF3 were significantly upregulated, accompanied by prominent fibrotic features such as thickened and disorganized collagen fibers. Clinical photographs and HE staining results showed significant differences in fibrosis between keloid and normal skin tissues. Furthermore, IHC analysis confirmed the high expression of EDN1 and NTF3 in fibrotic regions, suggesting their critical roles in fibroblast activity and ECM remodeling.

EDN1 is a potent vasoconstrictive peptide secreted by vascular endothelial cells, primarily regulating vascular tone and blood pressure through ET-A and ET-B receptors (Banecki and Dora, 2023). Under pathological conditions, overexpression of EDN1 is associated with various cardiovascular diseases, such as hypertension and atherosclerosis (Sutton et al., 2019). EDN1 promotes the generation of reactive oxygen species (ROS) by activating NADPH oxidase, leading to increased oxidative stress, endothelial dysfunction, and vascular remodeling (Cai and Harrison, 2000). This process suggests that EDN1 plays a critical role in ROS generation and fibrosis. In keloids, EDN1 likely exacerbates inflammatory responses and fibroblast overproliferation by promoting ROS generation and activating the MAPK signaling pathway, resulting in abnormal ECM deposition (Dagamajalu et al., 2020).

NTF3 belongs to the neurotrophin family and primarily participates in the survival, differentiation, and regeneration of neurons (Lin et al., 2018). By binding to the TrkC receptor, NTF3 promotes neural system development and functional maintenance (Yang et al., 2020). Although most studies on NTF3 have focused on the nervous system, its role in immune regulation has attracted increasing attention. NTF3 has been shown to influence immune cell activation and proliferation by modulating the MAPK and PI3K/Akt signaling pathways (Yang et al., 2022). Current studies on NTF3 are predominantly limited to the nervous system, and its interaction with ROS remains underexplored, necessitating further investigation.

This study found that the abnormal overexpression of key genes, such as EDN1 and NTF3, in keloid tissues is closely associated with significant fibrotic and inflammatory phenotypes. However, it remains controversial whether fibrosis and inflammation are the primary pathological drivers of keloid formation or secondary phenomena triggered by upstream events, such as genetic susceptibility, dysregulation of the TGF-β signaling pathway, or abnormal fibroblast proliferation (Ogawa and Regeneration, 2011). On the one hand, existing literature emphasizes that a positive feedback loop between persistent chronic inflammation and excessive extracellular matrix (ECM) deposition in keloids can sustain and even exacerbate excessive scar tissue overgrowth (McDougall et al., 2006). For instance, Ogawa et al. proposed that keloids can be considered the result of persistent chronic inflammation in the dermis, continuously stimulating fibroblast activation and inducing abnormal collagen deposition (Ogawa, 2017). On the other hand, some studies suggest that genetic factors or early gene mutations and signaling abnormalities are the initial pathogenic events, with fibrosis and inflammation potentially serving as amplifiers of these upstream pathological mechanisms, manifesting as secondary phenotypes at the tissue level (Zhang et al., 2009).

Functional enrichment analysis in this study revealed that key OSRDEGs (e.g., EDN1 and NTF3) were significantly enriched in signaling pathways closely associated with fibrosis and inflammation, including the MAPK signaling pathway, regulation of lymphocyte and monocyte proliferation, and cytokine activity regulation. ROS produced by OSRDEGs, as markers of oxidative stress, can activate key kinases in the MAPK pathway, such as ERK, JNK, and p38, through redox-sensitive signaling pathways (Son et al., 2011). These kinases further activate downstream transcription factors, such as AP-1 and NF-κB, which regulate gene expression to promote cell proliferation, inflammatory responses, and apoptosis (Geest and Coffer, 2009). In pathological conditions such as tumors and fibrosis, sustained activation of the MAPK pathway aggravates disease progression, induces persistent inflammatory responses, and promotes the increased secretion of cytokines such as TNF-α, IL-6, and IL-13 (Wagner and Nebreda, 2009). Persistently expressed inflammatory cytokines not only lead to overactivation of the immune system but also exacerbate tissue damage and fibrosis, such as abnormal fibroblast proliferation and excessive ECM deposition, which harden tissues and impair functionality (Kendall and Feghali-Bostwick, 2014).

Immune cell infiltration analysis revealed significant correlations between activated and resting mast cells and several key genes, including EDN1, NTF3, ADRB2, and VCAM1. Activated mast cells regulate local inflammation and fibrosis by releasing inflammatory mediators such as histamine, tryptase, and various cytokines (Galli et al., 2008a). Mast cell activation can be triggered by the binding of IgE to its high-affinity receptor (FcεRI) or by ROS generated during oxidative stress, thereby exacerbating inflammatory responses (Gilfillan and Beaven, 2011). In damaged tissues, mast cells release inflammatory factors that attract other immune cells and promote fibrosis (Artuc et al., 2002). Resting mast cells are typically inactive and maintain tissue homeostasis, initiating inflammatory responses only upon receiving specific signals, such as allergens (Theoharides, 1990). In this study, EDN1 was positively correlated with activated mast cells, suggesting that EDN1 may promote mast cell activation, drive fibrosis progression in keloid tissues, and amplify inflammatory responses. Chemokines released by mast cells can further attract monocytes and lymphocytes, enhancing local inflammation (Galli et al., 2008b). Conversely, NTF3 was negatively correlated with resting mast cells, suggesting that NTF3 may inhibit their resting state, promote mast cell activation, and increase inflammatory responses. NTF3 may also regulate MAPK and PI3K/Akt signaling pathways, influencing immune cell activation and proliferation (Liu et al., 2021). In this study, the high expression of EDN1 and NTF3 was significantly correlated with fibrosis markers (COL1A1, TGFB1) and inflammatory factors (IL6, TNFA), supporting their dual regulatory roles in fibrosis and inflammation. This suggests that fibrosis and inflammation may form a mutually reinforcing relationship in keloids. On the one hand, persistent chronic inflammation stimulates excessive fibroblast proliferation and collagen deposition; on the other hand, excessive ECM deposition further recruits and activates immune cells, amplifying the inflammatory response (Nangole and Agak, 2019). However, to fully elucidate the sequence and interplay between fibrosis and inflammation in keloid progression, further longitudinal studies and in vitro and in vivo functional experiments are required. For instance, dynamic monitoring of these molecules and signaling pathways during the early stages of keloid development (e.g., early wound healing phases) could be insightful. Meanwhile, gene knockout, overexpression, and pharmacological intervention studies targeting specific genes could help clarify their precise roles in keloid initiation and progression. Our findings suggest that although fibrosis and inflammation may not be the initial triggers, they play indispensable roles in the development and maintenance of keloids. Thus, early therapeutic interventions targeting the inflammation-fibrosis pathway may offer new approaches for the clinical management and prevention of keloids.

To enhance the representativeness and robustness of the study results, this research incorporated an independent external dataset, GSE7890, for validation. The results showed that although NTF3 was not significantly upregulated in the validation set, EDN1 and other genes (e.g., FLT1, FOXL2, VCAM1) displayed expression patterns consistent with the training set (GSE145725, GSE44270. The differences between the keloid and normal groups were statistically evaluated using the Mann-Whitney U test. Although some genes did not reach significance, likely due to the small sample size in the validation set, the overall results still support the importance of EDN1 and other genes in keloids. Furthermore, gene set enrichment analysis (GSEA) in both the training and validation datasets revealed significant enrichment of the IL-18 pathway, indicating the consistency of inflammatory signaling in keloids across different populations. WGCNA analysis identified EDN1 and NTF3 as key genes within the “brown” and “blue” modules, which were closely associated with fibrosis and inflammation. These modules were enriched with genes related to excessive ECM deposition and immune signaling. At the transcriptome-wide level, this further suggests that EDN1 and NTF3 may play dual regulatory roles in the fibrosis and inflammation processes of keloids. This analysis compensated for the limited statistical power of single-sample IHC studies and showed good concordance with immune cell-related heatmap findings. These findings further indicate that EDN1 and NTF3 may contribute to the development and progression of keloids by regulating fibrosis and inflammation.

The ROC analysis revealed that EDN1 and NTF3 exhibited high sensitivity and specificity in distinguishing keloid tissues from normal tissues, highlighting their potential value as biomarkers. This study also predicted miRNAs (e.g., hsa-miR-1264, mmu-miR-1933-5p), transcription factors (e.g., JUND, GATA3), and potential targeted drugs associated with EDN1 and NTF3. Previous studies have indicated that endothelin receptor antagonists, such as bosentan, can inhibit fibroblast activity and ECM deposition in various fibrotic diseases, warranting further evaluation in keloid models (Dhaun and Webb, 2019). Similarly, NTF3 may regulates MAPK/ERK and PI3K/Akt signaling pathways through TrkC receptors. Developing more specific TrkC inhibitors or antibodies, combined with localized drug delivery and multitarget approaches (e.g., p38/ERK inhibitors, anti-inflammatory agents), could potentially achieve better outcomes in suppressing the fibrosis-inflammation feedback loop (Kaplan and Miller, 2000).

Despite verifying the key roles of EDN1 and NTF3 through multiple analytical approaches, this study has some limitations. First, this study lacks functional experiments (e.g., gene knockdown or overexpression) to directly validate causal relationships. Second, the small sample size of the validation dataset might have affected the significance of some genes. Future studies will employ siRNA or CRISPR techniques in keloid-derived fibroblast or immune cell lines to observe their effects on ECM synthesis and inflammatory factor secretion. In addition, in vivo functional validation will be conducted as resources and conditions permit. Moreover, considering the ethnic and genetic susceptibilities of keloids, we provide demographic information (age, gender, lesion location, and ethnicity) of patients from the GSE145725 and GSE44270 datasets, presented in tabular form. Clinical keloid samples collected by our team were also summarized and presented in tabular format. Although these data do not fully encompass population variability, they enhance the extrapolative value and transparency of the study, laying the foundation for future multicenter or large-scale population studies.