The crystal structure of the ESPL1 protein was obtained from the Protein Data Bank (https://www.rcsb.org/) using PDB ID: 7NJ1. In parallel, the chemical structures of four ligands, namely, Hesperidin (CID: 10621), Quercetin (CID: 5280343), Paclitaxel (CID: 36314), and Docetaxel (CID: 148124) were retrieved from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/). Prior to docking, ligand structure minimization and protein structure preparation (including removal of heteroatoms and water molecules) were carried out using UCSF Chimera (https://www.cgl.ucsf.edu/chimera/). The AutoDock Tools 1.5.7 suite (https://ccsb.scripps.edu/mgltools/downloads/; https://autodock.scripps.edu/download-autodock4/) was then employed to evaluate the molecular binding interactions between ESPL1 and the selected ligands, which included both natural compounds and standard chemotherapeutic agents. For each ligand, ten docking conformations were generated. Among these, the top-scoring conformation, based on lowest binding energy, was selected for further analysis. The ligand–protein interaction profiles were subsequently visualized and interpreted using Discovery Studio Visualizer (http://163.15.166.20:9944/DS/), which enabled the generation of 2D interaction diagrams to identify key binding residues and hydrogen bond interactions.

To compare ESPL1 gene expressions between tumor and normal breast tissue samples, appropriate statistical tests were applied using curated datasets from online platforms. Log-rank tests were used to evaluate differences in survival outcomes, as well as to assess expression heterogeneity and functional enrichment patterns between high and low ESPL1 expression groups. P-values less than 0.05 (P < 0.05) were considered statistically significant. All analyses regarding prognostic performance, expression stratification, and interaction networks were performed using statistical tools integrated within the respective online databases, ensuring methodological reproducibility and robust data interpretation.

First, we analyzed the expression of ESPL1 across all cancers included in the TCGA, including BC, using publicly available web-based tools such as TIMER 2.0, UALCAN, FIREHOSE, and OncoMX. The results revealed that ESPL1 expression was significantly elevated in most cancer types compared to their corresponding normal tissues (Figures 1A–C; Supplementary Table S2). Next, focusing specifically on BC, we used the UALCAN database to compare ESPL1 expression in tumor vs. normal breast tissues. As shown in Figure 2A, ESPL1 was nearly 9-fold upregulated in BC tissues (P < 1.6e-12). This overexpression trend was consistently observed across multiple databases including ENCORI (P = 1.2e–93), GEPIA2 (P < 0.05), OncoDB (P = 2.9e–126), and the TCGA portal (P < 0.05) (Figures 2B–E).

Furthermore, ESPL1 expression was analyzed in an equal number of normal and tumor samples (n = 114 each), and the results confirmed significantly higher ESPL1 expression in tumor tissues (Figure 2F). To assess ESPL1 expression at the transcript level in BC we used the TCGAnalyzeR database, which is based on high-end single-cell RNA sequence transcriptomic data and found a 2.63 log fold change in its expression levels (Figure 2G). Given that certain reference or housekeeping gene, maintain stable expression under most conditions and serve as internal controls in cancer transcriptomic studies, we next compared ESPL1 expression with such genes using the CancerSEA database. This comparison further validated the upregulation of ESPL1 in BC samples (Supplementary Figure S1A). Moreover, to validate our In-silico findings, total RNA was isolated from tissue samples of all the 8 BC patients enrolled in the study. Following RNA quality assessment, only 4 matched pairs of tumor and adjacent normal tissues were deemed suitable for downstream analysis. Subsequent qRT-PCR-based gene expression profiling revealed a significant upregulation of ESPL1 mRNA in breast tumor tissues compared to their adjacent healthy counterparts (Figure 2H), with an overall average increase of 4-fold (Figure 2I).

To further investigate the expression pattern of ESPL1 across various breast cancer (BC) subtypes, we utilized the bc-GenExMiner v5.0 database, stratifying the data based on estrogen receptor (ER), progesterone receptor (PR), and triple-negative breast cancer (TNBC) status. ESPL1 expression was found to be significantly elevated in the more aggressive BC subtypes, including ER–versus ER+, PR–versus PR+, and TNBC versus non-TNBC, as well as in basal-like versus non-basal-like tumors (Supplementary Figures S1B–D). Further subgroup analysis revealed consistently higher ESPL1 expression in TNBC and basal-like patients compared to their non-TNBC and non-basal-like counterparts, respectively (Supplementary Figures S1E–G). We then assessed ESPL1 expression across molecular subtypes, such as, Luminal A, Luminal B, HER2-enriched, and Basal-like and observed notably higher expression in Luminal B and Basal-like subtypes (Supplementary Figures S1H,I). Moreover, using the UALCAN database, we evaluated ESPL1 expression across pathological stages and found it to be significantly upregulated in tumor samples across all stages compared to normal tissue (Supplementary Figure S1J). Collectively, these findings suggest that ESPL1 overexpression is associated with aggressive BC phenotypes, highlighting its potential role as a marker of poor prognosis.

To investigate the functional relevance of ESPL1 in breast cancer (BC), we explored its involvement in key oncogenic processes using the CancerSEA database. The results revealed that ESPL1 expression was significantly associated with critical cancer-related pathways including cell cycle (R = 0.33), proliferation (R = 0.31), metastasis (R = 0.31), and DNA damage response (R = 0.31) (Figures 3A–E). Then, we were interested in finding out the role of ESPL1 in metastasis. Using the TNMplot database based on RNA-Seq data, we compared ESPL1 expression across normal tissue, primary tumors, and metastatic lesions. We observed markedly higher ESPL1 expression in metastatic tissues, with P-value of 4.73e-30, compared to normal tissue. This finding was further validated using the DriverDBv4 platform, which confirmed a significant upregulation of ESPL1 in metastatic samples (P = 0.00078) (Figures 3F,G). To assess the prognostic potential of ESPL1 in BC, we performed Kaplan–Meier survival analyses using the KM Plotter database. Stratification of patients based on high and low ESPL1 expression revealed a strong correlation between elevated ESPL1 levels and poor clinical outcomes. Specifically, high ESPL1 expression was associated with reduced overall survival (OS) (HR = 1.35, 95% CI: 1.12–1.63, P = 0.0019), distant metastasis-free survival (DMFS) (HR = 1.58, 95% CI: 1.35–1.85, P = 6.7e-09), and relapse-free survival (RFS) (HR = 1.33, 95% CI: 1.2–1.48, P = 2.4e-08) (Figures 3H–J; Supplementary Table S3). These findings establish ESPL1 as a potential prognostic biomarker indicative of poor outcomes in BC patients.

Next, we examined the co-expressed genes of ESPL1 using the Enrichr database. Ten genes: AURKB, FOXM1, GTSE1, HJURP, KIF18B, KIF2C, KIFC1, PLK1, RRM2, and TROAP were found to be significantly co-expressed with ESPL1. Validation via the TIMER database demonstrated a strong positive correlation (P < 0.05) between ESPL1 and these genes, with correlation coefficients ranging from R = 0.754 to 0.880 (Table 1; Supplementary Figures S2A–J). Further analysis using the GSCA database showed that all ten co-expressed genes were overexpressed in BC, except KIF18B, HJURP, and GTSE1 (Supplementary Figure S2K). A combined gene set variation analysis (GSVA) revealed a significantly higher expression profile of the ESPL1 co-expression module in tumors compared to healthy tissue (Supplementary Figure S2L). Notably, these co-expressed genes also showed strong enrichment in TNBC and basal-like subtypes of BC, reinforcing their association with aggressive disease phenotypes (Supplementary Figure S2M).

Transcription factors (TFs) are critical regulators of gene expression and play essential roles in controlling diverse cellular processes such as cell cycle progression, DNA replication, apoptosis, and differentiation. Understanding the transcriptional regulation of oncogenic drivers like ESPL1 is therefore crucial for uncovering potential molecular mechanisms underlying cancer aggressiveness and for identifying novel therapeutic targets. To investigate the transcriptional regulators of ESPL1, Enrichr database was used to identify related pathways and found that E2F targets were significantly associated with ESPL1 as shown in Table 2. We focused on the E2F family of transcription factors, as they are known for their central role in controlling the G1/S transition of the cell cycle and promoting the transcription of genes involved in mitotic progression. This family comprises eight members (E2F1–E2F8), which function either as transcriptional activators or repressors depending on cellular context. Furthermore, Using TCGAnalyzeR, we analyzed the expression profiles of E2F genes in BC tumors and found a significant upregulation of E2F1, E2F2, E2F7, and E2F8, with log2 fold-change values of 2.09, 2.41, 2.94, and 3.49, respectively (Figures 4A–H). In contrast, E2F3, E2F4, E2F5, and E2F6 showed minimal or no significant change in expression (0.88, −0.11, 0.94, and 0.14, respectively). A combined expression plot (Figure 4I) demonstrated that E2F7 and E2F8 were in closest proximity to ESPL1, suggesting a potential regulatory axis between these TFs and ESPL1, which itself was upregulated with a log2 fold-change of 2.63.

To corroborate this potential relationship, we performed correlation analyses using the ENCORI and TIMER databases. Interestingly, ENCORI analysis confirmed that all E2Fs were significantly correlated with ESPL1 expression in BC (Figures 4J–Q). The TIMER analysis further supported these findings and highlighted only E2F7 (R = 0.741) and E2F8 (R = 0.814) as having the strongest positive correlations with ESPL1 among the E2F members (Figures 4R–Y). Together, these results suggest that E2F7 and E2F8 may act as critical upstream transcriptional regulators of ESPL1, potentially promoting its overexpression in aggressive subtypes of breast cancer.

While the transcriptional regulation of ESPL1 in breast cancer (BC) is increasingly understood, the mechanisms underlying its aberrant upregulation remain incompletely elucidated. Given the prominent role of microRNAs (miRNAs) in post-transcriptional gene silencing, we sought to identify miRNAs potentially involved in the negative regulation of ESPL1 expression in BC. We utilized the miRNet platform to construct a comprehensive miRNA–ESPL1 interaction network, identifying the top 11 miRNAs predicted to regulate ESPL1 (Supplementary Figure S3A). Further, the expression data from the UALCAN database revealed that four miRNAs, hsa-miR-10a-5p, hsa-let-7b-5p, hsa-miR-214-3p, and hsa-miR-1-3p were significantly downregulated in tumor tissues from BC patients (Supplementary Table S4). Among these, hsa-let-7b-5p and hsa-miR-10a-5p exhibited a strong negative correlation with ESPL1 expression, analyzed using the ENCORI database. To assess the prognostic value of these miRNAs, we analyzed survival plots using the CancerMIRNome database. Although hsa-miR-10a-5p showed a tumor suppressive expression pattern, it was not significantly associated with survival outcomes in BC patients (Supplementary Figure S3B). Intrestingly, hsa-let-7b-5p emerged as a key tumor suppressor miRNA, not only did it negatively correlate with ESPL1 expression, but it also demonstrated strong associations with transcriptional regulators, particularly E2F8, in the ENCORI database (Supplementary Figures S3C,D), suggesting a regulatory axis involving hsa-let-7b-5p, E2F8, and ESPL1.

We further validated the negative correlation between hsa-let-7b-5p and ESPL1 using ENCORI, TACCO, and CancerMIRNome datasets, all of which showed statistically significant negative relationships (Figures 5A–C). Additionally, Kaplan–Meier survival analyses using CancerMIRNome (HR = 0.58, P = 8.37e-04), KM Plotter (HR = 0.68, P = 0.00014), and ENCORI (HR = 0.56, P = 0.00014) consistently indicated that low expression of hsa-let-7b-5p is significantly associated with poor prognosis in BC patients (Figures 5D–F). Expression profiling of hsa-let-7b-5p across tumor and normal tissues in BC patients revealed significant downregulation in tumors using CancerMIRNome (P = 2.45e-07) and UALCAN (P = 4.6e-07) as shown in Figures 5G,H. Further stratification by clinical parameters assessed using UALCAN demonstrated that hsa-let-7b-5p expression was reduced in advanced pathological stages, nodal metastasis, and across both male and female patients (Figures 5I–K). Analysis through ExplORRnet confirmed its downregulation in metastatic BC patients (Figure 5L), while UALCAN analysis indicated lower expression in triple-negative breast cancer (TNBC) cases (Figure 5M), reinforcing its relevance in aggressive disease subtypes. Together, these findings implicate hsa-let-7b-5p as a critical post-transcriptional regulator of ESPL1 and highlight a potential ceRNA network involving hsa-let-7b-5p, E2F8, and ESPL1 that may drive tumor progression and poor prognosis in breast cancer.

In cellular systems, lncRNAs, miRNAs, and mRNAs can form ceRNA networks, collaboratively modulating gene expression and cellular processes. To identify lncRNAs that might modulate ESPL1 expression via interaction with hsa-let-7b-5p, we used the Enrichr database, and several candidate lncRNAs were identified to be associated with ESPL1, including LINC00618, SGO1-AS1, DEPDC1-AS1, LIX1L-AS1, and LINC01775 (Table 3). Further, using the UALCAN database, we found that DDX11-AS1, TMPO-AS1, DEPDC1-AS1, and CSRP3-AS1 were upregulated in BC tissues compared to normal controls. Among them, ENCORI analysis revealed that TMPO-AS1 had the strongest positive correlation with ESPL1 (R = 0.628, P = 2.39e-122) (Supplementary Table S5; Figure 6A). Notably, TMPO-AS1 showed a significant negative correlation with hsa-let-7b-5p (R = −0.204, P = 1.12e-11), suggesting a ceRNA-based regulatory relationship (Figure 6B). Given our earlier findings that E2F8 might regulate ESPL1, we further assessed its correlation with TMPO-AS1. ENCORI analysis showed a strong positive association between TMPO-AS1 and E2F8 (R = 0.514, P = 1.85e-75) (Figure 6C), further supporting their interaction within a regulatory axis. We validated the expression levels of TMPO-AS1 using multiple datasets. lncRNADisease v3.0 was used to assess the pan-cancer expression profile of TMPO-AS1 and it was found to be overexpressed amongst several malignancies with Breast Neoplasms at the top (Supplementary Table S6). Further, ENCORI, OncoDB, UALCAN, and TCGAnalyzeR databases consistently demonstrated significantly higher expression of TMPO-AS1 in BC tumors compared to normal tissues (Figures 6D–G). Moreover, UALCAN data showed TMPO-AS1 expression was elevated across different pathological stages and was also upregulated in TNBC compared to normal controls (Figures 6H,I). Additionally, using miRNet database and a direct TMPO-AS1/hsa-let-7b-5p/ESPL1 axis was revealed supporting a ceRNA model whereby TMPO-AS1 sponges hsa-let-7b-5p, thereby releasing repression on ESPL1 (Supplementary Figure S3E). To further illustrate these regulatory interactions, a lncRNA–miRNA–mRNA network was constructed using Cytoscape (3.10.3). The resulting visualization highlighted a potential regulatory cascade involving TMPO-AS1, hsa-let-7b-5p, E2F8, and ESPL1, with directional arrows indicating the predicted regulatory flow (Supplementary Figure S3F).

Furthermore, to assess this regulatory potential of hsa-let-7b-5p on key oncogenic drivers, interaction analyses were conducted with ESPL1 (NM_012291.5), E2F8 (NM_001256371), and the lncRNA TMPO-AS1 (NR_027157.1) (Table 4). miRWalk predicted strong binding affinities of hsa-let-7b-5p to ESPL1 and E2F8, with minimum free energy (MFE) values of −24.4 kcal/mol and −22.2 kcal/mol, respectively. RNA22v2 folding energy analysis corroborated these findings, showing energetically favorable heteroduplex formations with ESPL1 (−22.0 kcal/mol), E2F8 (−14.4 kcal/mol), and TMPO-AS1 (−12.4 kcal/mol). Additionally, the heteroduplex modeling revealed specific base pairing between hsa-let-7b-5p and each transcript, reinforcing the likelihood of direct post-transcriptional regulation. Collectively, these data suggest that hsa-let-7b-5p may suppress oncogenic mRNAs ESPL1 and E2F8, key cell cycle regulators, while its interaction with the oncogenic lncRNA TMPO-AS1 indicates a ceRNA mechanism. TMPO-AS1 likely sequesters hsa-let-7b-5p, relieving repression of ESPL1 and E2F8, thus, upregulation of hsa-let-7b-5p could disrupt this ceRNA network, attenuating cancer progression by simultaneously targeting oncogenic transcripts and their sponge lncRNA.

To develop effective therapeutic strategies for ER-/PR-breast cancers, it is essential to understand the molecular mechanisms underlying ESPL1 overexpression. Our analysis identified a significant regulatory axis involving TMPO-AS1, which showed a strong positive correlation with the proliferation marker MKI67 (Supplementary Table S6). MKI67 is widely used in pathological evaluations and is tightly associated with tumor cell proliferation and growth (Li et al., 2014). To establish the relation of MKI67 with our target genes ccorrelation analyses were performed using ENCORI database which revealed consistent co-expression patterns between MKI67 and ESPL1 (R = 0.827), E2F8 (R = 0.816), hsa-let-7b-5p (R = −0.153), and TMPO-AS1 (R = 0.540) in breast cancer and the same correlations were corroborated (Supplementary Figures S4A–G). Furthermore, we explored the relationship between these genes and the expression of estrogen receptors (ESR1) and progesterone receptors (PGR), both of which are key biomarkers and therapeutic targets in BC.

Using the ENCORI database, ESR1 and PGR exhibited significant negative correlations with ESPL1 (R = −0.200, P = 2.12e-11; R = −0.258, P = 4.70e-22) (Figures 7A,B), while E2F8 similarly negatively correlated with ESR1 (R = −0.266, P = 2.56e-19) and PGR (R = −0.329, P = 2.40e-29) (Figures 7C,D). In contrast, hsa-let-7b-5p positively correlated with both ESR1 (R = 0.267, P = 3.60e-19) and PGR (R = 0.274, P = 3.64e-20) (Figures 7E,F). TMPO-AS1 was inversely correlated with ESR1 (R = −0.167, P = 2.41e-08) and PGR (R = −0.232, P = 5.52e-15) (Figures 7G,H). MKI67 also showed strong negative correlations with ESR1 (R = −0.318, P = 2.38e-27) and PGR (R = −0.335, P = 2.30e-30) (Figures 7I,J). We did the cross-validation of the correlation values using OncoDB and found the correlation patterns (Supplementary Figures S4H–O). Consistent with these observations, expression profiling demonstrated that ESPL1, E2F8, and MKI67 were overexpressed across BC stages, while ESR1 and PGR were downregulated, particularly in more aggressive or hormone receptor-negative subtypes. Notably, ESPL1 expression patterns closely mirrored those of MKI67, as shown in Supplementary Figure S4P. To our knowledge, this is the first in-silico study to report ESPL1 overexpression mimicking MKI67 expression in breast cancer. These findings suggest that elevated ESPL1 expression could serve as a proxy for tumor cell proliferation, particularly in ER-/PR-subtypes and intermediate stages of disease. Thus, ESPL1 emerges as a promising marker of aggressive tumor growth and potential therapeutic vulnerability in hormone receptor-negative breast cancers.

Molecular docking studies were conducted to evaluate the binding affinities of natural and chemotherapeutic compounds with ESPL1. Among the tested ligands, Hesperidin exhibited the binding affinity with ESPL1 at −10.8 kcal/mol and Quercetin demonstrated the binding affinity of −7.8 kcal/mol indicating a strong and stable interaction. In comparison, Paclitaxel and Docetaxel showed binding affinities of −8.4 kcal/mol and −8.2 kcal/mol, respectively (Table 5). The highest docking score of Hesperidin suggests that it may serve as a potent inhibitor of ESPL1. These interactions were further supported by 2D molecular interaction diagrams generated in Discovery Studio Visualizer (Figures 8A–D). The combination of multiple hydrogen bonds, strong van der Waals contacts, and favorable π-alkyl interactions with key active site residues such as GLU432, GLN406, and HIS397 implies that Hesperidin not only fits tightly within the ESPL1 binding site but may also effectively inhibit its activity. In contrast, chemotherapeutic agents demonstrated fewer stabilizing interactions within the binding pocket. Notably, the extensive hydrogen bonding and π-π stacking of Hesperidin suggest a more stable and specific binding to ESPL1, highlighting its potential as a promising inhibitor.