Decoding the Tumor Microenvironment of Myoepithelial Cells in Triple-Negative Breast Cancer Through Single-Cell and Transcriptomic Sequencing and Establishing a Prognostic Model Based on Key Myoepithelial Cell Genes

Abstract

Background: Triple-negative breast cancer (TNBC) is an aggressive subtype with high malignancy, rapid progression, and a poor 5-year survival rate of ~77%. Due to the lack of targeted therapies, treatment options are limited, highlighting the urgent need for novel therapeutic strategies. Myoepithelial cells (MECs) in the tumor microenvironment may significantly influence TNBC development and progression. Methods: This study used single-cell RNA sequencing to compare the MEC gene expression in the normal versus TNBC tissues. TNBC-associated MECs showed increased expression of proliferation- and immune-related genes (e.g., FDCSP, KRT14, and KRT17) and decreased expression of inflammatory and extracellular matrix-related genes (e.g., CXCL8, SRGN, and DCN). Copy number variation and pseudotime analyses revealed genomic alterations and phenotypic dynamics in MECs. A CoxBoost-based prognostic model was developed and validated across 20 survival cohorts, integrating immune profiling, pathway enrichment, and drug sensitivity analyses. Mendelian randomization identified TPD52 as a TNBC risk-associated gene. siRNA knockdown of TPD52 was performed in TNBC cell lines to evaluate its effects on proliferation and migration. Results: TNBC MECs displayed significant changes in the gene expression and genomic integrity, impacting immune responses and tumor invasion. The prognostic model effectively predicted 1-, 3-, and 5-year survival outcomes, stratifying high-risk patients with enriched cell cycle and DNA replication pathways, reduced immune checkpoint expression, and chemotherapy resistance. TPD52 was identified as a tumor-promoting gene, and its knockdown suppressed TNBC cell proliferation and migration. Conclusion: This study highlights MECs' role in TNBC progression, provides a CoxBoost prognostic model for personalized treatment, and identifies TPD52 as a potential therapeutic target for TNBC intervention.

Keywords: TNBC; myoepithelial cells; prognosis.

Figure 1

Cellular landscape and molecular features of the seven identified cell subtypes in triple-negative breast cancer (TNBC). (a) Box plots showing the gene expression metrics (nCount, gene count, and mitochondrial gene ratio) for the seven identified cell subtypes in TNBC tissue. (b) The top three Gene Ontology (GO) enrichment pathways identified in each of the seven distinct cell subtypes. (c) UMAP visualization of the seven cell subtypes in TNBC and normal control (NC). (d) bar plots showing the differences in the proportions of the seven cell subtypes between TNBC and NC samples. (e) Differential gene expression analysis for the seven cell subtypes between normal and TNBC tissues, highlighting the top five significantly differentially expressed genes for each subtype.

Figure 2

Genomic alterations in myoepithelial cells from TNBC. (a) Heatmap depicting the landscape of copy number variations (CNVs) in myoepithelial cells, using T/NK cells as a reference. (b) K-means clustering based on CNVs, highlighting the similarities between TNBC myoepithelial cells in Cluster 4 and T/NK cells. (c) Monocle-based trajectory predictions of CRC tumor cell subpopulations. (d) Classification of pseudotime-dependent genes into four main categories and the associated pathway enrichment. (e) Heatmap showing the marker gene expression across different branches, annotated into four major clusters.

Figure 3

Identification of coexpression modules and prognostic hub genes in TNBC. (a) Weighted gene coexpression network analysis (WGCNA) of TNBC cells. (b) Visualization of the coexpression network structure across different modules. (c) Ranking of the top 10 eigengenes for each module based on module eigengene connectivity (kME). (d) GO enrichment analysis of highly enriched genes within the 16 TNBC coexpression modules.

Figure 4

Performance comparison of prognostic regression models. (a) Heatmap demonstrating the mean AUCs (1-year, 3-year, and 5-year) for Lasso regression, elastic net, ridge regression, stepwise Cox regression, and CoxBoost models across 15 cohorts. (b) Heatmap illustrating the regression coefficients of the input genes across various prognostic models. (c) Univariate Cox regression analysis of the risk scores combined with meta-analysis. (d) Kaplan–Meier survival analysis across 20 survival cohorts from 12 datasets.

Figure 5

Immune landscape analysis in high- and low-risk TNBC groups. (a) Heatmap comparing immune cell infiltration between the high- and low-risk groups using TIMER, CIBERSORT, CIBERSORT-ABS, QUANTISEQ, MCPCOUNTER, xCell, and EPIC algorithms. (b) Risk and subtype group comparisons in the TCGA dataset. (c) Differential expression analysis for immune checkpoint genes between the high- and low-risk groups (⁣^∗p < 0.05, ⁣^∗∗p < 0.01, and ⁣^∗∗∗p < 0.001). (d) Comparison of immune, stromal, and tumor purity scores between the two risk groups. (e) GSEA enrichment analysis results for the low-risk group. (f) GSEA enrichment analysis results for the high-risk group.

Figure 6

Drug sensitivity analysis. Boxplots showcasing the comparative drug sensitivity (IC50, half-maximal inhibitory concentration) across the high- and low-risk groups, revealing potential therapeutic responses.

Figure 7

Integration of Mendelian randomization analysis demonstrating the oncogenic role of TPD52 in breast cancer. (a) The upper panel shows results from Finngen_R10_C3_BREAST_ERPLUS_EXALLC dataset. (b) The lower panel shows results from Finngen_R10_C3_BREAST_EXALLC dataset. Each panel consists of three plots: Left plots display MR effect sizes for ‘TPD52' on ‘BRCA' with different SNPs (rs3117098, rs15126492, rs4739729, rs68042249, and rs79313939) and aggregated estimates (All-MR egger and All- inverse variance weighted). Middle plots show MR leave-one-out sensitivity analysis for ‘TPD52' on ‘BRCA'. Right plots illustrate SNP effects on TPD52 versus BRCA with different statistical methods (inverse variance weighted, MR egger, simple mode, weighted median, and weighted mode) represented by different colored lines. Error bars represent confidence intervals. The x-axis in all right plots represents SNP effect on TPD52, while the y-axis shows SNP effect on BRCA.

Figure 8

(a) IHC analysis shows the TPD52 protein expression levels in TNBC (Tumor) or the paired nontumor tissues (Normal), (scale bar : 50 μm). (b, c) Knocking out the TPD52 in breast carcinoma cell lines by siRNA. Assessing the effects of TPD52 by (d, e) CCK8 and (f, g) colony formation. (h, i) Assessment of migration by Transwell assay.