We obtained the TCGA BRCA dataset’s expression matrix using the TCGAbiolinks R package (10), which includes count sequencing data from 1222 breast invasive carcinoma samples. These consist of 1109 breast cancer samples and 113 adjacent non-cancerous tissue samples. The data was normalized to FPKM format for analysis. Additionally, we computed the TIDE Immunocore for patient samples using the TIDE (Tumor Immune Dysfunction and Exclusion) platform available at tide.dfci.harvard.edu (11). We accessed BRCA-related datasets (GSE7904, GSE29431, and GSE42568) from the Gene Expression Omnibus (GEO) database via the R package GEOquery (12–14).

We acquired 50 immune checkpoint genes from cited literature (15), precisely designated in Supplementary Table S2. We retrieved 32 HLA family genes from the GeneCards database using “Human Leukocyte Antigen” as the search term. Subsequently, we retained genes whose Gene Symbol commenced with “HLA,” as detailed in Tables 1 and 2. Somatic mutation data, including single-nucleotide (16) variations (SNVs), were obtained from the TCGA-BRCA dataset via the TCGA website. We visualized this data using the R package maftools. Copy Number Variation (CNV) data from the same dataset were procured using the TCGAbiolinks R package and analyzed with GISTIC 2.0 (17). We employed the cBioPortal database (18),accessed through http://www.cbioportal.org, to obtain data on microsatellite instability (MSI) and tumor mutation burden (TMB). The TIDE score was calculated using the TIDE algorithm from the TIDE website. Given patient privacy and data-sharing policies for the dataset, we used publicly available data from databases such as TCGA (The Cancer Genome Atlas) and GEO (Gene Expression Omnibus). To facilitate further access for researchers, we provide detailed instructions for data download and preprocessing on jianguoyun, along with specific guidelines on obtaining these data from public databases.

In this study, we employed the following criteria for sample inclusion and exclusion:

Inclusion Criteria:

Patients diagnosed with breast invasive carcinoma (BRCA), Availability of surgical tissue samples for analysis, Written informed consent was obtained from each patient.

Exclusion Criteria:

Patients with incomplete clinical data, Prior treatments that may influence tumor characteristics, and the Presence of other concurrent malignancies.

Normalization was performed using the limma R package, and dataset distribution was visualized via boxplots. Differentially expressed genes (DEGs) with logFC < 0 and p-value < 0.05 were identified and compared with Hedgehog-related genes (HRGs) using Venn diagrams, volcano plots, and differential rank maps generated with ggplot2.

First, we standardized the breast cancer datasets GSE7904, GSE29431 and GSE42568 using the R package limma package (Figures 1A–F). According to the figure, the dataset GSE7904 (Figures 1A, B), The expression profile data of GSE29431 (Figures 1C, D) and GSE42568 (Figures 1E, F) tended to be consistent among the samples after normalization.

To delineate the tumor immune microenvironment at single-cell resolution, we utilized the curated dataset GSE114725 from the GEO database, which comprises 23,688 single cells profiled using the 10x Genomics Chromium platform. Data preprocessing, quality control, and normalization were conducted using Cell Ranger and the Seurat R package, adhering to established best practices. Subsequent clustering analysis facilitated the identification of major cellular subsetsrcationntt tumor epithelial cells, T cells, B cells, macrophages, endothelial cells, and fibroblasts,onntt on canonical markers (e.g., EPCAM, CD3D, MS4A1, CD68, PECAM1, ACTA2). Cell type annotations were further validated against reference datasets such as the Human Primary Cell Atlas to ensure accuracy and reproducibility. Within this refined cellular landscape, we specifically examined the expression patterns of PSMB8 and BCL2 across immune cell subtypes relevant to breast cancer, leveraging the analytical framework provided by the IMMUcan database (18).We performed functional enrichment analysis using Gene Ontology (GO) and pathway enrichment analysis with the Kyoto Encyclopedia of Genes and Genomes (KEGG) for insights into gene functions at varying expression levels. Gene Ontology (GO) facilitates large-scale functional enrichment investigations, encompassing biological processes (BP), molecular function (MF), and cellular components (CC) (19).The KEGG database contains valuable information on genomes, biological pathways, diseases, and medications (20, 21). We ensured statistical significance using clusterProfiler (an R package) for GO annotation analysis of Hedgehog-Related Differentially Expressed Genes (HRDEGs), setting a significance threshold at p < 0.05 and a False Discovery Rate (FDR) value (q-value) below 0.05 (22).

We performed Gene Set Enrichment Analysis (GSEA) using the clusterProfiler software, prioritizing genes based on log fold change (logFC) and identifying enriched gene sets from the MSigDB database using the “c2.cp.v7.2.symbols.gmt” set. Significant enrichment was defined by a p-value < 0.05 and false discovery rate (FDR) < 0.25, which strikes a balance between biological discovery and error control (23). Additionally, we applied Gene Set Variation Analysis (GSVA), an unsupervised method independent of parameter settings. GSVA was conducted using the “h.all.v7.4.symbols.gmt” gene set from MSigDB on the TCGA-BRCA dataset. This allowed us to assess variations in functional enrichment between normal and BRCA patient samples (24).

To construct a prognostic model for Hedgehog-Related Genes (HRGs) in the Hh pathway and evaluate their predictive potential, we used univariate Cox regression analysis on the BRCA breast cancer dataset, combining Overall Survival (OS) and time. Genes with a p-value < 0.01 were further analyzed using multivariate Cox regression to develop the model and calculate RiskScore. Adopting a stricter p-value standard can help narrow the range of candidate genes or features, reduce the risk of false positives, and provide a more reliable basis for experimental validation and clinical application. The RiskScore, also referred to as the Hedgehog score, was constructed based on nine prognosis-related genes: TUBA3D, TUBA1C, PSME2, RPGRIP1L, PSMB8, BCL2, GLI1, GPC5, and BMP5. These genes were selected through multifactorial Cox regression analysis without any reduction in their number.

We used the STRING database to construct a protein-protein interaction network for differentially expressed genes, setting a minimum interaction score of 0.150 (25). The cytoHubba plugin, employing DEGREE Correlation, MNC, MCC, EPC, and DMNC algorithms, identified the top six HRDEGs (26). A correlation heatmap visualized connections between immune checkpoint genes, HLA family genes, and key genes. Access the starBase database (version 3.0) at https://rnasysu.com/chipbase3/index.php and the CHIPBase database (version 3.0) at https://discover.nci.nih.gov/cellminer/home.do (27, 28).We used the hTFtarget database to identify transcription factors interacting with essential genes (29).

We analyzed CNV, somatic mutations (SNPs), and functional similarity using TCGAbiolinks, GISTIC 2.0, VarScan, and maftools (30). The functional similarity was assessed using the GOSemSim R package.

We assessed the prevalence and abundance of immune-infiltrating cells using single-sample Gene Set Enrichment Analysis (ssGSEA) through the GSVA package in R. Enrichment scores were calculated to quantify immune cell infiltration across various categories (31). Additionally, following previously established guidelines (32), we estimated the presence of eight immune cell types and two non-immune stromal cell groups within the tumor microenvironment. Spearman’s correlation analysis evaluated the statistical relationship between MCPcounter abundance estimates, antigen gene expression, and immune infiltration, with statistical significance set at p-value < 0.05.

This section details an evaluation designed to determine how the expression stability of key genes influences the variation of other variables within our dataset. The GDSC database (cancerRxgene.org) contains extensive information on molecular markers associated with drug sensitivity and response in cancer cells. The CCLE repository (https://portals.broadinstitute.org/ccle/about) Genomic alterations in cancer significantly influence treatment response and often serve as biomarkers for drug therapy (33). The GDSC database (cancerRxgene.org) contains extensive information on molecular markers associated with drug sensitivity and response in cancer cells. The CCLE repository. The CellMiner database (https://discover.nci.nih.gov/cellminer/home) offers experimental data and vital genes for the initial dataset (34). These tools facilitate the comprehensive investigation of specific genes, cells, or test data, allowing for the storage, retrieval, and download of molecular profiling data for various cancer cell types. We analyzed drug sensitivity for essential genes using expression levels and drug information from the GDSC, CCLE, and CellMiner databases.

The Consensus Clustering algorithm uses resampling to determine the number of subgroups and evaluate clustering reliability. This technique aids in understanding cluster stability and parameter determination (35). To categorize distinct BRCA subtypes in the TCGA-BRCA dataset, we used the ConsensusClusterPlus R package, generating clusters ranging from 2 to 8 and conducting 50 iterations while randomly selecting 80% of the samples. The clustering algorithm was “km” with the “euclidean” distance metric.

Using ssGSEA, we determined the relative prevalence of individual genes in each dataset sample. Hs values for TCGA-BRCA samples were computed using the ssGSEA algorithm from the GSVA R package. We divided the BRCA samples into Hs_high and Hs-low groups based on the median Hs. Tumor purity was evaluated using the ESTIMATE algorithm, calculating immune and stromal scores from the expression data.

To investigate the predictive significance of Hedgehog scores (Hs) about immunotherapy, we obtained Hs data from the TCIA database (https://www.tcia.at/home) (36). Using the ggplot function in R, we generated a boxplot illustrating the differences between Hs high and low clusters of TCGA-BRCA patients across various Immune Phenotype Scores (IPS). We conducted a multivariate Cox regression analysis using TCGA-BRCA expression profile data, incorporating Hs and clinicopathological factors, to create a clinical prediction nomogram with the rms package in R. Calibration charts compared projected values of the nomogram with actual survival rates.

We collected 271 tissue samples from individuals who underwent surgical removal, with ethics committee approval from Nantong University (protocol code: 2020-L125) and informed consent from participants. RNA extraction from paraffin-embedded tissues was performed using a Thermo Fisher Scientific kit. qRT-PCR assessed mRNA levels of PSMB8, BCL2, BMP5, and PSME2, using GAPDH as the internal control. Primer sequences were as follows:

All data processing and analyses were conducted using R software (version 4.1.2). We first applied the Shapiro-Wilk test for each continuous variable to determine normal distribution, with a P-value below 0.05 indicating a non-normal distribution. An independent t-test assessed the statistical significance of constant variables for normally distributed data. At the same time, the Mann–Whitney U test (Wilcoxon rank-sum test) was used for non-normally distributed data. Categorical variables were analyzed using the chi-square or Fisher’s exact test, and the Kruskal–Wallis test compared three or more groups. Kaplan–Meier survival curves, generated using the R survival package, illustrated survival variations, with the Log-rank test determining significance. Statistical significance was defined as a P-value below 0.05; all P-values were two-sided.

Using the limma package, we analyzed the TCGA-BRCA dataset and identified 31,255 differentially expressed genes (DEGs) (Table 3), with 12,683 upregulated and 18,542 downregulated in breast cancer (BRCA) samples (Figure 2A).

To isolate Hedgehog-related DEGs (HRDEGs), we filtered for DEGs and Hedgehog-related genes (HRGs), identifying 207 HRDEGs (Figure 2B). A heatmap shows the expression of these HRDEGs across BRCA and normal samples (Figure 2C).

We performed enrichment analysis on 207 HRDEGs using GO and KEGG, with significant results (P < 0.05, FDR < 0.05) visualized in bubble charts, circular network diagrams, and bar graphs (Supplementary Figure S12). Key GO enrichments included Wnt, Notch, transcription regulation pathways, and components like the proteasome and peptidase complexes. KEGG analysis highlighted JAK-STAT, TGF-beta, Hedgehog, Wnt, and apoptosis pathways, providing insights into gene interactions and breast cancer progression.

In the context of GO and KEGG enrichment, significant items were selected based on their p-values and false discovery rate (FDR) values, with a threshold set at < 0.05. This analysis helps us understand the functional relevance and potential pathways associated with the Highly Differentially Expressed Genes, shedding light on their role in breast cancer development and progression.

We conducted Gene Set Enrichment Analysis (GSEA) and Gene Set Variation Analysis (GSVA) on the TCGA-BRCA dataset to explore the relationship between gene expression levels and associated biological processes. Gene Set Enrichment Analysis (GSEA): GSEA identified significant enrichments in the TP53, MAPK, Notch, Wnt, and Hedgehog pathways (Supplementary Figures S13A–F), suggesting their potential association with breast cancer development and progression. The enrichment results indicate the statistical clustering tendencies of the relevant gene sets; however, they do not directly imply that these pathways are necessarily in a higher activation state at the functional level. Detailed results are summarized in Table 4.

Gene Set Variation Analysis (GSVA): GSVA provided insights into differences between breast cancer samples (BRCA group) and standard samples (Normal group). This analysis identified significant variations in gene sets, particularly in the Notch and Wnt pathways (Supplementary Figure S12G). Further details are provided in, Table 5.

Overall, GSEA and GSVA elucidate gene expression’s impact on breast cancer and identify critical pathways involved in its progression.

To develop a predictive model for highly differentially expressed genes (HRDEGs) associated with the Hedgehog (Hh) pathway, we analyzed overall survival (OS) data from the TCGA-BRCA breast cancer dataset. Univariate Cox regression identified nine significant HRDEGs (TUBA3D, TUBA1C, PSME2, RPGRIP1L, PSMB8, BCL2, GLI1, GPC5, and BMP5) with p-values < 0.01 (Figure 3A).

These genes were incorporated into a multivariate Cox regression model, creating a risk score to categorize patients into high- and low-risk groups (Table 6). The model’s predictive ability was assessed using a nomogram (Figure 3B), indicating that PSME2 had the highest diagnostic efficacy while BMP5 had the lowest.

The nomogram was adapted for 1-, 3-, and 5-year predictions (Figures 3D–F), with calibration curves confirming alignment with actual survival outcomes. Decision Curve Analysis (DCA) demonstrated that the 1-year model offered greater clinical utility than the 3- and 5-year models (Figures 3G–I).

We also generated OS Kaplan-Meier (KM) curves using the risk scores from the TCGA-BRCA dataset (Supplementary Figure S2A). This analysis revealed that higher risk scores significantly predicted poorer OS (P < 0.001), as illustrated in Supplementary Figure S2B. Time-dependent AUC curves for 1, 3, and 5 years showed the RiskScore’s effectiveness, with AUCs of 0.7–0.9 for 1 and 3 years, while the 5-year prediction had lower accuracy (AUC 0.5–0.7).

In summary, the time-sensitive ROC curves (1, 3, and 5 years) for the RiskScore (Figure 3) demonstrated moderate accuracy for the first two years and relatively poor prediction for the fifth year.

The interaction relationships among the genes were visually represented using Cytoscape software (Figure 1A). We employed the cytoHubba plug-in to analyze five algorithms (MCC, DMNC, EPC, Degree, and Closeness), identifying the top four HRDEGs—PSMB8, BCL2, BMP5, and PSME2—as essential genes (Figures 1B, C). Using the GSCA online tool, we examined the expression patterns of these genes across 14 cancer types, which required a minimum of three tumor samples alongside their normal counterparts (Figure 1D and Table 1). The results indicated a consistent decrease in BCL2 and BMP5 expression in most tumor samples, while PSMB8 and PSME2 showed increased expression.

We utilized the ENCORI database to predict interactions between RNA-binding proteins (RBPs) and four essential genes: PSMB8, BCL2, BMP5, and PSME2. The interaction network shows connections between BCL2, BMP5, and PSME2 with 74 RBPs, totaling 94 connections (Figure 4A; Table 2).

We also analyzed the CHIPBase (version 3.0) and target databases to identify transcription factors (TFs) interacting with these genes. The network reveals interactions between PSMB8, BCL2, and PSME2 with 124 TFs, resulting in 211 connections (Figure 4B; Supplementary Table 5).

We used the IMMUcan database to analyze single-cell RNA data and examine the expression of PSMB8, BCL2, and PSME2 in the breast cancer immune microenvironment. Figure 5A shows the annotation of various immune cell types in breast cancer. PSMB8, BCL2, and PSME2 were expressed in tumor cells, stromal cells, and immune cell subsets (Figures 5B, C). PSMB8 was primarily expressed in macrophages, BCL2 in B cells, and PSME2 in monocytes (Figures 5D, E).

The TCGA-BRCA dataset demonstrated distinct accuracies for BMP5, PSMB8, PSME2, and BCL2, with varying AUC values across categories.

This figure compares the expression of key genes, including BMP5, PSMB8, PSME2, and BCL2, between BRCA and normal groups, where red represents BRCA and blue represents normal. A heat map visualizes the correlation among these significant genes within the TCGA-BRCA dataset. ROC curve analysis was conducted for BMP5, PSMB8, PSME2, and BCL2, assessing their diagnostic performance. The area under the curve (AUC) values indicates varying levels of precision: values above 0.9 suggest high precision, values between 0.7 and 0.9 indicate moderate precision, and values between 0.5 and 0.7 reflect relatively low precision. Additionally, PCR results reveal that PSME2 and PSMB8 are highly expressed in cancer tissues, while BMP5 and BCL2 are more highly expressed in normal tissues (Figures 6A–C).

We analyzed the expression of four essential genes (PSMB8, BCL2, BMP5, PSME2) using the GSE7904 dataset (Figure 4A) and examined their differential expression between BRCA and normal groups in GSE29431 (Figure 4B) and GSE42568 (Figure 4C). After excluding genes not present in these datasets, we focused on PSMB8, BCL2, and BMP5. In GSE7904, PSMB8 and BCL2 showed significant differences between BRCA and normal groups (p < 0.05). In GSE29431, BMP5 (p < 0.001) and BCL2 (p < 0.05) exhibited significant differences, and in GSE42568, both BMP5 and BCL2 (p < 0.001) and PSMB8 (p < 0.01) were significantly different. ROC curves for the three genes in these datasets were plotted, excluding results with AUC values below 0.6. In GSE7904, PSMB8 and BCL2 had low accuracy (AUC 0.5–0.7). BMP5 and BCL2 in GSE29431 and BMP5, BCL2, and PSMB8 in GSE42568 demonstrated moderate accuracy (AUC 0.7–0.9).

The TCGA-BRCA breast cancer dataset exhibits separate Kaplan–Meier survival curves (KM) for four notable genes (PSMB8, BCL2, BMP5, and PSME2). Statistical significance was set below 0.05, indicating that the four crucial genes in the TCGA-BRCA breast cancer dataset were statistically significant (Figures 7A–D).

Figure 8 depicts the prognostic analysis of key genes using Kaplan–Meier (KM) curves for overall survival (OS). The KM curves for PSME2 (A), PSMB8 (B), BCL2 (C), and BMP5 (D) were examined for their prognostic value. In these curves, the Low-expression group of breast cancer patient samples is represented in blue, while the High-expression group is shown in red. Notable distinctions are indicated by P < 0.05, while differences with a high significance level are marked as P < 0.01.

We examined the clinical correlation of four key genes (PSMB8, BCL2, BMP5, and PSME2) with breast cancer prognosis. A violin plot (Figure 9A) showed gene expression variations by age (≤65 and >65) in the TCGA-BRCA dataset, with BMP5 demonstrating a significant difference (p < 0.001). We also analyzed the relationship between gene expression and clinicopathological factors, including disease-specific survival (DSS) and progression-free interval (PFI) (Supplementary Figures 8B-D).

Overall survival (OS) analysis revealed significant differences (p < 0.01), particularly for PSMB8, which differed notably between Alive and Dead groups (p < 0.05). BCL2 showed significant variation in DSS (p < 0.01), while PSME2 significantly influenced both DSS (p < 0.05) and PFI (p < 0.01). Additionally, PSMB8 exhibited significant variations in PFI between the Alive and Dead groups (p < 0.05).

A comprehensive analysis of copy number variations (CNVs) in the TCGA-BRCA cohort, utilizing GISTIC 2.0, identified four key genesifiediv BCL2, BMP5, and PSME2ified recurrent genomic alterations. Among these, PSMB8 exhibited the highest frequency of amplifications, while BCL2 demonstrated the most frequent deletions (Supplementary Figure S5A). Notably, PSME2 displayed a substantial number of both amplification and deletion events (Supplementary Figure S5A), indicating underlying genomic instability at this locus. Somatic mutation profiling conducted using the maftools R package further delineated nine major mutation categories, with missense mutations being predominant. Single nucleotide variants (SNVs), particularly C-to-T transitions, were the most frequently observed (Supplementary Figure S5B). Functional similarity (Friends) analysis underscored a potentially central role for PSME2 in breast cancer pathogenesis, as it approached the threshold of statistical significance (Supplementary Figure S5C). Consistent with this finding, subsequent immune infiltration and survival analyses revealed significant associations between PSME2 alterations and the levels of various tumor-infiltrating immune cells, as well as overall survival (OS). These results suggest that CNV aberrations in PSME2 may contribute to breast tumorigenesis by modulating gene expression and shaping the tumor immune microenvironment.

We applied the ssGSEA method to assess the relationship between the expression profiles of 28 immune cells in the TCGA-BRCA dataset across different categories (Normal/BRCA). Boxplots (Supplementary Figure 6A) revealed significant differences (p < 0.05) in the expression of 21 immune cell types, including activated B cells, CD4 T cells, CD8 T cells, dendritic cells, and natural killer cells. A heatmap (Supplementary Figure 6B) showed the correlation between the expression of four key genes (PSMB8, BCL2, BMP5, PSME2) and immune cell infiltration, with significant associations (p < 0.05). PSMB8 correlated positively with CD56dim NK cells and PSME2 with Type 1 Th cells, while negative correlations were observed for PSME2 with plasmacytoid dendritic cells and BCL2 with CD56bright NK cells. Using the MCPCounter algorithm, we identified 10 distinct immune cell categories in the TCGA-BRCA dataset, with significant correlations, such as PSMB8 with T cells/CD8 T cells and inverse relationships between PSME2 and endothelial cells and BCL2 with monocytic lineage cells (Supplementary Figure 6C).

This section presents a drug sensitivity analysis using data from the CellMiner, GDSC, and CCLE databases. We developed a ridge regression model using the prophetic algorithm to predict gene responsiveness to anticancer drugs based on IC50 measurements by integrating mRNA expression profiles with drug activity data for key genes. In the CellMiner database, BCL2 was associated with 6-chloro-2-methoxy-9-[(5-piperidinopentyl)amino]acridine and enhydrin A, while PSMB8 correlated with gw772405x and negatively with zimelidine hydrochloride (Supplementary Figure 7A). The GDSC database showed positive correlations between BCL2, Tanespimycin, and BMP5 with T0901317, though BCL2 had negative correlations with most small drug molecules (Supplementary Figure 7B). In the CCLE database, BCL2 was inversely correlated with L-685458 and panobinostat (Supplementary Figure 7C).

We identified 29 immune checkpoint genes by intersecting a set of 50 known genes with the TCGA-BRCA dataset. We compared the expression levels of four key genes (PSMB8, BCL2, BMP5, and PSME2) with these 29 immune checkpoint genes (Supplementary Figure 8A). Notably, PSMB8 showed a significant correlation with CD27 (R = 0.502, p-value = 0).

In addition, we integrated 19 genes from the HLA family by merging 32 HLA genes with the TCGA-BRCA dataset. Our correlation analysis revealed specific relationships between the four main genes and HLA genes (Supplementary Figure S8B). For instance, PSMB8 strongly correlated with HLA-A (R = 0.812, p-value = 0), while PSME2 moderately correlated with HLA-E and HLA-A (R = 0.514).

Using the ‘ConsensusClusterPlus’ R package, we analyzed the expression of four key genes (PSMB8, BCL2, BMP5, PSME2) in BRCA samples from the TCGA-BRCA dataset. Consensus clustering identified two BRCA subtypes: Cluster 1 (894 samples) and Cluster 2 (215 samples) (Figure 8A). Principal Component Analysis (PCA) revealed significant differences between these subtypes (Figure 8B). The delta chart and CDF diagram confirmed optimal clustering with two clusters (k=2) (Supplementary Figures 12C, D). Expression levels of PSMB8, BCL2, and PSME2 differed significantly (P < 0.001) between clusters, while BMP5 showed significant differences (P < 0.01) (Figure 8E). A Kaplan-Meier curve was created based on the overall survival (OS) and OS time data (Supplementary Figure S1). ROC curves showed that PSMB8 and PSME2 had high precision (AUC > 0.9) in distinguishing BRCA subtypes, while BCL2 had limited precision (AUC 0.5-0.7) (Supplementary Figures 12F–H).

In the TCGA-BRCA dataset, we developed a prognostic model based on the Hs score. A forest plot (Supplementary Figure 2A) summarizes the results of the Cox regression analysis for Hs. A predictive chart combining Hs, age, and stage was created to assess breast cancer patients’ overall survival (OS) (Supplementary Figure 2B).

We performed a GSEA enrichment analysis on the TCGA-BRCA dataset to assess gene expression impacts on breast cancer across high and low Hs score groups (Figure 8A). This analysis revealed significant associations with various biological pathways, including Cornified Envelope formation (Figure 8B), Keratinization (Figure 8C), Fcgr Activation (Figure 8D), and Fatty Acid Metabolism (Figure 8F), among others (Table 6).

This section examines genetic alterations, focusing on single nucleotide polymorphisms (SNPs) in the TCGA-BRCA dataset. Using R and map tools, mutation analysis identified nine primary types of somatic alterations, with missense mutations being the most common (Supplementary Figure S10). SNPs, particularly C to-T transitions, were predominant (Figures 10A, C). PIK3CA had the highest mutation rate in the high HR score group (38%), while TP53 accounted for 51% of mutations in the low HR score group (Figures 10B, D).

We evaluated drug sensitivity using GDSC database data to predict responses of breast cancer patients with varying HR scores to anti-tumor medications. The Mann–Whitney U test identified significant differences in sensitivity to the top 20 drugs (Figure 12). Notably, distinct sensitivity profiles emerged for these medications between high- and low-HR score categories, indicating substantial disparities.