Abstract Background M2-like tumor-associated macrophages (M2-like TAMs) function crucially in the tumor microenvironment (TME) and cancer development. This study developed a prognostic signature based on M2-like TAM-related genes for breast cancer (BRCA) applying transcriptome and scRNA-seq analysis. Methods TCGA-BRCA, [32]GSE20685, and [33]GSE176078 datasets were downloaded from UCSC xena and GEO databases. AUCell score of immune-related genes (IRGs) was calculated using R package. Genes related to M2-like TAMs were screened by WGCNA. Prognostic genes were further identified by univariate Cox and LASSO regression analyses to form a RiskScore model, which was validated in external dataset. Furthermore, a nomogram was established by integrating RiskScore and clinical characteristics, and correlation analysis between the RiskScore and TME or chemotherapeutic drugs was conducted. Finally, the mRNA expression levels of the key genes identified were verified using quantitative real time polymerase chain reaction (qRT-PCR). Results As macrophages exhibited the highest AUCell score of IRGs in single-cell transcriptomic atlas of BRCA, the cells were further classified into Macrophages C1 and C2 subtypes, with the C1 subtype showing a high expression of M2 macrophage marker genes. ARHGAP26, RILP, KLRB1, CSTA, KLHDC7B, PSMB8, KYNU, RNASE1, LONRF3, and TRPM2 were screened as the prognostic signature genes from a total of 903 M2-like TAM-related genes to establish a robust RiskScore model. Furthermore, a nomogram with a strong predictive performance was constructed combining stage, Age, and RiskScore, and we found that most immune cells showed a negative correlation with RiskScore. Multiple drugs were closely associated with the RiskScore, notably, Ribociclib_1632 had higher a half-maximal inhibitory concentration (IC50) value in high-risk group. Finally, qRT-PCR demonstrated that the mRNA expression levels of the 10 genes were significantly different in control and BRCA cell lines. Conclusion We identified 10 M2-like TAM-related prognostic signature genes for BRCA, providing potential therapeutic targets for the treatment of the cancer. Supplementary Information The online version contains supplementary material available at 10.1007/s12672-025-02161-7. Keywords: Breast cancer, M2 macrophages, Tumor microenvironment, Immune-related genes, Prognosis, Therapeutic targets Introduction Breast cancer (BRCA) is a prevalent female malignancy [[34]1–[35]3] that shows increasing incidence and mortality in recent years [[36]4–[37]6]. BRCA is characterized by specific somatic mutations and variations in gene and protein expressions [[38]7], based on which the cancer is categorized into different subtypes of estrogen receptor (ER), progesterone receptor (PR), and human epidermal growth factor receptor 2 (HER2) [[39]8]. Noticeably, different molecular subtypes of BRCA are distinct in clinical features and prognosis [[40]9]. For example, triple-negative breast cancer (TNBC) is characterized by high metastatic potential and invasiveness, poor prognosis, and a lack of efficacious therapeutic targets [[41]10, [42]11]. In clinical practice, common treatment methods for BRCA are breast-conserving therapy, mastectomy, drug therapy, radiotherapy, chemotherapy, etc. [[43]12]. Despite significant advances in early diagnosis and treatment, which have improved the survival of BRCA patients, the 5-year overall survival (OS) of some patients still remains dismal due to recurrence, metastasis or drug resistance [[44]13, [45]14]. Therefore, the identification of more accurate and reliable prognostic traits and models is essential to stratify risk and provide potential therapeutic targets for BRCA patients. The tumor microenvironment (TME) of BRCA is a complex micro-ecosystem that significantly influences the development and prognosis of BRCA through the interactions between tumor cells and stromal cells or immune cells in the TME [[46]15]. In the TME of breast cancer, TAMs are the most abundant immune cells that appear different phenotypic polarization to perform specific functions [[47]16]. During tumor progression, macrophages can produce opposite effects depending on the polarization of M1 and M2 [[48]17]. Generally, M1 macrophages contribute to a tumor-suppressive TME through pro-inflammatory, cytotoxic, and anti-tumor activities [[49]18], while M2 macrophages are tumor-promoting cells that stimulate angiogenesis by secreting VEGF, provide nutrition for tumor cells to facilitate the progression of BRCA [[50]19, [51]20]. Moreover, M2-like TAMs also play an indispensable role in promoting tumor metastasis and enhancing tumor drug resistance [[52]21]. It has been proven that high infiltration of TAMs is related to worse prognosis and a low survival rate in many cancers [[53]22, [54]23]. Recent research suggested that glycyrrhetinic acid inhibits the M2-like polarization of macrophages by activating JNK1/2 signaling pathway and downregulates VEGF, thereby suppressing the development and metastasis of BRCA [[55]24, [56]25]. Hence, targeting M2-like TAMs may become a new direction for controlling tumor metastasis and drug resistance of BRCA. However, to the best of our knowledge, reliable gene signatures related to M2-like TAMs have not been established for predicting the prognosis of BRCA patients. The single-cell RNA sequencing (scRNA-seq) can be applied to reveal the molecular characteristics of various cell clusters in TME [[57]26, [58]27]. In this study, we downloaded the Cancer Genome Atlas (TCGA)-BRCA, [59]GSE20685, [60]GSE176078 datasets from the public databases. Analysis of the single-cell transcriptomic atlas of BRCA and heterogeneity of macrophages were conducted using [61]GSE176078 dataset. Then, M2-like TAM-related genes in TCGA-BRCA dataset were identified by weighted gene co-expression network analysis (WGCNA), prognostic signature genes were selected to construct a RiskScore model for BRCA, and the robustness of the model was validated in [62]GSE20685 dataset. Furthermore, a nomogram was formed to verify whether the RiskScore was independent of other clinicopathological features. Finally, the correlation between RiskScore and TME and chemotherapeutic drugs was analyzed, aiming to contribute to personalized treatment for BRCA. Material and methods Data collection and preprocessing The TCGA-BRCA dataset in FPKM format were acquired from the University of California, Santa Cruz (UCSC) xena Database ([63]https://xena.ucsc.edu/). After excluding samples without clinical follow-up data or with a survival time shorter than 10 days, a total of 1053 BRCA samples and 113 normal samples were obtained. Then, the Ensembl was converted to Gene symbol, and the maximum value was taken to indicate the expressions of multiple gene symbols. The expression matrix was transformed to Trusted Platform Module (TPM) format and log2-converted. The TCGA-BRCA dataset served as the training set. The scRNA-seq data and clinical data of [64]GSE20685 dataset containing 327 BRCA samples were obtained from the Gene Expression Omnibus (GEO, [65]https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE20685) database and served as the validation set. The probe was mapped to the genes based on the annotation information, and the maximum value was taken to indicate the expression level of the gene when multiple probes matched one gene. Information on the characteristics of the patients was shown in Table [66]S1. The scRNA-seq data of [67]GSE176078 dataset, which contained 26 BRCA samples (10 TNBC, 11 ER+ , and 5 HER2+), were collected from the GEO database ([68]https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE176078) [[69]28]. Only genes expressed in at least 3 cells and a cell expressing at least 200 genes were selected, and cells with mitochondrial gene ratios < 15% and gene numbers between 500–6000 were retained [[70]29]. The SCTransform function was applied to normalize the data [[71]30], and then the RunPCA function was employed for principal component analysis (PCA) [[72]31]. Next, the batch effects between different samples were eliminated using the Harmony R package [[73]32]. Cell clustering analysis and identification of marker genes The top 30 PCs in [74]GSE176078 dataset were subjected to the Uniform Manifold Approximation and Projection for Dimension Reduction (UMAP) for dimensionality reduction [[75]33, [76]34]. The FindNeighbors and FindClusters functions were utilized for cell clustering analysis at a resolution of 0.1 for all BRCA cells and at 0.3 for macrophages. Then, cell types were annotated using the marker genes provided by the CellMarker2.0 database ([77]http://bio-bigdata.hrbmu.edu.cn/CellMarker/). The marker genes with a high expression in each cell subpopulation were identified using the FindAllMarkers function under the conditions of logfc.threshold = 0.5 and min.pct = 0.25. Calculation of AUCell enrichment scores A total of 2483 IRGs were collected from the ImmPort database ([78]https://www.immport.org/shared/). Then, the AUCell enrichment scores of IRGs for each cell type were computed using the AUCell R package [[79]35]. Subsequently, cell clusters containing active "gene sets" in the single-cell data were identified. The AUC score reflected the proportion of top-performing genes identified from a selection of pathway genes in each cell and indicated the activity of particular gene sets within each cell. Pathway enrichment analysis High-expressed marker genes in macrophages were subjected to Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis using the clusterProfiler R package [[80]36]. The pathways with p < 0.05 signified a significant enrichment. Single sample gene set enrichment analysis (ssGSEA) The M2-like TAM scores between normal and tumor samples in TCGA-BRCA dataset were compared by ssGSEA [[81]37], with specific high-expressed genes of M2-like TAMs as the background gene set. According to the optimal critical value of M2-like TAM scores, the patients were split into high- and low- risk groups, and the OS probability of the two risk groups was assessed. WGCNA The key genes relevant to M2-like TAMs of BRCA patients in TCGA dataset were screened using the WGCNA R package [[82]38]. The optimal soft threshold (power) was determined by the pickSoftThreshold function, and the minimum gene number in each module was set as 80 to obtain co-expression modules. Next, the module-trait relationships between the modules and M2-like TAM scores were analyzed using Spearman’s correlation method, and the key module with the highest correlation coefficient was determined. Then, the correlation analysis between gene significance (GS) for M2-like TAM scores and module membership (MM) in the key module was conducted to identify M2-like TAM-related genes under the criteria of |GS|> 0.3 and |MM|> 0.3. Construction and validation of a RiskScore model The prognostic genes for BRCA were identified from the M2-like TAM-related genes by univariate Cox regression analysis (p < 0.05). Further LASSO Cox regression analysis was conducted using the glmnet R package to reduce the total gene number [[83]39]. Then, a prognostic signature based on M2-like TAMs were formulated by stepwise regression analysis, and a RiskScore model was established using threefold cross validation. The RiskScore developed based on the TCGA-BRCA dataset was calculated as follow [[84]40]: [MATH: RiskScore=∑βi∗ExPi :MATH] βi represents the coefficient of a gene in Cox regression model, and ExPi indicates the expression level of the gene. The BRCA patients were separated into high- and low-risk groups based on the optimal critical value of RiskScore, and the OS rate of the patients in different risk groups was compared. Meanwhile, the receiver operating characteristic (ROC) curve was plotted using the timeROC R package [[85]41]. Finally, the robustness of RiskScore model was validated in the [86]GSE20685 validation set. Establishment and evaluation of a nomogram Univariate and multivariate Cox regression analysis was conducted for the TCGA-BRCA dataset to further investigate whether the RiskScore was independent of other clinicopathological features (Age, Gender, stage, pathologic_M, pathologic_N, pathologic_T). Subsequently, a nomogram was established by integrating stage, Age, and RiskScore using the rms R package [[87]42]. The prediction accuracy and reliability of the nomogram were assessed by the calibration curve and decision curve analysis (DCA), respectively. Correlation analysis between the RiskScore and TME The relationship between the RiskScore and TME was analyzed using ESTIMATE and MCP-counter algorithms. The ImmuneScore, StromalScore, and ESTIMATEScore of high- and low-risk groups in TCGA-BRCA dataset were calculated using the ‘ESTIMATE’ R package and ssGSEA [[88]43]. Additionally, the infiltration of 10 immune cells was calculated by the MCP-counter algorithm [[89]44], followed by analyzing the correlation between the RiskScore and infiltration of the 10 types of immune cells applying Spearman method. Correlation analysis between the RiskScore and chemotherapeutic drugs The oncoPredict R package was utilized to calculate the IC50 value for each drug in the TCGA-BRCA [[90]45]. Then, the correlation between the RiskScore, key model genes and 8 chemotherapeutic drugs was evaluated by Pearson correlation analysis. P < 0.05 and |cor|> 0.3 denoted a significant correlation. Cell acquisition and culture Breast cancer cell lines (MCF-7, MDA-MB-231, and SK-BR-3) and normal human breast epithelial cells MCF-10A were sourced from the American Type Culture Collection (ATCC). SK-BR-3 cells were grown in DMEM that contained 10% (v/v) heat-inactivated fetal bovine serum (FBS), 2 mM L-glutamine, and a combination of 100 U/mL penicillin along with 100 µg/mL streptomycin. And other cells were cultivated in DMEM (Gibco Laboratories, USA) supplemented with 10% FBS (Gibco Laboratories, USA) and 1% penicillin–streptomycin, and incubated at 37 °C in a 5% CO[2] atmosphere. Quantitative real-time polymerase chain reaction (qRT-PCR) Total RNA was isolated using TRIzol (Invitrogen, USA) and used to synthesize cDNA with a SuperScript™ III reverse transcriptase kit (Invitrogen, USA). Quantitative RT-PCR (qRT-PCR) was conducted on an ABI Prism 7500 utilizing SYBR Green qPCR Master Mix (Takara, Shanghai, China). The qPCR conditions included an initial denaturation at 94 °C for 30 s, followed by 40 cycles of 94 °C for 5 s and 60 °C for 30 s. The relative expression levels of ARHGAP26, RILP, KLRB1, CSTA, KLHDC7B, PSMB8, KYNU, RNASE1, LONRF3, and TRPM2 were measured by the 2^−ΔΔCt method, with β-actin serving as an internal control. All the experiments were performed in triplicate. Primer sequences are shown in Table [91]1. Table 1. Sequence of primers Gene Primers (5ʹ–3 ʹ) Forward (F) Reverse (R) ARHGAP26 GACTGCTGCCTCGATAGTCC CTTCCGCTTCGCTGAAGACAA RILP TTGGCCTATGGTATCGGGGTA CCCCAGACAAAGGTGTTCGT KLRB1 TGGCATCAATTTGCCCTGAAA TCCAAGGGTTGACAGTGTGAG CSTA AAACCCGCCACTCCAGAAATC CACCTGCTCGTACCTTAATGTAG KLHDC7B GCACCATGCACAACTACCTGT ATTCGCCACCGATGGCATAG PSMB8 GCAGGCTGTACTATCTGCGAA AGAGCCGAGTCCCATGTTCAT KYNU GGCTCTCCACCTAGATGAGGA GCTGCTATTTTGGCCCACTTAT RNASE1 ACTGTAACCAAATGATGAGGCG GTACCTGGAGCCGTTTGTCA LONRF3 GCTCCCACATTGTTCTAGTCAG AGAGGTAGCATCCCACTTCTC TRPM2 TTCGTGGATTCCTGAAAACATCA CCAGCATCAGACAGTTTGGAAC β-actin GGACATCCGCAAAGACCTGTA GCTCAGGAGGAGCAATGATCT [92]Open in a new tab Statistical analysis All statistical analyses were performed in R language (version 4.2.0) and GraphPad Prism (version 8.0.2). Comparisons were performed using a two-way ANOVA combined with Turkey or Sidak test. The OS rate was visualized by the Kaplan–Meier method. We used log-rank test to compare survival differences between the two groups. In addition, differences in continuous variables between the two groups were compared by Wilcoxon rank-sum test. The Spearman analysis was performed to assess the correlation between risk score and the level of immune cell infiltration. A P < 0.05 was considered to be statistically significant. Results Macrophages exhibited the highest AUCell score of IRGs in single-cell transcriptomic atlas of BRCA After data normalization, downscaling and cell clustering analysis, a total of 13 cell clusters were identified based on the scRNA-seq data of BRCA samples in the [93]GSE176078 dataset (Fig. [94]1A). According to the expressions of classical maker genes in each cell cluster, the following 6 cell types were classified (Fig. [95]1B): T cells (marked with CD2, CD3D, CD3E, and CD69), epithelial cells (marked with EPCAM, KRT18, and KRT8), macrophages (marked with CD68, CD163, AIF1, FCER1G, LYZ, and TYROBP), endothelial cells (marked with PECAM1, PLVAP, and VWF), fibroblasts (marked with COL1A1, COL1A2, DCN, and LUM), and B cells (marked with CD79A, DERL3, JCHAIN, MZB1, and MS4A1) (Fig. [96]1C). Furthermore, quantification of each cell type in 26 BRCA samples showed that T cells accounted for a higher proportion in the samples (Fig. [97]1D). In addition, macrophages exhibited the highest AUCell score of IRGs in comparison to the rest 5 cell types (Fig. [98]1E). Fig. 1. [99]Fig. 1 [100]Open in a new tab The single-cell transcriptomic atlas of breast cancer (BRCA) in [101]GSE176078 dataset. A UMAP plot showing the 13 BRCA cell clusters from the [102]GSE176078 dataset. B UMAP plot showing the 6 cell types of the scRNA-seq dataset from BRCA after annotation. C Bubble plot of classical marker genes in each cell type. D The proportion of cell types in 26 BRCA samples. E AUCell enrichment scores of immune-related genes (IRGs) in different cell types The marker genes of M2 macrophages showed a high expression in Macrophage C1 The macrophages were re-clustered by UMAP clustering into two subpopulations of Macrophages C1 and C2 (Fig. [103]2A), which were marked as Macrophage C1 (CD163, FGL2, and MSR1) and Macrophage C2 (CIB1 and KYNU), respectively (Fig. [104]2B). Macrophage C1 showed higher AUCell score of IRGs than that in Macrophage C2 (Fig. [105]2C). Moreover, Macrophage C1 had higher expression of the marker genes of M2 macrophages such as CD163, FGL2, and MSR1 (Fig. [106]2D). Thus, Macrophage C1 was selected to represent M2 macrophages for further study. KEGG enrichment analysis indicated that specifically high-expressed genes of M2 macrophages were notably enriched in Phagosome, IL-17 signaling pathway, Chemokine signaling pathway, NK-kappa B signaling pathway, Antigen processing and presentation, TNF signaling pathway, Apoptosis, etc. (Fig. [107]2E). Fig. 2. [108]Fig. 2 [109]Open in a new tab The heterogeneity of macrophages. A UMAP plot showing two macrophage subtypes (C1 and C2) in BRCA. B The marker genes with a high expression in Macrophage C1 and Macrophage C2. C AUCell scores of IRGs in Macrophage C1 and Macrophage C2. D Violin plot of the expression levels of marker genes in M2 macrophages. E KEGG pathway enrichment analysis of the marker genes in M2 macrophages A total of 903 genes closely related to M2-like TAMs were identified by WGCNA The ssGSEA showed that the M2-like TAM score of the tumor group was noticeably lower than that of normal group in TCGA-BRCA dataset (Fig. [110]3A), indicating that M2-like TAMs were dysregulated in the tumor samples. According to the Kaplan–Meier survival curve, high M2-like TAM score group showed a higher 10-year OS rate of BRCA patients than that in low M2-like TAMs score group, however, the trend was reversed after 10 years (Fig. [111]3B). Further, WGCNA was applied to identify key module correlated with M2-like TAMs of BRCA patients in TCGA-BRCA dataset. The optimal soft threshold (power) was determined to be 6 by pickSoftThreshold function (Fig. [112]S1) and 9 modules were acquired (Fig. [113]S2). The number of genes contained in each module was shown by a lollipop plot (Fig. [114]3C). The heatmap of module-trait relationships presented that the black module was significantly positively correlated with M2-like TAM score (cor:0.88, p:0.00) (Fig. [115]3D). Then, 903 genes closely linked to M2-like TAMs were screened based on the |GS|> 0.3 and |MM|> 0.3 (Fig. [116]3E). Fig. 3. [117]Fig. 3 [118]Open in a new tab Identification of M2-like tumor-associated macrophage (TAM)-related genes in TCGA-BRCA dataset by WGCNA. A M2-like TAM scores in tumor and normal samples by ssGSEA. **indicates p < 0.01. B OS survival curves between high and low M2-like TAM score groups. C The number of genes included in each module identified by WGCNA. D The heatmap of module-trait relationships. E Scatter plot of module membership (MM) in black module versus gene significance (GS) for M2-like TAM score A robust RiskScore model was developed based on 10 M2-like TAM-related prognostic genes Firstly, the prognostic genes were screened from a total of 903 M2-like TAM-related genes by univariate and LASSO Cox regression analysis, and 17 genes with lambda = 0.0172 were selected for subsequent study (Fig. [119]4A, B). Next, 10 prognostic signature genes associated with M2-like TAMs were obtained, including 6 protective genes (ARHGAP26, RILP, KLRB1, CSTA, KLHDC7B, and PSMB8) and 4 risk genes (KYNU, RNASE1, LONRF3, and TRPM2) (Fig. [120]4C). The RiskScore model was formulated: RiskScore = 0.42*TRPM2 + 0.295*RNASE1-0.266*RILP-0.146*PSMB8 + 0.353*LO NRF3 + 0.197*KYNU-0.252*KLRB1-0.153*KLHDC7B-0.166*CSTA-0.287*ARHGAP26. Notably, high-risk group had a lower OS rate and shorter survival time than the low-risk group in TCGA-BRCA dataset (Fig. [121]4D). The area under ROC curve (AUC) analysis showed that the Riskscore model was reliable, with 1-year AUC of 0.73, 3-year AUC of 0.76, 5-year AUC of 0.72, and 7-year AUC of 0.73 (Fig. [122]4D). In addition, validation of the Riskscore model using the [123]GSE20685 dataset demonstrated that the 1-, 3-, 5-, and 7- year AUC value reached 0.77, 0.67, 0.66, and 0.61, respectively (Fig. [124]4E). Fig. 4. [125]Fig. 4 [126]Open in a new tab Establishment and validation of a RiskScore model. A The coefficients of each independent variable changed with lambda. B The confidence interval under lambda. C The prognostic signature of RiskScore model. D The ROC and OS survival curves of RiskScore model in TCGA-BRCA dataset. E The ROC and OS survival curves of RiskScore model in validation set [127]GSE20685 Development of a nomogram with strong performance for predicting the prognosis of BRCA The results of univariate and multivariate Cox regression analyses showed that stage, Age, and RiskScore were the independent prognostic factors for BRCA (Fig. [128]5A, B). Then, a nomogram was developed by integrating stage, Age, and the RiskScore (Fig. [129]5C). It was observed that the predictive calibration points for 1-, 3-, 5-, and 7 year calibration curve were close to the standard curves (Fig. [130]5D), indicating that the nomogram had a strong predictive performance. The net benefit of the nomogram was also higher than the extreme curve by DCA curve (Fig. [131]5E), which supported that the nomogram was reliable. Fig. 5. [132]Fig. 5 [133]Open in a new tab Construction and validation of a nomogram. A Univariate Cox regression analysis on the RiskScore and clinicopathological features in TCGA-BRCA dataset. B Multivariate Cox regression analysis on the RiskScore and clinicopathological features in TCGA-BRCA dataset. C Construction of a nomogram by combining stage, Age, and RiskScore. D The calibration curve of the nomogram. E The DCA curve of the anomogram BRCA patients with a low RiskScore were more sensitive to chemotherapeutic drugs The StromalScore, ImmuneScore, and ESTIMATEScore of high-risk group were considerably lower than the low-risk group in TCGA-BRCA dataset (Fig. [134]6A), suggesting that a lower immune score was related to tumor progression and a worse prognosis of BRCA patients. Moreover, correlation analysis between 10 types of immune cells and the RiskScore was assessed by MCP-counter algorithm, and it was found that most immune cells exhibited a negative correlation with the RiskScore (Fig. [135]6B). In addition, there was a close correlation between the RiskScore and the IC50 of 8 chemotherapeutic drugs (FDR < 0.05 and |cor|> 0.3) (Fig. [136]6C), such as Ribociclib_1632 (R = 0.34) (Fig. [137]6D). The IC50 of Ribociclib_1632 in high-risk group was higher than that in low-risk group (Fig. [138]6E), indicating that the BRCA patients with a low RiskScore were more sensitive to Ribociclib_1632. Fig. 6. [139]Fig. 6 [140]Open in a new tab Correlation analysis between the RiskScore and TME or chemotherapeutic drugs in TCGA-BRCA dataset. A StromalScore, ImmuneScore, and ESTIMATEScore of high- and low-risk groups assessed by ESTIMATE algorithm. B Correlation analysis between the RiskScore and infiltration of 10 types of immune cells calculated by MCP-counter algorithm. C Correlation analysis between the RiskScore, key model genes and 8 chemotherapeutic drugs. D Correlation analysis between the RiskScore and Ribociclib_1632. E IC50 values of Ribociclib_1632 between high- and low-risk groups. ****indicates p < 0.0001; ***indicates p < 0.001; **indicates p < 0.01; *indicates p < 0.05 Cell-based experiments to validate the expression levels of key genes To further validate the key genes, we measured the mRNA levels of these 10 genes using breast cancer cell lines MCF-7, MDA-MB-231, and SK-BR-3 and normal breast epithelial cell line MCF-10A. It was observed that the four risk genes were significantly overexpressed in breast cancer cells (Figs. [141]7A, [142]S3A), while the mRNA expression levels of the six protective genes were significantly downregulated in breast cancer cells when compared to normal breast epithelial cells MCF-10A (Figs. [143]7B, [144]S3B). Fig. 7. [145]Fig. 7 [146]Open in a new tab Using qRT-PCR to measure the mRNA expression levels of the 10 key genes. A Differential mRNA expression levels of the four risk genes (KYNU, RNASE1, LONRF3, and TRPM2) in breast cancer cell lines MCF-7 and MDA-MB-231 and in normal breast epithelial cells MCF-10A. B Differential mRNA expression levels of the six protective genes (ARHGAP26, RILP, KLRB1, CSTA, KLHDC7B, and PSMB8) in breast cancer cell lines MCF-7 and MDA-MB-231 and in normal breast epithelial cells MCF-10A Discussion M2-like TAMs play a critical role in the development and metastasis of BRCA via regulating immune response and TME, which is found to be correlated with a worse prognosis of BRCA patients [[147]46, [148]47]. This study constructed the single-cell transcriptomic atlas of BRCA, and found that macrophages exhibited the highest AUCell score of IRGs. After further division of macrophages into C1 and C2 subtypes, it was found that Macrophage C1 showed the characteristics of M2 macrophages. Based on WGCNA, a total of 903 M2-like TAM-related genes were obtained. Then, 10 prognostic signature genes were identified, including 6 protective genes (ARHGAP26, RILP, KLRB1, CSTA, KLHDC7B, and PSMB8) and 4 risk genes (KYNU, RNASE1, LONRF3, and TRPM2). We further developed a RiskScore and a nomogram, which all exhibited a strong performance for predicting the prognosis of BRCA. Moreover, the RiskScore model was closely linked to immune cell infiltration and drug sensitivity in BRCA. Numerous studies showed that TAM-related genes can serve as prognostic genes for BRCA [[149]48]. Previous research identified 45 prognostic genes from TAM-regulated BRCA cell MCF7 transcriptome. Further verification showed that these genes were almost relevant to all the clinical features of BRCA, showing a highly accurate prediction of the OS and prognostic outcomes for BRCA patients [[150]49]. Additionally, another study discovered 37 genes associated with human breast TAMs, which were significantly enriched in the aggressive BRCA subtypes and were related to shorter disease-specific survival [[151]50]. M2-like TAMs were dominant in the breast TME, markedly influencing the progression and metastasis of BRCA [[152]51]. Chang identified 722 M2-like TAM genes using WGCNA, classified two modular subtypes of M2-like TAMs in BRCA, and constructed a prognostic signature for evaluating the prognosis and immunotherapy efficacy of BRCA [[153]52]. In our present study, since macrophages exhibited the highest AUCell score of IRGs, the cells were further separated into C1 and C2 subtypes, with Macrophage C1 showing the characteristics of M2 macrophages. Then, WGCNA identified 903 M2-like TAM-related genes, based on which 10 prognostic genes were screened to establish a RiskScore model. According to the RiskScore, BRCA patients were divided into high- and low-risk groups, with high-risk patients having lower OS rate and worse prognosis. Moreover, the RiskScore was strongly correlated with immune cell infiltration and drug sensitivity. The 10 M2-like TAM-related prognostic genes for BRCA contained 6 protective genes (ARHGAP26, RILP, KLRB1, CSTA, KLHDC7B, and PSMB8) and 4 risk genes (KYNU, RNASE1, LONRF3, and TRPM2). ARHGAP26 belongs to a Rho GTPase-activating protein family that can transform GTPases into an inactive form, thereby suppressing the signaling transduction [[154]53]. ARHGAP26 is involved in the progression of numerous cancers, including ovarian cancer, and is therefore considered as a putative tumor suppressor [[155]54]. RILP, a Rab-interacting lysosomal protein, fulfills a crucial functions in regulating the transport of lysosomal [[156]55]. Study showed that RILP inhibits the proliferation and invasion of BRCA [[157]56], prostate cancer, and lung cancer cell lines [[158]57]. KLRB1 is a vital member of killer cell lectin-like receptor family [[159]58]. It has been reported that KLRB1 is an independent factor for BRCA prognosis and low-expressed KLRB1 in BRCA patients is predictive of a poor survival outcome [[160]59, [161]60]. CSTA is a cysteine protease inhibitor and an important estrogen-regulated gene in BRCA cells [[162]61]. The deletion of CSTA in breast tumors can facilitate extracellular matrix remodeling and induce the invasion and metastasis of BRCA [[163]62]. KLHDC7B is a Kelch domain-containing protein that plays crucial roles in protein degradation, gene expression, and extracellular communication [[164]63]. KLHDC7B has been recognized as an epigenetic marker, which is highly methylated at the promoter but shows an upregulated expression in BRCA [[165]64]. PSMB8 encodes a critical subunit of the special immunoproteasome complex [[166]65]. Previous study showed that a higher mRNA expression of PSMB8 is indicative of a better prognosis of basal-like and TNBC patients [[167]66]. In addition, KYNU serves as an enzyme involved in tryptophan catabolism [[168]67]. The expression of KYNU is associated with the malignancy grade of BRCA and usually indicates poor outcomes [[169]68]. RNASE1 belongs to the pancreatic RNase A superfamily and is involved in hemostasis, anti-inflammatory, and innate immunity [[170]69]. The overexpression of RNASE1 in BRCA predicts a worse poor survival [[171]70]. LONRF3 encodes a ring-finger domain-containing protein. Study revealed that LONRF3 is linked to unfavorable prognosis in pancreatic cancer [[172]71], but the role of LONRF3 in BRCA has not been reported. TRPM2 is a key channel of plasmalemma and lysosome in immune and tumor cells [[173]72]. TRPM2 is high-expressed in many cancers, such as ovarian cancer, neuroblastoma, and BRCA, suggesting that TRPM2 may facilitate the development of cancers [[174]73]. Therefore, the above results demonstrated that M2-like TAM-related prognostic genes can be considered as potential targets for BRCA treatment. However, our current study also had some limitations. Firstly, the data used in the present work had a small samples size and were only collected from the public database, therefore further validation is needed by adding clinical samples and detailed experimental data. Secondly, the prognostic genes involved in the progression of BRCA, such as ARHGAP26 and LONRF3, are rarely reported in BRCA. Hence, the biological functions and molecular mechanisms of the 10 prognostic genes in BRCA should be further investigated using basic experiments and large-scale clinical trials. Conclusion In summary, by integrating transcriptome and single-cell sequencing analysis, this study identified a total of 10 M2-like TAM-related prognostic genes for BRCA, including ARHGAP26, RILP, KLRB1, CSTA, KLHDC7B, PSMB8, KYNU, RNASE1, LONRF3, and TRPM2. A robust RiskScore model was established as an effective tool for the prediction of BRCA prognosis. Further analysis showed that the RiskScore was closely linked to immune cell infiltration and drug sensitivity in BRCA. To conclude, the present work discovered potential targets and predicted potential drugs for the clinical treatment of BRCA. Supplementary Information [175]Additional file1 (DOCX 379 KB)^ (381.9KB, docx) [176]Additional file2 (ZIP 21 KB)^ (20.7KB, zip) Abbreviations AUC Area under ROC curve BRCA Breast cancer DCA Decision curve analysis ER Estrogen-receptor ESTIMATE Estimation of STromal and Immune cells in MAlignant tumour tissues using expression data FPKM Fragments Per Kilobase of exon model per million mapped fragments GEO Gene expression omnibus GS Gene significance HER2 Human epidermal growth factor receptor 2 IC50 Half maximal inhibitory concentration IRGs Immune-related genes KEGG Kyoto encyclopedia of genes and genomes LASSO Least absolute shrinkage and selection operator M2-like TAMs M2-like tumor-associated macrophages MCP-counter Microenvironment cell populations-counter MM Module membership OS Overall survival PCA Principal component analysis PR Progesterone-receptor ROC Receiver operating characteristic scRNA-seq Single-cell RNA-sequencing ssGSEA Single sample gene set enrichment analysis TAMs Tumor-associated macrophages TCGA The cancer genome atlas TME Tumor microenvironment TNBC Triple-negative breast cancer TPM Trusted platform module UCSC University of Cingifornia Sisha Cruz UMAP Uniform manifold approximation and projection for dimension reduction VEGF Vascular endothelial growth factor WGCNA Weighted gene co-expression network analysis Author contributions Jinhong Zhou and Qing Zhang contributed to the study conception and design. Material preparation, data collection and analysis were performed by Yanghaochen Xu and Peiyan Lin. The first draft of the manuscript was written by Peiyan Lin and then Ye Zhu modified it. All authors discussed the results and revised the manuscript. All authors reviewed and approved the submission of final manuscripts. Funding The authors declare that no funds, grants, or other support were received during the preparation of this manuscript. Data availability The datasets used and analyzed in this study are publicly available through the Gene Expression Omnibus (GEO) database ([177]http://www.ncbi.nlm.nih.gov/geo) and UCSC xena database ([178]https://xena.ucsc.edu/). The TCGA-BRCA dataset in FPKM format were acquired from UCSC xena Database. The scRNA-seq data were acquired from GEO database, with the accession numbers: [179]GSE20685, [180]GSE176078. Declarations Ethics approval and consent to participate This manuscript is not a clinical trial, hence the ethics approval and Consent to participation are not applicable. Competing interests The authors have no relevant financial or non-financial interests to disclose. Consent for publication Not applicable. Footnotes Publisher's Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Yanghaochen Xu, Peiyan Lin have contributed equally to this work . References