Comparison Bioinformatic Analysis of Extracellular Vesicles-Related Genes and MicroRNAs in Breast Cancer

TL;DR

This study uses bioinformatics to analyze genes and microRNAs related to extracellular vesicles in breast cancer, identifying potential biomarkers and therapeutic targets.

Contribution

The study identifies specific EV-related genes and microRNAs with significant expression patterns and survival correlations in breast cancer subtypes.

Findings

CD36, DAB2, and CXCL14 were downregulated in breast cancer tissues, while TIMP1 was upregulated.

Higher levels of TIMP1, DAB2, and CXCL14 were associated with improved overall survival in breast cancer patients.

miR-222-3p was downregulated and positively correlated with TIMP1 and DAB2, while miR-181b-5p was upregulated and negatively correlated with CXCL14.

Abstract

Breast cancer (BC) remains a leading cause of cancer-related mortality in women, with treatment challenges due to the lack of targeted therapies. Extracellular vesicles (EVs) play a crucial role in BC progression by carrying bioactive molecules. This study analyzed EV-associated molecules (ENPEP, TIMP1, CD36, MARCKS, DAB2, CXCL14, miR-181b-5p, miR-222-3p) using bioinformatics tools. We used GEPIA2; Human Protein Atlas (HPA) 24.0; bc-GenExMiner v5.1; UALCAN 2022; Kaplan–Meier plotter 2025; ENCORI database v2.0; Enrichr-KG web tool 2021; Cancer Hallmark Enrichment tool 2025; Tumor, Normal, and Metastatic (TNM) plot database 2025; MicroRNA Target Prediction Database 6.0; TargetScan 8.0; and STRING database 12.0. CD36, DAB2, and CXCL14 were significantly downregulated, while TIMP1 was upregulated in BC tissues (p < 0.05). CD36, CXCL14, and DAB2 were predominantly low in triple-negative and basal-like subtypes, whereas TIMP1 was higher in HER2+, ER+, and PR+ tumors (p < 0.01). These changes correlated with promoter methylation patterns. Higher TIMP1, DAB2, and CXCL14 levels were associated with improved overall survival (p < 0.05). miR-222-3p was downregulated and positively correlated with TIMP1 and DAB2, while miR-181b-5p was upregulated and negatively correlated with CXCL14. TNM analysis confirmed these expression changes. Functional enrichment linked these molecules to key cancer hallmarks, including proliferation and angiogenesis. CD36, DAB2, CXCL14, TIMP1, miR-222-3p, and miR-181b-5p may serve as biomarkers for BC pathogenesis and potential therapeutic targets. Further studies are needed to validate these findings.

Figures

Full text

1. Introduction

Breast cancer (BC) is one of the most common malignancies among women worldwide and remains a leading cause of cancer-related death, despite advances in detection and treatment strategies [1]. The heterogeneity of BC and the development of resistance to therapy continue to pose major clinical challenges [2,3]. Therefore, identifying new molecular markers and understanding the mechanisms that drive disease progression are crucial for improving patient outcomes.

Extracellular vesicles (EVs), especially exosomes, have emerged as key players in cancer biology. These lipid-bound vesicles are secreted by cells and serve as carriers of bioactive molecules, including proteins, lipids, microRNAs (miRNAs), and long non-coding RNAs (lncRNAs) [4]. EVs play a pivotal role in intercellular communication, influencing tumor proliferation, angiogenesis, immune modulation, and metastasis [5].

Recent evidence indicates that EVs reflect the molecular characteristics of their cell of origin. In BC, the cargo of EVs—particularly regulatory RNAs and proteins—has been implicated in promoting tumor aggressiveness and therapy resistance [6]. For example, exosomal miRNAs can activate pro-metastatic signaling pathways or remodel the tumor microenvironment. Thus, EV content profiling may represent a promising source of non-invasive biomarkers for early detection, prognosis, and therapeutic monitoring [7].

Although several studies have examined the role of EVs in BC, most have focused on limited sets of molecules or lacked comprehensive bioinformatic integration. The molecular mechanisms regulating EV cargo selection and secretion also remain insufficiently understood. These gaps limit our ability to fully exploit EVs as diagnostic or therapeutic tools [8,9,10,11].

EVs, including exosomes and microvesicles, are lipid bilayer-enclosed particles secreted by virtually all cell types and play critical roles in intercellular communication. They facilitate the horizontal transfer of bioactive molecules such as proteins, mRNAs, and lncRNAs that can modulate the tumor microenvironment and influence cancer progression [12]. Several molecules identified in our study, including TIMP1, CD36, DAB2, and miR-181b-5p, have been previously reported to be enriched in EVs derived from BC cells, as documented in EV-specific repositories such as Vesiclepedia and ExoCarta [13,14]. These EV-associated molecules are known to contribute to hallmark cancer processes such as proliferation, angiogenesis, metastasis, and immune modulation [15]. Therefore, we focused our analysis on these EV-linked targets to explore their expression patterns, prognostic relevance, and potential functional interactions in BC subtypes.

In recent years, bioinformatics-based approaches have emerged as indispensable tools for unraveling the complex molecular networks underlying cancer progression. These analyses in silico facilitate the integration and interpretation of large-scale omics data from various public repositories, such as TCGA and GEO. By leveraging such datasets, researchers can systematically identify differentially expressed genes, non-coding RNAs, and signaling pathways that may serve as potential biomarkers or therapeutic targets. Importantly, these computational frameworks provide a cost-effective and scalable alternative to wet-lab experiments, allowing for hypothesis generation and validation without the immediate need for direct experimental assays [16]. This study is entirely based on bioinformatics analyses, with no experimental procedures conducted. The authors relied exclusively on previously published, publicly available datasets to investigate EV-associated molecules in BC. In this study, we aimed to address this gap by systematically evaluating a panel of EV-associated molecules in BC using integrated bioinformatics tools. We focused on ENPEP, TIMP1, CD36, MARCKS, DAB2, CXCL14, miR-181b-5p, and miR-222-3p—genes and miRNAs previously implicated in cancer but not comprehensively analyzed in the context of EVs in BC. Our approach contrasts with prior work by combining multiple datasets and tools to assess gene expression, prognostic significance, and potential regulatory interactions.

In order to investigate the expression and regulation of EV-related genes in BC, we systematically selected a panel of molecules (ENPEP, TIMP1, CD36, MARCKS, DAB2, CXCL14, miR-181b-5p, and miR-222-3p) based on their documented significance in the literature and their relevance to EV biology and BC pathogenesis. We utilized comprehensive public databases, including GEPIA2; Human Protein Atlas (HPA) 24.0; bc-GenExMiner v5.1; UALCAN 2022; Kaplan–Meier plotter 2025; ENCORI database v2.0; Enrichr-KG web tool 2021; Cancer Hallmark Enrichment tool 2025; Tumor, Normal, and Metastatic (TNM) plot database 2025; MicroRNA Target Prediction Database 6.0; TargetScan 8.0; and STRING database 12.0. These resources facilitated the assessment of gene/microRNA expression patterns, methylation status, prognostic significance, and molecular interactions across various BC subtypes.

2. Results

To explore the relevance of EV-associated molecules in BC, we performed a comprehensive multi-omics analysis integrating transcriptomic expression, subtype-specific profiling, promoter methylation, prognostic impact, and functional enrichment. We began by assessing the differential expression of selected EV-associated genes and miRNAs between tumor and normal breast tissues. This was followed by correlation analysis to investigate potential regulatory interactions. We then examined their prognostic significance using survival analysis, explored epigenetic regulation via promoter methylation status, and assessed functional implications through enrichment of biological pathways and cancer hallmarks. This sequential approach enabled us to characterize the expression behavior, clinical relevance, and mechanistic roles of these candidate molecules in the BC context.

2.1. Expression Profile Results

CD36, DAB2, and CXCL14 were statistically downregulated while TIMP1 was upregulated in the BC cohort compared to normal adjacent tissue (p < 0.05). However, there was no statistical change in ENPEP and MARCKS (p > 0.05) (Figure 1A). Additionally, Immunohistochemistry (IHC) evaluation results of the genes examined in HPA are shown in Figure 1B. These findings highlight significant alterations in EV-associated genes, whose differential expression in BC may reflect functional implications in EV-mediated communication.

2.2. Gene Correlation and Targeted Expression Analysis Results

Additionally, the correlations of the examined genes with each other are shown in Figure 2. Accordingly, while no correlation was found between TIMP1 and MARCKS or between CXCL14 and MARCKS (p = 0.13 and p = 0.36, respectively), a positive significant correlation was observed between all the others (p < 0.001).

Five criteria are suggested for categorizing tumors: estrogen receptor (ER) via IHC, progesterone receptor (PR) via IHC (all RNA-seq data), HER2 status via IHC (all RNA-seq data), intrinsic molecular subtypes determined by PAM50 (all RNA-seq data), triple-negative breast cancer (TNBC) identified through IHC (all RNA-seq data), and basal-like subtypes identified through PAM50. The analysis revealed that the genes CD36, CXCL14, DAB2, MARCKS, and ENPEP exhibited low expression levels, particularly in TNBC and basal-like subtypes (p < 0.01). Conversely, the TIMP1 gene demonstrated high expression in HER2+, ER+, and PR+ tumors (p < 0.01). Specifically, the CD36 gene showed significantly reduced expression in ER-, PR-, HER2+, TNBC, and basal-like groups (p < 0.001, p < 0.001, p = 0.0190, p < 0.001, p < 0.001, respectively) as illustrated in Figure 3A. The CXCL14 gene exhibited low expression levels in ER-, PR-, HER2+, TNBC, and basal-like tumors (p < 0.001) (Figure 3B). The DAB2 gene exhibited low expression in ER- and TNBC tumors (p = 0.0139, p = 0.0008, respectively), while it showed high expression levels in HER2+ tumors (p = 0.0373) (Figure 3C). The ENPEP gene was detected at significantly higher levels in PR-, HER2+, and basal-like tumors (p = 0.0199, p = 0.0016, p = 0.0124, respectively) (Figure 3D). The MARCKS gene demonstrated high expression in ER-, PR-, HER2+, TNBC, and basal-like tumors (p < 0.001) (Figure 3E). The TIMP1 gene exhibited low levels of expression in HER2+, ER-, and PR- tumors (p < 0.001, p = 0.0097, p < 0.001, respectively) (Figure 3F).

Five criteria are suggested for categorizing tumors: ER via IHC, PR via IHC (all DNA microarrays), HER2 status via IHC (all DNA microarrays), intrinsic molecular subtypes determined by PAM50 (all DNA microarrays), TNBC identified through IHC (all DNA microarrays), and basal-like subtypes identified through PAM50. Specifically, the CD36 gene showed significantly reduced expression in ER-, TNBC, and basal-like groups (p < 0.001, p < 0.001, p < 0.001, respectively) as illustrated in Figure 4A. The CXCL14 gene showed decreased expression levels in ER-, PR-, HER2+, TNBC, and basal-like tumors (p < 0.001) (Figure 4B). The DAB2 gene exhibited low expression in ER-, PR-, HER2+, TNBC, and basal-like tumors (p < 0.001) (Figure 4C). The ENPEP gene was detected at significantly lower levels in HER2, TNBC, and basal-like tumors (p = 0.0075, p < 0.001, p < 0.001, respectively) (Figure 4D). The MARCKS gene demonstrated high expression in ER-, PR-, HER2+, TNBC, and basal-like tumors (p < 0.001) (Figure 4E). The TIMP1 gene exhibited low levels of expression in PR- tumors (p = 0.0438, respectively) (Figure 4F).

2.3. The Promotor Methylation Analysis

CD36, CXCL14, and DAB2 were statistically hypermethylated, while TIMP1 was hypomethylated in BC compared to normal (p < 0.05). However, there was no statistical change in ENPEP or MARCKS (p > 0.05) (Figure 5). As a result of the examination conducted in BC subgroups, downregulation in expression profiles generally occurred, while upregulation in TIMP1 expression levels was detected compared to normal adjacent tissue (Figure 6).

2.4. The Overall Survival Status

The elevated expression levels of TIMP1, DAB2, and CXCL14 were statistically associated with a longer OS (p < 0.05). However, no significant association was observed between the expression of ENPEP, CD36, and MARCKS and overall survival (p > 0.05) (Figure 7). A decreased expression level of miR-222-3p was associated with a longer OS (p = 0.011). However, no significant association was observed between miR-181b-5p expression and overall survival (p > 0.05) (Figure 8).

2.5. The Relationship Between MicroRNA vs. BC and MicroRNA vs. RNA

miR-222-3p was statistically downregulated while miR181b-5p was upregulated in the BC cohort compared to normal adjacent tissue (p < 0.05). miR181b-5p was negatively correlated with CXCL14 while miR-222-3p was positively correlated with TIMP1, CD36, and DAB2 (p < 0.05). miR181b-5p was upregulated and showed a statistically significant negative correlation with CXCL14 (p = 0.00000719), suggesting an inverse expression pattern, while it was positively correlated with MARCKS (p = 0.00000040), indicating a possible co-expression. miR-222-3p was positively correlated with TIMP1, MARCKS, and TIMP1 (p = 0.000208, p = 0.000363, p = 0.0183, respectively), while it was negatively correlated with CD36 and CXCL14 (p = 0.0315, p = 0.000934, respectively) (Figure 9).

2.6. Gene Set Enrichment Analysis of ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14

According to KEGG pathway analyses, CD36 is a component of the cholesterol metabolism pathway, ENPEP is part of the renin-angiotensin system, and DAB2 is involved in the fat digestion and absorption pathway. Data from the Jensen Lab indicate that the MARCKS gene is associated with autistic disorder. ENPEP is linked to vaccinia, antiphospholipid syndrome, and bladder cancer. TIMP1 is related to ureteral diseases. DisGeNET data reveal associations between myelodysplastic syndrome and the genes CD36, DAB2, ENPEP, and TIMP1. Complete hydatidiform mole is associated with the DAB2 and ENPEP genes. Squamous cell carcinoma is linked to the expression of MARCKS, CD36, DAB2, ENPEP, and TIMP1. Choriocarcinoma is associated with MARCKS, DAB2, ENPEP, and TIMP1. Infectious skin diseases are associated with CD36 expression. Enrichment analysis results of ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14 are shown in Figure 10A,B. The enrichment analysis of these EV-linked genes suggests that their altered expression patterns may influence cancer-related pathways via extracellular vesicle transport mechanisms.

2.7. Tumor, Normal, and Metastatic (TNM) Analysis

The expression levels of the genes CD36, CXCL14, DAB2, MARCKS, ENPEP, and TIMP1 in normal, tumor, and metastatic tissues were analyzed using various statistical methods. Violin and box plots were employed to illustrate the expression levels of each gene across these tissue types. The Dunn test results indicated significant expression changes between normal and tumor tissues for all genes (p < 0.05). Notably, significant differences were observed for the CD36 (p = 2.22 × 10^−58^) and DAB2 (p = 4.36 × 10^−25^) genes. The differences between tumor and metastatic tissues were gene-specific; for instance, no significant changes were observed in some genes, such as CD36 and ENPEP (p = 0.42), whereas significant changes were detected in others, such as MARCKS (p = 5.40 × 10^−7^). A multivariate Cox regression analysis was conducted to assess the potential roles of these genes in cancer progression. Fold change (FC) values were used to reflect gene expression changes during the transition from normal tissue to tumor and metastatic stages. The DAB2 gene exhibited a 0.59-fold decrease between tumor and normal tissues, with a measurement of 0.67 in metastatic tissues. TIMP1 showed a 2.56-fold increase in tumors and a 1.57-fold increase in the metastatic stage, resulting in a total change of 4.02-fold. CD36 demonstrated a significant decrease (0.24) during the transition from normal to tumor, with an even more pronounced decrease in the metastatic stage (FC_MvsN = 0.15) (Figure 11).

2.8. Cancer Hallmark Enrichment Analysis

Genes such as TIMP1, ENPEP, CD36, and CXCL14 have been identified as contributors to the mechanism of “sustaining proliferative signaling” (adj. p = 0.06175). Additionally, the genes CXCL14, DAB2, ENPEP, and MARCKS are significant in the “evading growth suppressors” category (adj. p = 0.01012). Furthermore, the genes TIMP1, CD36, and CXCL14 are implicated in the “sustained angiogenesis” mechanism (adj. p = 0.01012) (Figure 12 and Table 1).

2.9. Result of miRNAs Associated with Genes

The miRNAs that influence the CD36, CXCL14, DAB2, MARCKS, ENPEP, and TIMP1 gene groups are detailed in Table 2 and illustrated in the Venn diagram in Figure 13. The CD36 and CXCL14 are affected by the largest number of common miRNAs (n = 13). This is followed by the common miRNAs (n = 11 impacting the ENPEP and MARCKS as well as the ENPEP and DAB2 pairs). Among the triple gene groups, the combination of CD36, DAB2, and CXCL14 is influenced by the highest number of common miRNAs (n = 9).

2.10. Gene–Gene Interaction

String analysis results indicate that the top ten nodes associated with CXCL14 are CXCL17, CCL21, CXCL12, CXCL11, CCL13, CXCL9, CCL17, CCL11, CXCL5, and PPBP (Figure 14A). Similarly, the top ten nodes linked to CD36 include GP2, GP5, LYN, FYN, THBS1, TRL2, TRL4, THBS2, CD47, and TRL6 (Figure 14B). For TIMP1, the top ten associated nodes are IL10, IL6, MMP10, MMP2, FN1, CD63, MMP3, MMP9, MMP1, and MMP14 (Figure 14C). The nodes most closely related to MARCKS are BASP1, CALML5, CALML3, CALML2, CALML6, CALML4, MARCKS1, PRKCD, and CALML1 (Figure 14D). In the case of DAB2, the top ten nodes are LDLR, FCHO2, AP2B1, MYO6, LRP2, GRB2, DVL3, AXIN1, DAB2IP, and TGFBR1 (Figure 14E). Lastly, the nodes associated with ENPEP are C4B, C4A, ITIH4, CD209, DPP4, CXCR4, CD4, CCR5, AGT, and CALML6 (Figure 14F).

3. Discussion

To the best of our knowledge, this study is among the first to integrate multiple bioinformatics platforms to investigate the expression and regulation of ENPEP, TIMP1, CD36, MARCKS, DAB2, CXCL14, and miRNAs (miR-222-3p, miR-181b-5p) in BC. Our findings underscore their potential roles as prognostic biomarkers and therapeutic targets, particularly in aggressive BC subtypes such as TNBC. These molecules were selected based on their documented relevance to BC biology and prior evidence suggesting their presence in EVs according to literature and databases such as 5.1 and ExoCarta v1.0. However, we did not analyze EV-specific datasets nor apply vesicular cargo prediction tools; thus, our findings may be interpreted as suggestive of possible EV involvement, rather than conclusive.

EVs have emerged as pivotal mediators in cancer biology, reflecting the molecular profile of their cells of origin. Specific EV-associated miRNAs, such as miR-21, miR-155, miR-222-3p, and miR-181b-5p, have been implicated in metastasis and chemoresistance in BC. These miRNAs contribute to oncogenic signaling and immune modulation and are detectable in circulation, rendering them promising candidates for liquid biopsy-based diagnosis [17,18,19]. miR-181b-5p was upregulated and showed a negative correlation with CXCL14, while miR-222-3p exhibited positive correlations with TIMP1, DAB2, and CD36. Although these associations are consistent with predicted regulatory interactions, they may be interpreted as hypothesis-generating findings that require further experimental validation.

CD36, a scavenger receptor implicated in fatty acid metabolism and angiogenesis, has been extensively investigated in the context of cancer. While it is commonly linked to metastasis and poor prognosis in various malignancies, including BC [20], our analysis revealed that CD36 was significantly downregulated in the BC cohort, particularly in ER-, PR-, HER2+, TNBC, and basal-like subtypes. Notably, CD36 expression did not exhibit a significant correlation with OS in our cohort, suggesting a potentially context-specific or subtype-limited prognostic role. The observed downregulation may be attributed to promoter hypermethylation [21], as demonstrated by our methylation analysis. It is also noteworthy that CD36 expression is often diminished in tumor-associated stroma, contributing to metabolic rewiring and tumor progression. Previous studies have indicated that reduced stromal CD36 is associated with increased tumor aggressiveness [22,23]. However, in patients with HER2+ breast cancer, elevated CD36 expression correlates with a poor prognosis [20]. Interestingly, in our BC cohort, CD36 expression levels were not associated with OS. According to an IHC investigation, the significant downregulation of CD36 in ER-, PR-, HER2+, TNBC, and basal-like breast cancer may suggest a shift in metabolic and tumor microenvironmental dynamics. These findings suggest that CD36 may function as a context-dependent modulator of fatty acid metabolism and tumor–stroma interactions in BC. Its downregulation in aggressive subtypes, particularly TNBC and basal-like tumors, supports its potential role as a subtype-specific biomarker and warrants further functional investigation.

CXCL14 expression was similarly downregulated in BC tissues, particularly in TNBC and basal-like subtypes, known for their aggressiveness and limited treatment options. Previous research indicates that stromal CXCL14 expression serves as an independent predictor of decreased BC-specific and recurrence-free survival, especially in alpha ER- tumors [24,25]. Both in vitro and in vivo studies suggest that re-expressing CXCL14 can suppress tumor growth and metastasis by enhancing immune cell recruitment and modulating the tumor microenvironment [26]. In our dataset, this downregulation may be attributed to promoter hypermethylation and post-transcriptional regulation via miR-181b-5p. Furthermore, IHC findings suggest that CXCL14 downregulation in ER-, PR-, human epidermal growth HER2+, TNBC, and basal-like subtypes may lead to immune evasion, reduced tumor immunogenicity, and increased aggressiveness. Taken together, the consistent downregulation of CXCL14 across multiple aggressive BC subtypes and its association with shorter OS indicate that it may serve as a prognostic marker and therapeutic target, particularly through its role in immune modulation and tumor microenvironment regulation.

TIMP1, a tissue inhibitor of metalloproteinases, plays a context-dependent role in BC. While it has been implicated in extracellular matrix remodeling, proliferation, and chemoresistance, its prognostic value remains debated [27,28,29,30]. In contrast to these studies, another study suggested that low TIMP1 expression levels may serve as an independent prognostic and predictive biomarker for breast cancer [31]. Our analysis showed TIMP1 upregulation in luminal and TNBC subtypes and an association with improved OS, contrasting with some studies linking high TIMP1 levels to poor outcomes. This discrepancy may reflect the dual role of TIMP1: it may act as a tumor suppressor by stabilizing the ECM in certain subtypes while promoting immune evasion or resistance in others. IHC data further revealed lower TIMP1 expression in HER2+, ER–, and PR– tumors. These findings highlight the need for subtype-stratified evaluation to accurately assess TIMP1’s biomarker potential.

DAB2, known as a tumor suppressor and adaptor protein, is involved in TGF-β signaling, cytoskeletal regulation, and immune modulation [32]. Consistent with its documented role in various cancers [32], our analysis revealed a significant downregulation of DAB2 in BC tissues, especially in TNBC. This reduction is associated with decreased OS. Promoter hypermethylation may play a crucial role in silencing DAB2 expression. Notably, DAB2 affects the endocytosis of TGF-β receptors and contributes to immune tolerance regulation through Treg cell induction [32,33]. Its loss could lead to unregulated TGF-β signaling, promoting a tumor-friendly microenvironment characterized by immune suppression and metastatic potential [34,35]. The low expression of DAB2 in HER2-, ER-, and TNBC tumors further highlights its potential role as a tumor suppressor in aggressive BC subtypes. Interestingly, our findings also showed a positive correlation between DAB2 and miR-222-3p expression, indicating complex regulatory dynamics that merit further investigation. The downregulation of DAB2, along with its known tumor-suppressive and immunomodulatory functions, suggests a mechanistic contribution to tumor progression in TNBC. Its epigenetic silencing may thus represent a key event in disrupting TGF-β signaling and immune regulation in aggressive BC phenotypes.

Although miR-222-3p is often identified as an oncomiR in various cancers [36], it was found to be downregulated in our BC cohort. This observation contradicts previous reports linking elevated miR-222-3p expression with tumor proliferation, metastasis, and drug resistance. However, low serum levels of miR-222-3p have also been associated with poor prognosis in BC [37]. Our study confirmed its upregulation in BC, particularly in TNBC, which aligns with its known roles in promoting proliferation and chemoresistance. Given its association with EVs and expression in immune cells, miR-222-3p may serve as both a biomarker and an effector in tumor-immune interactions. These findings underscore its significance in the development of targeted therapies. They suggest a complex and potentially subtype-specific association pattern for miR-222-3p in BC. Nevertheless, due to the absence of experimental confirmation, its role remains speculative and requires functional validation. Although traditionally characterized as an oncomiR, its downregulation and correlation with tumor suppressor genes in this cohort may indicate a subtype-specific functional shift. The upregulation of miR-181b-5p and its negative correlation with CXCL14 expressions suggest a potential role in the post-transcriptional suppression of tumor-suppressive pathways. These observations support its investigation as a candidate mediator of immune evasion and metastatic progression in BC.

In our study, we identified promoter hypermethylation as a significant regulatory mechanism, particularly influencing the expression of CD36, CXCL14, and DAB2. The epigenetic silence of these tumor suppressors may facilitate tumor progression, immune evasion, and metastasis. Methylation-based repression patterns can differ among tumor subtypes, suggesting their potential as biomarkers for prognosis and therapeutic stratification [38,39]. The identification of promoter hypermethylation in CD36, CXCL14, and DAB2 underscores the role of epigenetic silencing in BC pathogenesis. These findings also emphasize the potential of methylation profiling as a complementary biomarker strategy for stratifying tumor subtypes.

STRING network analysis revealed that our target genes are involved in critical cancer-related processes, including immune modulation, inflammation, cell migration, and extracellular matrix remodeling. Enrichment analyses showed that TIMP1, ENPEP, CD36, and CXCL14 contribute to sustaining proliferative signaling. Additionally, CXCL14, DAB2, ENPEP, and MARCKS were implicated in evading growth suppressors, while TIMP1, CD36, and CXCL14 were linked to sustained angiogenesis. These processes are hallmarks of cancer and highlight the significance of our findings.

Despite the increasing interest in EVs as mediators of intercellular communication and cancer progression, the field encounters several unresolved challenges. A significant limitation is the heterogeneity of EV populations, which vary considerably in size, content, and cellular origin [40]. Currently, there is no universally accepted set of molecular markers that definitively distinguish EV subtypes, complicating efforts to attribute specific biological functions to EV cargos [41]. Moreover, it is increasingly recognized that molecules from the extracellular environment may adsorb onto EV surfaces, forming a “protein corona.” This phenomenon can influence the detection, function, and downstream interpretation of EV-associated proteins. These factors pose challenges to the consistency and reproducibility of EV studies, including those involving bioinformatics analyses [42,43]. While our study focuses on molecules previously reported to be enriched in EVs, we acknowledge that some may represent dynamic or context-dependent associations rather than core vesicular components. Further experimental validation using EV isolation and proteomic/miRNA profiling will be essential to confirm the vesicular origin and functionality of the candidate biomarkers reported here. These molecules have been implicated in EV-mediated communication that enhances metastatic potential, modulates immune responses, or facilitates resistance to therapy. For instance, TIMP1 has been shown to be packaged into EVs and promote tumor invasiveness via extracellular matrix remodeling [44,45]. Similarly, CD36 and DAB2, when carried by EVs, contribute to lipid uptake and signal transduction within the tumor stroma [34,46]. The miRNAs analyzed in our study, miR-181b-5p and miR-222-3p, are also well-characterized EV components derived from cancer cells and function as post-transcriptional regulators in recipient cells [17,47]. Thus, while our results reflect differential expression patterns and regulatory associations, the implication of these molecules in EV-mediated BC pathogenesis remains speculative and warrants experimental validation through EV isolation and profiling studies.

Overall, our findings suggest that these genes and miRNAs are central to BC pathogenesis and progression, with relevance to aggressive subtypes such as TNBC. Their involvement in essential biological processes and their potential regulation via methylation and miRNAs highlight their utility as diagnostic and therapeutic targets. Future studies should incorporate functional assays, in vitro validation, and stratified patient cohorts to further elucidate the mechanistic pathways underlying these associations and assess their clinical applicability. Importantly, the inclusion of these molecules was grounded in their consistent identification as EV cargos. Their altered expression across BC subtypes, along with prognostic implications, strengthens the hypothesis that they may participate in EV-mediated transfer of oncogenic signals or immune evasion mechanisms. Future studies should investigate whether these molecules directly modulate recipient cell behavior via EV transport.

Limitation

This study encountered several limitations. The data were sourced from specific databases, which restricted access to comprehensive data on BC patients. Consequently, it was not feasible to perform OS analysis individually for each BC subtype, resulting in incomplete data in this regard. Additionally, the disparity in the number of normal adjacent tissue samples compared to tumor tissue samples may introduce bias. Therefore, the findings derived from these databases require validation through further studies involving homogeneous patient groups. While several of the molecules investigated in this study exhibited statistically significant differential expressions between tumor and normal breast tissues, it is important to emphasize that statistical significance does not inherently imply clinical utility as a biomarker. This is particularly relevant in large data sets such as TCGA, where even small differences may yield extremely low p-values. However, as illustrated in Figure 1, the expression distributions of some molecules (e.g., TIMP1, miR-181b-5p, miR-222-3p) show considerable overlap between the cancer and normal groups. This overlap suggests that these markers may yield a non-negligible rate of false positives or false negatives if used diagnostically without additional filtering criteria. For a molecule to be considered a robust biomarker, beyond statistical significance, it should demonstrate a high discriminatory power as measured by sensitivity, specificity, and preferably an AUC (Area Under the Curve) greater than 0.80 in receiver operating characteristic (ROC) analysis. Additionally, effect size (e.g., fold-change or Cohen’s d) and reproducibility across independent cohorts are essential metrics. These considerations highlight the need for further validation and classification modeling, which should be addressed in future studies incorporating cross-validation and independent datasets.

4. Materials and Methods

This is a bioinformatic-based study. In the TCGA-BRCA BC cohort analyzed in this study, which includes RNA-sequencing and clinical data from 1085 primary breast tumor samples and 291 normal breast tissue samples, patients were further stratified by molecular subtype (luminal A, luminal B, HER2-enriched, and triple-negative) according to PAM50 annotations provided in the TCGA database.

We used the following bioinformatic databases while evaluating prognostic features of ENPEP, TIMP1, CD36, MARCKS, DAB2, CXCL14, miR181b-5p and miR-222-3p in BC. The selection of ENPEP, TIMP1, CD36, MARCKS, DAB2, CXCL14, miR-181b-5p, and miR-222-3p was based on a dual strategy combining literature-based prioritization and bioinformatics screening. Initially, we conducted a comprehensive review of EV-related studies and BC literature to identify genes and microRNAs commonly reported to be enriched in EVs and implicated in tumor progression, immune modulation, metastasis, or angiogenesis. These molecules were further screened using public datasets and EV-focused databases to confirm their detectable expression in breast cancer tissues and their availability across multiple omics platforms. Differential expression patterns and previously reported roles in EV-mediated intercellular communication (especially for miR-181b-5p and miR-222-3p) further supported their inclusion in the current study. By focusing on these specific candidates, we aim to ensure both biological relevance and analytical feasibility for downstream functional and clinical correlation analyses.

Expression profile analysis of ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14

The access date: 15 February 2025

The GEPIA2 (http://gepia2.cancer-pku.cn/, accessed on 15 February 2025) web database was used for expression profile analysis. GEPIA2 is an enhanced portion of the web server that offers more detailed gene expression analysis, including transcription-level quantification, analysis of specific cancer subtypes, and the ability for users to upload their own RNA-seq fragments to separate with existing datasets [48]. With the GEPIA2 web server, it is possible to compare the expression levels of genes examined in both tumoral tissue and adjacent healthy tissue. It is also possible to compare the expression levels of the gene or genes examined in cancer subtypes according to the type of cancer. In our BC cohort, we compared the expression levels of ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14 genes with the healthy tissue adjacent (n = 291) to the tumor tissue (n = 1085). We also used the Human Protein Atlas (HPA) (https://www.proteinatlas.org/, accessed on 15 February 2025) to confirm the expression results obtained from GEPIA2. This section of the Human Protein Atlas highlights the expression patterns of genes in human tissues, analyzing both mRNA and protein levels. Protein expression data from 44 normal human tissue types is obtained through antibody-based protein profiling, utilizing conventional and multiplex immunohistochemistry techniques [49].

Targeted expression analysis and gene correlation targeted analysis

The access date: 15 February 2025

We utilized bc-GenExMiner v5.1 web tool for analysis of targeted expression and gene correlation. bc-GenExMiner v5.1 is a statistical analysis tool designed for examining annotated breast cancer transcriptomic data, including DNA microarrays (n = 11,552) and RNA-seq datasets (n = 5023). It enables users to investigate the expression of specific genes of interest in the context of breast cancer [50,51].

The methylation status analysis of ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14

The access date: 15 February 2025

UALCAN (https://ualcan.path.uab.edu/, accessed on 15 February 2025) is a comprehensive, user-friendly, interactive web resource for analyzing cancer transcriptomic data. It aims to facilitate the exploration of TCGA (The Cancer Genome Atlas) data and enable researchers to analyze gene expression profiles and perform in-depth analyses of various cancer types. UALCAN provides insight into the methylation status of gene promoters. Users can compare methylation levels between tumor and normal samples. It allows us to study the correlation between expression levels of two genes. The server is useful for identifying potential gene interactions and pathways [52]. We performed promoter region methylation analyses (tumor tissue n = 793), (normal tissue n = 97) with the UALCAN web server and examined the expression profiles of other input genes in BC, including ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14. Additionally, we also compared the BC subtypes (luminal (n = 566), HER2+ (n = 37), triple negative breast cancer (TNBC) (n = 116)) expression profiles, including ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14, with each other and normal adjacent tissue (n = 114).

The survival analysis of ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14

The access date: 15 February 2025

KM Plotter (Kaplan–Meier Plotter) (https://kmplot.com/analysis/index.php?p=home, accessed on 15 February 2025) is an online tool designed to assess the effect of gene survival on various cancers using clinical data. It enables researchers to conduct meta-analyses of gene expression data to identify potential prognostic biomarkers. Users can generate Kaplan–Meier survival plots to visualize the relationship between gene expression levels and patient survival outcomes. It integrates data from several high-quality sources such as GEO (Gene Expression Omnibus), EGA (European Genome-phenome Archive), and TCGA and combines clinical data with gene expression profiles to provide a comprehensive analysis [53].

The relationship between MicroRNA vs. BC and MicroRNA vs. RNA

The access date: 15 February 2025

The ENCORI Pan-Cancer Analysis Platform (https://rnasysu.com/encori/index.php, accessed on 15 February 2025) is developed to decode Pan-Cancer Networks involving lncRNAs, miRNAs, pseudogenes, snoRNAs, RNA-binding proteins (RBPs), and all protein-coding genes. It achieves this by analyzing their expression profiles across 32 cancer types, utilizing data from approximately 10,000 RNA-seq and 9900 miRNA-seq samples integrated from the TCGA project [54].

Enrichment analysis of ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14

The access date: 15 February 2025

Differentially expressed gene (DEG) analysis was conducted using the Enrichr-KG web tool (https://maayanlab.cloud/enrichr-kg, accessed on 15 February 2025) across several categories: Gene Ontology (GO) biological process, Kyoto Encyclopedia of Genes and Genomes (KEGG), Jensen_DISEASES for exploring disease-gene associations, and DisGeNET for compiling data on genes and variants linked to diseases. The analysis was limited to the top 20 terms, with statistical significance defined as p < 0.05 [55].

Cancer Hallmark Enrichment plot and

The access date: 15 February 2025

The Cancer Hallmark Enrichment tool (https://cancerhallmarks.com/, accessed on 15 February 2025) compiles a consensus list of cancer hallmark genes by integrating 6763 genes from various mapping resources. CancerHallmarks.com leverages the hallmark concept as a powerful organizational framework, linking genes to their corresponding biological functions [56].

Tumor, Normal, and Metastatic (TNM) Plot Analysis

The access date: 15 February 2025

The Tumor, Normal, and Metastatic (TNM) plot database provides users with a real-time comparison of gene expression changes across tumor, normal, and metastatic tissues for all genes, utilizing multiple platforms. The multiple gene analysis page provides an overview of the selected gene set in the selected tissue using gene chip-based data. Accessible without registration at www.tnmplot.com (accessed on 15 February 2025) the portal features three distinct analysis options. Among these, the pan-cancer analysis tool allows for a simultaneous comparison of normal and tumor samples across 22 different tissue types [57].

Prediction of miRNAs associated with genes

The access date: 15 February 2025

We used the MicroRNA Target Prediction Database (miRDB) (https://mirdb.org/, accessed on 15 February 2025) [58] and TargetScan 8.0 (https://www.targetscan.org/vert_80/, accessed on 15 February 2025) to identify and predict the target genes of differentially expressed miRNAs [58,59].

Gene–gene interaction analysis

The access date: 15 February 2025

The STRING database (https://string-db.org/, accessed on 15 February 2025) was employed to elucidate the interactions and potential mechanisms involving ENPEP, TIMP1, CD36, MARCKS, DAB2, and CXCL14 along with other related proteins. STRING systematically predicts protein–protein interactions by considering both functional and physical associations. The data were sourced from various origins, including automated text mining of the scientific literature, computational predictions based on co-expression patterns and conserved genomic contexts, as well as curated interaction experiment databases. These interactions undergo rigorous evaluation, scoring, and hierarchical orthology-based transfer to organisms with limited prior research [60].

5. Conclusions

This integrative bioinformatics study underscores the pivotal roles of EV-associated genes and miRNAs, specifically ENPEP, TIMP1, CD36, MARCKS, DAB2, CXCL14, miR-222-3p, and miR-181b-5p, in the pathogenesis and progression of BC. Our findings indicate that these molecules are not only differentially expressed across molecular subtypes of BC, particularly in aggressive forms such as TNBC, but also engage in key oncogenic processes, including immune modulation, extracellular matrix remodeling, angiogenesis, and metastasis.

We observed significant expression changes in CD36, CXCL14, DAB2, and TIMP1, often associated with promoter methylation or regulatory miRNAs. Notably, the expression levels of TIMP1, DAB2, and CXCL14 were positively correlated with patient OS, whereas miR-181b-5p may serve as a negative prognostic factor by downregulating tumor suppressor genes. The downregulation of CD36 and CXCL14, alongside the upregulation of TIMP1, underscores the complex and context-dependent roles these genes play in BC biology.

Our enrichment and network analyses further corroborate that these genes and miRNAs—selected due to their reported association with extracellular vesicles—are functionally implicated in core hallmarks of cancer, providing a compelling rationale for their consideration as potential biomarkers or therapeutic targets. This is particularly crucial in TNBC, where effective targets remain scarce. Their vesicular context suggests that they may exert systemic effects by modulating recipient cells through EV-mediated pathways, highlighting their translational potential in diagnosis and therapy.

Collectively, these findings lay the groundwork for future functional studies and clinical validation efforts to explore the diagnostic, prognostic, and therapeutic potential of these EV-related molecules in BC. A deeper understanding of their regulatory mechanisms, especially those involving promoter methylation and miRNA interaction, could ultimately contribute to more personalized and effective treatment strategies.