Association of MPO Expression with the Immune Microenvironment in Breast Cancer: Insights from Bioinformatics and Single-Cell Analyses

Summary

This article presents a reproducible bioinformatics and single-cell workflow for exploring associations between myeloperoxidase (MPO) expression and immune/myeloid features in breast cancer. Because the analyses are based on public datasets and in silico methods, the findings are interpreted as exploratory and hypothesis-generating.

Abstract

Breast cancer remains a major cause of cancer-related mortality, and exploratory computational workflows can help prioritize immune-associated markers for further investigation. Here, we used the cancer genome atlas breast invasive carcinoma (TCGA-BRCA) bulk transcriptomic data and the public single-cell dataset GSE161529 to examine associations between myeloperoxidase (MPO) expression, clinical outcomes, immune infiltration, methylation, upstream-regulator annotations, single-cell expression patterns, virtual knockdown sensitivity outputs, drug–gene interaction retrieval, and absorption, distribution, metabolism, excretion, and toxicity (ADMET) annotation. MPO expression was lower in breast cancer tissues than in adjacent non-tumor tissues. Higher MPO expression was associated with a longer progression-free interval, whereas its associations with overall survival and disease-specific survival were not statistically significant. Receiver operating characteristic (ROC) analysis suggested tumor–normal separation within the analyzed public dataset, but this should not be interpreted as clinical diagnostic validation. Immune deconvolution and enrichment analyses indicated that MPO expression mainly tracked with immune- and myeloid-related transcriptional features, rather than establishing tumor-intrinsic regulation of the immune microenvironment. At single-cell resolution, the MPO signal was sparse, with only 85 MPO-positive cells detected before k-nearest neighbor (KNN)-based neighborhood expansion. Detectable MPO signal and MPO-associated scores were interpreted cautiously because they may be influenced by sparse expression, cell-type annotation uncertainty, dropout, doublets, or ambient RNA. In silico virtual knockdown suggested candidate immune- and inflammatory-related transcriptional changes, but these results were considered exploratory and require validation. Drug-gene interaction database (DGIdb)-based drug-gene retrieval and ADMET annotation were used only as preliminary chemical annotations and were not interpreted as therapeutic evidence. Overall, this study provides a reproducible in silico workflow for generating hypotheses about MPO-associated immune/myeloid features in breast cancer, which require external cohort validation and experimental confirmation.

Introduction

Breast cancer is a highly heterogeneous immune-related malignancy¹. Its disease progression, risk of recurrence and metastasis, and treatment response are closely associated with the composition and functional status of the tumor immune microenvironment (TIME)². Despite ongoing optimization of comprehensive treatment strategies, some patients still experience progression or recurrence, underscoring the urgent need to identify molecular biomarkers that characterize TIME status and support risk stratification, while elucidating their underlying mechanisms.

Myeloperoxidase (MPO) is a heme-containing peroxidase predominantly expressed in neutrophils and, to a lesser extent, monocytes and macrophages. Through the generation of hypochlorous acid and other reactive oxidants, MPO contributes to antimicrobial defense but may also promote oxidative tissue injury and chronic inflammation. In cancer, the biological significance of MPO appears to be context-dependent³. On the one hand, MPO-mediated oxidative stress has been implicated in carcinogenesis and tumor progression through DNA damage, lipid and protein oxidation, inflammatory signaling, and remodeling of the tumor microenvironment^4,5,6. On the other hand, MPO-positive innate immune or myeloid cell infiltration has been associated with a favorable prognosis or antitumor immune activity in certain tumor contexts^7,8,9. These apparently conflicting findings suggest that the clinical and biological significance of MPO may depend on tumor type, disease stage, cellular source of MPO, and the immune composition of the tumor microenvironment. However, the expression pattern and prognostic relevance of MPO in breast cancer, particularly at the single-cell level, remain incompletely characterized.

The tumor immune microenvironment (TIME) contains heterogeneous myeloid, lymphoid, stromal, and epithelial compartments¹⁰. MPO is classically associated with neutrophils and other myeloid-lineage cells, and MPO-related signals in bulk tumor profiles may therefore reflect immune-cell composition rather than tumor-cell-intrinsic activity¹⁰. In breast cancer, the distribution of MPO signal across bulk and single-cell datasets, its association with immune-infiltration estimates, and the reproducibility limits of downstream computational analyses remain insufficiently characterized. This study, therefore, treats MPO as an immune-associated marker for exploratory workflow development, not as a proven causal regulator of the TIME or a validated therapeutic target. Compared with single-cohort differential expression analyses or single-platform immune-infiltration estimates, an integrated workflow combining bulk transcriptomics, immune deconvolution, methylation annotation, single-cell mapping, and computational perturbation can provide a broader exploratory view of gene-associated immune context. This approach is useful for prioritizing candidate markers and generating testable hypotheses, especially when experimental datasets are not yet available. However, such computational integration cannot by itself determine cellular source, causality, pharmacological activity, or clinical utility. With the advancement of large-scale public cancer cohorts and single-cell transcriptomic technologies, bioinformatics approaches can be used to explore associations among gene expression, clinical outcomes, immune-cell composition, and transcriptional states at both population and single-cell levels¹¹. Computational perturbation methods based on single-cell gene regulatory networks may further provide hypothesis-generating information regarding gene-associated transcriptional sensitivity^12,13. Therefore, this study aimed to characterize the expression pattern, survival association, immune/myeloid context, methylation profile, single-cell distribution, and exploratory computational perturbation profile of MPO in breast cancer. The overall workflow is shown in Figure 1.

Protocol

Acquisition from the TCGA database

RNA sequencing data and clinical information for the TCGA breast invasive carcinoma (TCGA-BRCA) cohort were obtained from the genomic data commons portal¹⁴. STAR workflow RNA-seq data in transcripts per million (TPM) format were extracted together with matched clinical annotations. RNA-seq samples lacking corresponding clinical information were excluded. For expression-based analyses, TPM values were transformed as log2(TPM + 1). MPO expression was extracted using the gene symbol MPO and Ensembl gene ID ENSG00000005381.8. For analyses requiring MPO-high and MPO-low grouping, only TCGA-BRCA tumor samples were included, and adjacent normal samples were excluded from group assignment. Tumor samples were divided according to the median value of log2(TPM + 1)-transformed MPO expression among TCGA-BRCA tumor samples. Samples with MPO expression greater than or equal to the median were assigned to the MPO-high group, whereas samples below the median were assigned to the MPO-low group. This median-based grouping strategy was used for survival analysis, differential expression analysis, enrichment analysis, methylation grouping, and immune-cell enrichment comparisons unless otherwise specified. Clinical pathological characteristics, including sex, age, ethnicity, pathological T stage, histological grade, PAM50 subtype, pathological stage, tumor status, and survival endpoints, including overall survival (OS), progression-free interval (PFI), and disease-specific survival (DSS), were analyzed using R version 4.2.1.

Public Immunohistochemistry Image Retrieval

Representative MPO immunohistochemistry (IHC) images of adjacent normal breast tissue and breast cancer tissue were used as qualitative protein-level references. These images were not included in quantitative morphometric or statistical analyses. The boxed areas indicate regions shown at higher magnification. Scale bars indicate 100 µm in the 20× images and 50 µm in the 40× images.

Expression correlation analysis

The TCGA-BRCA dataset was used to examine genes co-varying with MPO expression in breast cancer. Genome-wide Pearson correlation coefficients were calculated between MPO and protein-coding genes, and the top 30 positively and top 30 negatively correlated genes were selected for visualization. For correlation analyses involving multiple tested genes, nominal p-values were adjusted using the Benjamini-Hochberg false discovery rate method. The MPO-associated protein-protein interaction (PPI) network was constructed using the search tool for the retrieval of interacting genes/proteins (STRING) database, with protein pairs displaying interaction scores greater than 0.40 retained for visualization¹⁵.

Functional enrichment analysis

Differentially expressed genes (DEGs) were identified by comparing the MPO-high and MPO-low TCGA-BRCA tumor groups using thresholds of |log2FC| > 1 and a Benjamini-Hochberg-adjusted p-value < 0.05. Functional enrichment analysis of DEGs was performed using the R package clusterProfiler version 4.4.4, including gene ontology (GO) biological process, cellular component, molecular function, and Kyoto encyclopedia of genes and genomes (KEGG) pathway analyses^{16,17,18,19,20}. Enriched GO and KEGG terms were considered significant when the adjusted p-value was < 0.05.

Gene set enrichment analysis (GSEA) was performed using a pre-ranked gene list based on differential expression statistics between MPO-high and MPO-low groups. The MSigDB C2 Canonical Pathways collection c2.cp.all.v2022.1.Hs.symbols.gmt, corresponding to MSigDB v2022.1.Hs and containing 3,050 gene sets, was used^21,22. Enriched terms were considered significant according to Benjamini–Hochberg adjusted p-value < 0.05, FDR q-value < 0.25, and |normalized enrichment score| > 1. Where applicable, Z-scores for significantly enriched terms were computed using the GOplot package for visualization.

Analysis of immune-cell enrichment in tumors

Immune and stromal components in the TCGA-BRCA cohort were evaluated using the ESTIMATE algorithm implemented in the R package estimate version 1.0.13. Log2(TPM + 1)-transformed expression data were used as input, and immune score, stromal score, and ESTIMATE score were calculated for each tumor sample. TIMER/TIMER2.0 was used to evaluate associations between MPO expression and estimated infiltration levels of major immune-cell populations in the TCGA-BRCA cohort, including B cells, CD8+ T cells, CD4+ T cells, macrophages, neutrophils, and dendritic cells^23,24,25. TIMER-based results were interpreted as immune-infiltration estimates derived from the corresponding online resource. For immune-cell enrichment analysis across 24 immune-cell types, single-sample gene set enrichment analysis (ssGSEA) was implemented using the R package GSVA version 1.46.0²⁶. The LM22 immune-cell signature matrix used for CIBERSORT-based deconvolution of 22 immune-cell types is provided in Supplementary Table 1. Correlations between MPO expression and immune-cell enrichment scores were evaluated using Spearman's rank correlation. Differences in immune-cell enrichment scores between the median-defined MPO-high and MPO-low tumor groups were compared using the Wilcoxon rank-sum test. For analyses involving multiple immune-cell types, p-values were adjusted using the Benjamini–Hochberg false discovery rate method.

DNA methylation of the MPO gene

DNA methylation patterns within the MPO locus were evaluated using MethSurv. CpG methylation beta values and survival associations for TCGA-BRCA were obtained from the MethSurv platform. Selected MPO-related CpG sites were visualized, and their associations with survival outcomes were evaluated using the survival analysis outputs provided by MethSurv²⁷. For analyses involving multiple CpG sites, p-values were adjusted across the tested MPO-related CpG sites using the Benjamini-Hochberg false discovery rate method. These methylation analyses were interpreted as exploratory epigenetic annotations.

PPI network construction and correlation analysis of neutrophil-related genes

To examine the association between MPO and neutrophil-related biology, a systematic network analysis was performed. A gene set comprising established mediators of neutrophil activation and associated inflammatory processes was curated from the current literature. The complete neutrophil-related gene list is provided in Supplementary Table 2. Gene symbols were harmonized to official gene symbols, duplicate entries were removed, and available genes were intersected with the TCGA-BRCA expression matrix before STRING/PPI analysis, hub-gene prioritization, and MPO–hub gene correlation analysis. The PPI network among these genes was constructed using the STRING database (version 11.5) with a medium confidence interaction score threshold (>0.40). Hub genes within this network were algorithmically prioritized based on degree centrality, which quantifies the number of direct interactions per node. The top 20 genes with the highest degree scores were selected for downstream correlation analysis.

Subsequently, the expression profiles of these hub genes and MPO were extracted from the TCGA-BRCA transcriptomic dataset. The association between MPO and each hub gene was statistically evaluated using Spearman's rank correlation. To characterize correlation patterns among the hub genes themselves, a pairwise Spearman correlation matrix was computed across all tumor samples. These correlation analyses provided the quantitative basis for subsequent visualizations, including the lollipop plot of MPO-hub gene correlations and the chord diagram/heatmap depicting hub-gene correlation patterns.

Prediction of upstream transcription factors and miRNAs targeting MPO

The KnockTF database (https://bio.liclab.net/KnockTF/index.php)^28,29, the ChIP database (http://chip-atlas.org/)^30,31, and the GTRD database^32,33 (https://gtrd.biouml.org/#!) were used to predict MPO's target TFs. In addition, the TargetScan database (https://www.targetscan.org/vert_80/) was utilized to predict potential miRNA binding sites targeting MPO. Venn diagrams were generated using the MicroBioinformatics website (https://www.bioinformatics.com.cn/static/others/jvenn/example.html)³⁴.

Single-cell analysis of MPO

The specific dataset GSE161529 originates from the Gene Expression Omnibus (GEO). Data preprocessing first performed cell-level filtering to exclude low-quality cells—those meeting any of the following criteria: mitochondrial gene expression exceeding 25%, total unique molecular identifier (UMI) count below 5000, or fewer than 2500 detected genes. Subsequently, ambient RNA contamination and technical batch effects were corrected³⁵. Principal component analysis (PCA) was performed for dimensionality reduction to assess cellular similarity, followed by UMAP for cell clustering and visualization. Then, according to the typical marker genes of cells, different clusters were annotated to cell types¹¹. The MPO-associated gene set used for single-cell signature scoring is provided in Supplementary File 1. Before scoring, gene symbols were harmonized to official gene symbols, duplicate entries were removed, and available genes were intersected with the GSE161529 expression matrix. AUCell, Seurat AddModuleScore, and ssGSEA were used to calculate per-cell MPO-associated scores. Scores from the three methods were Z-score normalized, scaled to a comparable range, and integrated to generate a composite MPO-associated score for downstream descriptive analyses. Cell-cell interaction networks were interrogated to compare inferred ligand–receptor communication patterns involving epithelial tumor cells stratified by MPO-associated signal and diverse partner cell types. These outputs were interpreted as descriptive communication patterns rather than evidence that MPO-expressing cells directly mediate intercellular communication.

Single-cell virtual knockdown of MPO and pathway enrichment analysis using scTenifoldKnk

Single-cell virtual knockdown of MPO was performed by integrating Seurat and scTenifoldKnk. Following standard quality control (200–6,000 genes per cell; mitochondrial fraction < 10%), data were log-normalized, and 2,000 highly variable genes were selected for dimensionality reduction and clustering. To enrich for MPO-relevant contexts, cells scoring in the top 50% for a myeloid/neutrophil gene module were retained. From these cells, an MPO-neighborhood subset was defined by expanding from MPO-positive seeds using k = 40 nearest neighbors in PCA space. The expanded subset was not treated as a pure MPO-positive population, and no cell-type proportion conclusion was drawn from this KNN expansion step. This subset was subjected to virtual knockdown analysis via scTenifoldKnk, using the union of highly variable genes and MPO (expressed in ≥25 cells) as the gene set. Significantly perturbed genes were identified (FDR < 0.05, BH-adjusted). Resultant genes were further analyzed for functional enrichment in GO Biological Processes and KEGG pathways (q < 0.05).

Exploratory drug-gene retrieval and ADMET annotation

DGIdb was queried to obtain preliminary MPO-associated drug–gene or chemical–gene interaction records. Because database-derived interaction lists may include entries supported by heterogeneous evidence types and may not directly correspond to clinically actionable therapeutic agents, the retrieved compounds were treated as exploratory annotations rather than prioritized treatment candidates. SwissADME and ADMETlab were subsequently used to summarize predicted physicochemical, pharmacokinetic, and toxicological properties. These in silico annotations were used to provide preliminary context for compound-level interpretation and to highlight the need for further pharmacological, toxicological, and clinical curation before any therapeutic relevance can be considered³⁶.

Results

MPO expression patterns and exploratory survival associations in breast cancer

To describe MPO expression patterns across cancer datasets, we analyzed MPO RNA-seq data from the TCGA pan-cancer dataset and observed lower MPO expression in tumor tissues from bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), glioblastoma multiforme (GBM), head and neck squamous cell carcinoma (HNSC), kidney chromophobe (KICH), liver hepatocellular carcinoma (LIHC), lung adenocarcinoma (LUAD), lung squamous cell carcinoma (LUSC), pancreatic adenocarcinoma (PAAD), prostate adenocarcinoma (PRAD), and thyroid carcinoma (THCA), and higher MPO expression in colon adenocarcinoma (COAD), kidney renal papillary cell carcinoma (KIRP), and other tissues (Figure 2A). We then evaluated associations between MPO expression and clinical outcomes in each cancer type. In the TCGA-BRCA cohort, both unpaired and paired comparisons showed lower MPO expression in tumor tissue than in normal/adjacent tissue (Figure 2B,C). After stratifying TCGA-BRCA tumor samples using the median tumor MPO expression cutoff, Kaplan-Meier analysis showed that patients with higher MPO expression had a longer progression-free interval (Hazard Ratio (HR) = 0.67, p = 0.028) (Figure 2D). Overall survival (OS) (p = 0.296; Supplementary Figure 1A) and disease-specific survival (DSS) (p = 0.18; Supplementary Figure 1B) were not statistically significant. The tumor-versus-normal ROC curve suggested separation between tissue groups in this dataset (Figure 2E), but this analysis should not be interpreted as clinical diagnostic validation. This exploratory discrimination may be influenced by normal-sample source, batch effects, tumor purity, and tissue-composition differences. MPO expression was also associated with pathological T stage (Figure 2F) and PAM50 subtype distribution (Figure 2G). Representative MPO immunohistochemistry (IHC) images of adjacent normal breast tissue and breast cancer tissue were included as qualitative protein-level references (Figure 2H). The boxed areas indicate regions shown at higher magnification. The 20× overview images include 100 µm scale bars, whereas the 40× higher-magnification images include 50 µm scale bars.

Correlation and enrichment analysis of MPO in the TCGA-BRCA cohort

Pearson correlation analysis identified the top 30 genes positively correlated with MPO, which showed coordinated upregulation along the MPO expression gradient (Figure 3A), whereas the top 30 negatively correlated genes exhibited an inverse expression pattern (Figure 3B). At the pathway level, MPO expression was significantly and positively associated with multiple tumor-related signature scores, including the inflammatory response signature (r = 0.41; Figure 3C), EMT markers (r = 0.264; Figure 3D), and the reactive oxygen species (ROS)-related gene set score (r = 0.415; Figure 3E), suggesting that MPO expression tracks with inflammatory/oxidative and mesenchymal-like transcriptional states in the TCGA-BRCA cohort.

Unsupervised clustering of MPO-associated genes further stratified tumors into expression patterns that aligned with clinical annotations, including pathological T stage and PAM50 intrinsic subtypes (Figure 3F). To explore possible connectivity among MPO-associated genes, we constructed a protein-protein interaction (PPI) network using STRING, revealing an interconnected module among several MPO-correlated genes (Figure 3G). In the PPI network, ESR1, FOXA1, XBP1, GATA3, and KRT18 had high network connectivity within this correlation-derived module. These results identify genes co-varying with MPO expression but do not establish MPO-related pathogenesis or directionality. Differential expression analysis between the MPO-high and MPO-low groups revealed transcriptomic differences summarized in the volcano plot (Figure 3H). A total of 1,159 upregulated and 854 downregulated genes were identified, providing input for subsequent enrichment analyses.

We next interrogated the functional relevance of the differentially expressed genes (DEGs) between the MPO-high and MPO-low groups using the clusterProfiler package in R. Gene Ontology (GO) enrichment analysis indicated that these DEGs were predominantly involved in immune-related biological processes, including regulation of immune response-related cell surface receptor signaling and lymphocyte-mediated immunity, with enrichment also observed in cellular components such as the T-cell receptor complex and molecular functions related to receptor activator activity (Figure 4A). Consistently, KEGG pathway analysis highlighted immune and inflammation-associated pathways, including cytokine–cytokine receptor interaction, chemokine signaling, T-cell receptor signaling, natural killer cell–mediated cytotoxicity, Th1/Th2 and Th17 differentiation, NF-κB signaling, primary immunodeficiency, and the intestinal immune network for IgA production (Figure 4B).

To further integrate expression directionality with functional terms, the GO plot was used to compute term-level Z-scores based on DEG |log2FC| values, which again highlighted immune-enriched transcriptional programs such as humoral immune response, leukocyte/lymphocyte-mediated immunity, immune response activation, and signaling transduction (Figure 4C). Gene set enrichment analysis (GSEA) based on the ranked gene list also showed enrichment of immune system pathways, including adaptive immune system, cytokine–cytokine receptor interaction, and neutrophil degranulation (Figure 4D–G). Because MPO is a myeloid/neutrophil-associated gene, these enrichments are interpreted as evidence that MPO-high samples exhibit stronger immune/myeloid transcriptional signals, rather than as evidence that MPO itself remodels the immune microenvironment.

Correlation between MPO expression and immune cell infiltration in breast cancer

We assessed the relationship between MPO expression and tumor microenvironment characteristics in the TCGA-BRCA cohort. Application of the ESTIMATE algorithm revealed significant positive correlations between MPO expression and the ESTIMATE score (R = 0.347, p < 0.001), immune score (R = 0.361, p < 0.001), and stromal score (R = 0.232, p < 0.001) (Figure 5A). The distribution of these scores across samples is shown in Figure 5B. Analysis using the TIMER/TIMER2.0 resource indicated that MPO expression was associated with estimated infiltration levels of major immune-cell populations, including B cells, CD8+ T cells, neutrophils, CD4+ T cells, macrophages, and dendritic cells in the TCGA-BRCA cohort (Figure 5C). This association pattern was further evaluated using ssGSEA-based immune-cell enrichment scores for 24 immune-cell types. After Benjamini–Hochberg false discovery rate correction, MPO expression showed positive associations with multiple immune-cell enrichment scores, including T cells, B cells, cytotoxic cells, dendritic-cell subsets, macrophages, T helper subsets, regulatory T cells, CD8+ T cells, NK cells, mast cells, and neutrophils (Figure 5D). These findings are interpreted as immune-composition associations rather than evidence that MPO directly controls immune-cell infiltration. A heatmap was generated to visualize sample-level immune-cell enrichment patterns across the TCGA-BRCA cohort (Figure 5E). We then compared ssGSEA-estimated immune-cell enrichment scores between median-defined MPO-high and MPO-low tumor groups. Several immune-cell enrichment scores differed between the two groups, including activated dendritic cells (aDC), B cells, CD8+ T cells, cytotoxic cells, neutrophils, T cells, Tregs, Th1 cells, Th2 cells, Th17 cells, γδ T cells, follicular helper T cells (TFH) highly variable gene (HVG) cells, effector memory T cells, central memory T cells, and T helper cells (Figure 5F,G). In addition, CIBERSORT-based deconvolution using the LM22 signature matrix was performed to estimate the relative fractions of 22 immune-cell types, and the resulting immune-cell composition patterns are shown in Figure 5H.

DNA methylation analysis of MPO in the TCGA-BRCA cohort

Using the same median tumor MPO expression cutoff, TCGA-BRCA samples were grouped into MPO-high and MPO-low groups, and DNA methylation patterns were visualized for each group (Figure 6A). Selected CpG sites within the MPO locus showed survival associations in the MethSurv analysis, including cg22331200, cg14619064, and cg11151395 (Figure 6B–G). These methylation-related results were interpreted as exploratory epigenetic annotations and require independent validation before prognostic or mechanistic conclusions can be drawn.

Association between MPO expression and neutrophil-related gene networks in breast cancer

The TCGA-BRCA cohort was used to examine the association between MPO expression and neutrophil-related genes. A STRING-based PPI network was constructed for neutrophil-associated genes, and hub genes were prioritized according to network topology (Figure 7A). The top 20 hub genes were subsequently evaluated for their correlation with MPO expression. As shown in the lollipop plot, MPO exhibited predominantly positive correlations with multiple neutrophil-related mediators, with stronger associations observed for chemokine/innate immune signaling components such as CCL5, CCL2, and TLR2, as well as TLR4, CXCR4, TNF, and MMP9 (Figure 7B).

To further characterize the co-regulation pattern among these hub genes, we visualized their pairwise relationships using a chord diagram and a correlation heatmap, which revealed extensive positive inter-gene correlations across the hub module, consistent with a coordinated inflammatory/neutrophil-associated transcriptional program (Figure 7C,D). Collectively, these results indicate that higher MPO expression is accompanied by coordinated expression of a neutrophil-related gene network in breast cancer.

Candidate transcription-factor annotation for MPO

To explore candidate transcription factors potentially associated with MPO, public transcription-factor resources, including KnockTF, ChIP-Atlas, and GTRD, were queried and intersected. Candidate transcription factors were further summarized using network-based prioritization and correlation analysis. A graphical summary is shown in Supplementary Figure 2, and the full tabular results are provided in Supplementary File 2. Because these databases integrate evidence from heterogeneous experimental contexts, database overlap and network degree were used only for candidate annotation and prioritization. These results were not interpreted as functional evidence of direct transcriptional regulation of MPO in breast cancer. Candidate factors, including MYC, are therefore presented as supplementary exploratory annotations rather than as validated upstream regulators.

Single-cell clustering and descriptive cell-cell communication analysis stratified by MPO signal

To annotate cell types, we first performed cluster-specific expression analysis based on canonical markers for each lineage. The average expression levels and the percentage of cells expressing these key genes across clusters are shown, supporting subsequent annotation (Figure 8A). Accordingly, the annotated cell clusters are visualized in a uniform manifold approximation and projection (UMAP) plot, in which each population is color-coded according to its identified type, including plasmacytoid dendritic cells, endothelial cells, myoepithelial cells, cycling epithelial cells, plasma cells, cytotoxic T cells, epithelial tumor cells, B cells, activated CD4 T cells, monocytes–macrophages, fibroblasts, and conventional T cells (Figure 8B). The heatmap displays expression levels of selected genes across cell clusters (C1-C8). Each row represents a gene, and each column represents a cell cluster. The color gradient indicates expression levels, with red denoting high expression and blue denoting low expression. The left dendrogram clusters genes with similar expression patterns (Figure 8C). The MPO-associated score was computed per cell using the MPO-associated gene set provided in Supplementary File 1. AUCell, Seurat AddModuleScore, and single-sample gene set enrichment analysis (ssGSEA) were used to calculate per-cell scores. Scores from the three methods were Z-score normalized, scaled to a comparable range, and integrated to yield a composite MPO-associated score for downstream descriptive analysis (Figure 8D).

This cell–cell interaction analysis compared inferred ligand–receptor communication patterns between cell groups stratified by MPO-associated signal, including the interaction network, signaling-pattern heatmaps, outgoing signaling bubble plot, and incoming signaling bubble plot (Figure 8E–H). Because MPO signal was sparse at the single-cell level and the apparent distribution across annotated cell types may be affected by dropout, ambient RNA, doublets, and annotation uncertainty, these communication plots should be interpreted as descriptive workflow outputs. They do not demonstrate that MPO-expressing cells mediate or control intercellular communication. Detectable MPO signal was observed in a limited number of cells, including annotated epithelial tumor cells and monocytes–macrophages (Figure 8I). Given that MPO is canonically associated with neutrophil/myeloid lineages, this pattern requires validation in independent single-cell datasets or by orthogonal experimental methods.

Exploratory scTenifoldKnk sensitivity analysis based on sparse MPO-positive cells

Multiple 10x Genomics samples were integrated, followed by normalization and HVG selection, PCA-based dimensionality reduction, construction of a k-nearest neighbor graph, and Louvain clustering. Canonical marker-gene expression patterns across clusters were summarized using a DotPlot, supporting subsequent cell-type annotation (Figure 9A). UMAP visualization showed the annotated single-cell populations in the integrated dataset (Figure 9B). Canonical lineage markers (e.g., EPCAM and KRT8/KRT18 for epithelial cells; PTPRC for immune cells; MS4A1 for B cells; LST1/S100A8/S100A9 for myeloid cells; PECAM1 for endothelial cells; and COL1A1 for fibroblast/smooth muscle lineages) showed cluster-specific expression patterns, supporting cell-type annotation (Figure 9C). Sample-stratified stacked bar plots indicated that each sample contained multiple clusters with limited overall batch-to-batch variation (Figure 9D).

MPO expression was relatively sparse in the single-cell dataset, with only 85 MPO-positive cells detected initially (Figure 9E). Given this limited number, KNN-based neighborhood expansion was used only to define a local MPO-neighborhood subset for exploratory sensitivity analysis. This expanded subset should not be interpreted as a pure MPO-positive population because it may include neighboring cells with low or undetectable MPO expression. Within this MPO-neighborhood subset, virtual knockdown of MPO was performed using scTenifoldKnk as a computational sensitivity analysis. The resulting volcano plot, manifold displacement analysis, manifold alignment visualization, GO/KEGG enrichment results, and top-displacement genes (Figure 9F-N) highlighted candidate transcriptional programs related to antigen presentation, myeloid/lymphocyte activation, cytokine production, and phagosome-related pathways. These results should be interpreted as exploratory transcriptional sensitivity outputs rather than direct evidence that MPO mechanistically regulates these pathways in breast cancer. Independent single-cell datasets and orthogonal experimental validation, such as immunohistochemistry, flow cytometry, qPCR, or functional assays, will be required to substantiate these observations.

Exploratory drug–gene interaction retrieval and ADMET annotation

As an exploratory extension of the MPO-centered analysis, drug–gene interaction information was retrieved from DGIdb. A graphical summary is shown in Supplementary Figure 3, and compound-level results are provided in Supplementary Table 3. The DGIdb query returned a heterogeneous set of MPO-associated chemical entries, including compounds with limited clinical plausibility or unfavorable toxicological profiles. Therefore, these database-derived compounds were not considered therapeutic candidates for breast cancer based on the present analysis. ADMET-related information was summarized to provide preliminary annotation of predicted physicochemical, pharmacokinetic, and toxicological properties. Database-based compound retrieval and ADMET annotation are not equivalent to clinically curated drug prioritization. These results therefore serve only as screening-level chemical annotations and illustrate the need for careful pharmacological, toxicological, and clinical filtering before any compound can be considered for therapeutic investigation. The main findings of this study focus on the association between MPO expression and immune/myeloid-related transcriptional features.

DATA AVAILABILITY:

TCGA-BRCA transcriptomic and clinical data were obtained from the Genomic Data Commons portal (https://portal.gdc.cancer.gov; downloaded on 26 August 2025; data release/version 202208). The single-cell dataset GSE161529 was obtained from the Gene Expression Omnibus (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE161529). No new sequencing data were generated in this study. Analysis scripts are publicly available at https://github.com/tengfeitcm/MPO.

Figure 1: Flowchart of the data collection and analysis process. Please click here to view a larger version of this figure.

Figure 2: MPO expression patterns and exploratory survival associations in breast cancer. (A) MPO expression levels were analyzed in 33 distinct cancer types and their adjacent normal tissues using the TCGA database. (B) Unpaired samples were selected from the TCGA-BRCA dataset to analyze MPO mRNA expression in breast cancer and normal tissues. (C) Paired samples were selected from the TCGA-BRCA dataset to analyze MPO mRNA expression in breast cancer and normal tissues. (D) Kaplan-Meier analysis of PFI in patients stratified by the median tumor MPO expression cutoff in the TCGA-BRCA cohort. (E) Exploratory ROC curve evaluating tumor-normal discrimination based on MPO expression in the analyzed public transcriptomic dataset. (F) MPO expression across different pathological T stages. (G) MPO expression across PAM50 molecular subtypes, with subtype labels shown. (H) Representative MPO immunohistochemistry (IHC) images of adjacent normal breast tissue and breast cancer tissue. The boxed areas indicate regions shown at higher magnification. The 20× overview images include 100 µm scale bars, whereas the 40× higher-magnification images include 50 µm scale bars. These images are shown as qualitative protein-level references and were not used for quantitative morphometric or statistical analysis. Please click here to view a larger version of this figure.

Figure 3: MPO-associated correlation and differential-expression analysis in breast cancer. (A) Top 30 of coding genes positively correlated with MPO expression at the mRNA level based on Pearson correlation coefficients from the TCGA database. (B) Top 30 coding genes negatively correlated with MPO expression at the mRNA level based on Pearson correlation coefficients. (C) Scatter plots illustrating Spearman correlations between MPO and genes upregulated by inflammatory response. (D) Scatter plots illustrating Spearman correlations between MPO and genes upregulated by EMT markers. (E) Scatter plots illustrating Spearman correlations between MPO and genes upregulated by ROS. (F) Heatmap of MPO-associated gene clusters based on clinical significance (T stage and PAM50). (G) PPI network predicted using the STRING database for MPO-associated proteins. (H) Volcano plot of differentially expressed genes between median-defined MPO-high and MPO-low tumor groups in the TCGA-BRCA cohort. Please click here to view a larger version of this figure.

Figure 4: Enrichment analysis of MPO in breast cancer. (A) Gene Ontology enrichment analysis of the 2,013 differentially expressed genes between MPO-high and MPO-low groups. (B) Kyoto Encyclopedia of Genes and Genomes pathway enrichment analysis of the 2,013 differentially expressed genes. (C) Combined Gene Ontology enrichment visualization integrating enriched terms with differential-expression direction and |log2FC| values. (D) Representative GSEA enrichment plot for an MPO-associated immune-related gene set; the gene-set name, normalized enrichment score, and FDR q-value are shown in the panel. (E) Representative GSEA enrichment plot for an additional MPO-associated immune-related gene set; the gene-set name, normalized enrichment score, and FDR q-value are shown in the panel. (F) Representative GSEA enrichment plot for an additional MPO-associated immune-related gene set; the gene-set name, normalized enrichment score, and FDR q-value are shown in the panel. (G) Representative GSEA enrichment plot for an additional MPO-associated immune-related gene set; the gene-set name, normalized enrichment score, and FDR q-value are shown in the panel. Please click here to view a larger version of this figure.

Figure 5: Correlation between immune-cell enrichment and MPO expression in breast cancer. (A) Scatter plots showing correlations between MPO expression and ESTIMATE score, immune score, and stromal score. (B) Box plots showing differences in ESTIMATE score, immune score, and stromal score between median-defined MPO-high and MPO-low tumor groups. (C) TIMER/TIMER2.0-based analysis showing associations between MPO expression and estimated infiltration of major immune-cell populations. (D) Lollipop plot showing Spearman correlations between MPO expression and ssGSEA-estimated enrichment scores for 24 immune-cell types. P-values from multiple immune-cell correlations were adjusted using the Benjamini–Hochberg false discovery rate method. (E) Heatmap illustrating sample-level immune-cell enrichment patterns across the TCGA-BRCA cohort. (F) Box plots showing the first set of ssGSEA-estimated immune-cell enrichment score differences between median-defined MPO-high and MPO-low tumor groups; group comparisons were performed using the Wilcoxon rank-sum test with Benjamini–Hochberg correction. (G) Box plots showing the second set of ssGSEA-estimated immune-cell enrichment score differences between median-defined MPO-high and MPO-low tumor groups; group comparisons were performed using the Wilcoxon rank-sum test with Benjamini–Hochberg correction. (H) Stacked bar plot showing CIBERSORT-estimated immune-cell fractions based on the LM22 signature matrix for 22 immune-cell types in median-defined MPO-low and MPO-high tumor groups. Please click here to view a larger version of this figure.

Figure 6: DNA methylation analysis of the MPO gene in breast cancer. (A) Heatmap showing MPO methylation patterns in median-defined MPO-high and MPO-low groups. (B) Kaplan-Meier survival curve demonstrating the prognostic significance of methylation at the cg27456487 site. (C) Kaplan-Meier survival curve demonstrating the prognostic significance of methylation at the cg02668773 site. (D) Kaplan-Meier survival curve demonstrating the prognostic significance of methylation at the cg07110356 site. (E) Kaplan-Meier survival curve demonstrating the prognostic significance of methylation at the cg11151395 site. (F) Kaplan-Meier survival curve demonstrating the prognostic significance of methylation at the cg14619064 site. (G) Kaplan-Meier survival curve demonstrating the prognostic significance of methylation at the cg22331200 site. Please click here to view a larger version of this figure.

Figure 7: Analysis of MPO and neutrophil-related gene correlations at the mRNA level using the TCGA database. (A) Visualization of the protein interaction network, depicting interactions between the core protein and other proteins. (B) Correlation analysis of the top 20 neutrophil-related genes with MPO, showing correlation coefficients and P-value distributions for different genes. (C) Chord diagram of correlations among the top 20 neutrophil-related genes, visually representing the strength and direction of gene associations. (D) Correlation heatmap of the top 20 neutrophil-related genes, displaying correlation coefficients and significance levels through color gradients and statistical markers. Please click here to view a larger version of this figure.

Figure 8: Single-cell clustering and MPO-associated cell–cell communication analysis in the breast cancer single-cell dataset. (A) DotPlot of canonical marker genes across clusters for cell-type annotation. (B) UMAP visualization of annotated cell populations. (C) Heatmap of selected marker genes across cell clusters. (D) DotPlot summarizing MPO-associated scores across annotated cell types computed using AUCell, ssGSEA, and Seurat AddModuleScore based on the gene set provided in Supplementary File 1. (E) Cell-cell interaction network depicting communication between epithelial tumor cells stratified by MPO-associated signal and other cell types; edge width represents interaction strength and node size reflects overall interaction activity. (F) Heatmaps showing outgoing and incoming signaling patterns across cell types. (G) Bubble plot of outgoing signaling pathways from epithelial tumor cells stratified by MPO-associated signal to other cell types. (H) Bubble plot of incoming signaling pathways from other cell types toward epithelial tumor cells stratified by MPO-associated signal. (I) MPO expression distribution across annotated cell types. Please click here to view a larger version of this figure.

Figure 9: Single-cell atlas analysis and exploratory virtual knockdown of MPO sensitivity output. (A) DotPlot showing the expression of canonical marker genes across single-cell clusters; dot size represents the percentage of cells expressing each marker, and color intensity represents average expression level. (B) UMAP visualization of annotated single-cell populations, with each color representing a distinct cell type or cluster. (C) UMAP visualization of key marker gene expression, showing expression distribution of marker genes for cell types, including myeloid cells. (D) Stacked bar chart of cell cluster proportions across samples. (E) UMAP visualization of MPO gene expression. (F) Violin plot showing single-cell sequencing QC metrics. (G) Clustering plot of key marker genes. (H) DotPlot of canonical marker genes at the cluster level. (I) Volcano plot of genes changed in the virtual knockdown sensitivity analysis. (J) Scatter plot of displacement versus significance. (K) Manifold alignment arrow plot. (L) GO BP enrichment analysis of genes from the virtual knockdown output. (M) KEGG pathway enrichment analysis of genes from the virtual knockdown output. (N) Top 20 genes with the highest manifold displacement after excluding MPO. Please click here to view a larger version of this figure.

Supplementary Figure 1: Additional survival analyses for MPO in the TCGA-BRCA cohort. (A,B) This file contains supplementary Kaplan-Meier survival analyses for (A) overall survival and (B) disease-specific survival stratified by the median tumor MPO expression cutoff. These analyses are provided as supplementary outcome analyses to Figure 2D and were not statistically significant in the current cohort.Please click here to download this file.

Supplementary Figure 2: Exploratory candidate transcription-factor annotation for MPO. (A) Venn diagram showing the intersection of candidate transcription factors from three public transcription-factor resources. (B) MYC expression comparison output. (C) Transcription-factor correlation heatmap with row and column labels. (D) MPO–MYC correlation output. (E) MYC survival-analysis output. (F) MYC ROC output. MYC-related outputs are shown only as supplementary candidate-transcription-factor annotations and are not used to support mechanistic upstream-regulator conclusions.Please click here to download this file.

Supplementary Figure 3: Exploratory DGIdb drug-gene retrieval output for MPO. Gray nodes represent the MPO gene, orange nodes represent retrieved small-molecule entries, and connecting lines indicate database-predicted drug–gene relationships.Please click here to download this file.

Supplementary Table 1: The LM22 immune-cell signature matrix used for CIBERSORT-based immune-cell deconvolution analysis of 22 immune-cell types. Gene symbols were harmonized, duplicate entries were removed, and available genes were intersected with the corresponding TCGA-BRCA or GSE161529 expression matrices before downstream analysis.Please click here to download this file.

Supplementary Table 2: The neutrophil-related gene list used for STRING/PPI analysis, hub-gene prioritization, and MPO–hub gene correlation analysis. Please click here to download this file.

Supplementary Table 3: Exploratory DGIdb drug–gene retrieval and ADMET annotation outputs for MPO. This file contains DGIdb-retrieved MPO-associated chemical-gene interaction records and compound-level predicted physicochemical, pharmacokinetic, and toxicity-related annotations. These outputs are provided only as preliminary chemical annotations and should not be interpreted as therapeutic candidate lists. They do not establish MPO inhibition, target engagement, ligand specificity, selectivity, safety, therapeutic efficacy, or clinical suitability. The values in this table represent predicted physicochemical and drug-likeness parameters for the listed compounds. Molecular weight is expressed in grams per mole (g/mol). Hydrogen-bonding acceptor and hydrogen-bonding donor values indicate the predicted numbers of hydrogen-bond acceptors and donors, respectively. The Moriguchi octanol–water partition coefficient indicates predicted lipophilicity. Lipinski's violations indicate the number of Lipinski rule-of-five criteria not satisfied by each compound. The bioavailability score represents the predicted oral bioavailability-related score, and topological surface area refers to the predicted topological polar surface area.Please click here to download this file.

Supplementary File 1: The MPO-associated gene list used for single-cell signature scoring with AUCell, Seurat AddModuleScore, and ssGSEA. Please click here to download this file.

Supplementary File 2: Exploratory candidate transcription-factor and miRNA annotation outputs for MPO. This file contains database-derived candidate transcription-factor and miRNA annotation results based on public resources, including KnockTF, ChIP-Atlas, GTRD, and TargetScan. These annotations are provided only for exploratory candidate prioritization and are not interpreted as functional evidence of upstream regulation of MPO in breast cancer.Please click here to download this file.

Discussion

This study presents an exploratory public-dataset and in silico workflow for examining associations between MPO expression and immune/myeloid features in breast cancer. The TCGA-BRCA analyses showed that MPO expression was lower in tumor tissues than in adjacent non-tumor tissues and that higher MPO expression was associated with a longer progression-free interval. However, overall survival and disease-specific survival were not statistically significant. Therefore, MPO should not be interpreted as a robust or established prognostic biomarker based on the current evidence. Future studies should evaluate MPO using multivariable Cox regression models adjusted for established clinicopathological variables, independent validation cohorts, and subtype-stratified analyses.

Compared with conventional single-cohort differential-expression analysis or single-platform immune-infiltration estimation, this MPO-centered workflow integrates bulk transcriptomics, immune enrichment, methylation annotation, single-cell mapping, and virtual perturbation to provide a broader exploratory view of MPO-associated immune/myeloid features. However, this workflow remains complementary to, rather than a replacement for, external cohort validation, spatial or protein-level validation, and experimental perturbation assays.

The immune-infiltration and enrichment results should be interpreted as MPO-associated immune context rather than MPO-driven immune remodeling. MPO is predominantly expressed in neutrophils and other myeloid-lineage cells³⁷. Therefore, positive correlations between MPO expression and ESTIMATE scores, immune scores, immune-cell enrichment scores, neutrophil-related genes, cytokine pathways, antigen-presentation signatures, or neutrophil degranulation pathways are biologically plausible and may largely reflect differences in immune/myeloid cell abundance within bulk tumor samples. This interpretation is consistent with prior studies showing that MPO-positive neutrophil infiltration is associated with a favorable prognosis in breast cancer and that MPO has been implicated in dendritic-cell function and T-cell-driven tissue inflammation^38,39. Bulk RNA-seq data cannot determine whether MPO has tumor-cell-intrinsic activity or whether the observed signal primarily reflects infiltrating immune cells. Independent single-cell datasets, spatial profiling, immunohistochemistry, flow cytometry, or perturbation-based experimental models would be needed to clarify the cellular source and function.

The single-cell analysis provides additional descriptive information but remains limited by sparse MPO detection. Only 85 MPO-positive cells were initially detected before KNN-based neighborhood expansion. Although KNN expansion allowed a sensitivity analysis of cells in the local transcriptional neighborhood of MPO-positive cells, this procedure may include cells that do not directly express MPO. Accordingly, the scTenifoldKnk virtual knockdown output should be interpreted as an exploratory computational sensitivity analysis rather than evidence of MPO-mediated pathway regulation⁴⁰. Validation in independent single-cell breast cancer datasets and orthogonal experimental assays will be required before mechanistic conclusions can be drawn.

The transcription-factor analysis should also be interpreted cautiously. The overlap of KnockTF, GTRD, and ChIP-Atlas predictions, followed by degree-based prioritization, can nominate candidate transcription factors but cannot establish functional transcriptional regulation of MPO in breast cancer. MYC and other candidate factors were therefore retained only as exploratory annotations. Because transcription-factor activity is highly context-dependent and may vary by tumor subtype, cellular composition, assay platform, and preprocessing strategy, context-specific validation would be required before any candidate factor can be assigned an upstream regulatory role. Such validation should include ChIP-qPCR or ChIP-seq, promoter reporter assays, and transcription-factor perturbation followed by measurement of MPO expression.

The drug–gene retrieval and ADMET annotation should also be interpreted cautiously. DGIdb can return heterogeneous chemical–gene associations, including compounds that are not selective MPO ligands and may have limited clinical plausibility or unfavorable toxicological properties³⁶. ADMET predictions provide preliminary chemical annotations but do not establish target engagement, potency, selectivity, safety, or therapeutic efficacy⁴¹. Therefore, the current compound-level results should not be used to infer therapeutic potential. A meaningful translational evaluation would require a curated set of pharmacologically relevant MPO inhibitors or probes, comparison with established MPO-targeting compounds, and validation using biochemical, cellular, and pharmacological assays. The context-dependent and potentially dual role of MPO in cancer should also be considered when interpreting drug-related findings. MPO may contribute to tumor-promoting processes through oxidative stress, reactive oxidant generation, DNA damage, chronic inflammation, and remodeling of the tumor microenvironment. At the same time, MPO expression in bulk tumor datasets may reflect infiltration by neutrophils or other myeloid immune cells, which in some contexts may be associated with an immune-active microenvironment and more favorable clinical outcomes^39,42. The biological interpretation of MPO, therefore, depends on tumor type, cellular source, disease stage, and immune-cell composition.

DNA methylation findings within the MPO locus were also considered exploratory. Selected CpG sites showed survival associations in the MethSurv analysis, but these results require independent validation before prognostic or mechanistic conclusions can be drawn. Epigenetic regulation of MPO may interact with transcription-factor binding and chromatin-level regulation, but such interactions remain speculative without functional chromatin or perturbation data⁴³.

The clinical implications of MPO expression should be interpreted cautiously. The present results do not establish MPO as a clinically actionable biomarker or as a marker that can currently guide immunotherapy decisions in breast cancer. Rather, MPO may reflect a myeloid/neutrophil-related immune context within the tumor microenvironment. In future studies, MPO expression could be evaluated together with established immunotherapy-related markers, including tumor-infiltrating lymphocytes, PD-L1 expression, immune checkpoint gene expression, molecular subtype, and validated immune signatures. Such analyses should include independent cohorts, multivariable models, and treatment-response datasets before MPO can be considered for patient stratification or immunotherapy decision-making.

Several workflow steps are critical for reproducibility, including consistent TCGA-BRCA data preprocessing, tumor-only median-based MPO-high/MPO-low grouping, predefined statistical thresholds and multiple-testing correction, immune-cell enrichment algorithms and signature sets, single-cell quality control and annotation, KNN-based MPO-neighborhood expansion, and the exploratory handling of virtual knockdown and DGIdb/ADMET outputs. Changes in these parameters may affect downstream results and interpretation; therefore, they should be reported and reproduced carefully. For troubleshooting, inconsistent outputs should be addressed by checking sample source, expression normalization, grouping cutoff, multiple-testing correction, immune-cell signature sets, single-cell QC and annotation, KNN-neighborhood definition, virtual knockdown thresholds, enrichment cutoffs, and heterogeneous DGIdb/ADMET compound records.

Several limitations of this study should be acknowledged. First, this study was based on retrospective analyses of public databases using TCGA-BRCA and publicly available single-cell data, and therefore may be affected by cohort heterogeneity, sample-source differences, batch effects, incomplete clinical annotation, and tumor-composition differences. Second, the present findings are derived mainly from transcriptomic associations and in silico analyses and lack direct experimental validation. Therefore, the observed associations between MPO expression, immune/myeloid features, methylation annotations, transcription-factor candidates, and virtual knockdown outputs should not be interpreted as causal mechanisms. Third, MPO detection in the single-cell dataset was sparse, with only 85 MPO-positive cells detected before KNN-based neighborhood expansion. The expanded MPO-neighborhood subset may include cells with low or undetectable MPO expression and should not be considered a pure MPO-positive population. Fourth, because MPO is predominantly associated with neutrophils and other myeloid-lineage cells, MPO-related signals in bulk RNA-seq data may be confounded by immune-cell abundance, tumor purity, and cellular composition rather than reflecting tumor-cell-intrinsic activity. Finally, the DGIdb-based drug–gene retrieval and ADMET annotation were exploratory chemical annotations only. These outputs do not establish MPO inhibition, target engagement, selectivity, safety, therapeutic efficacy, or clinical suitability. Future studies using independent cohorts, spatial or protein-level validation, and functional experiments are required to confirm the biological and clinical relevance of these findings.

In summary, the current public-dataset analyses support an association between MPO expression and immune/myeloid transcriptional features in breast cancer, with higher MPO expression associated with longer progression-free interval in the analyzed cohort. The findings remain exploratory and hypothesis-generating. They do not establish that MPO causally regulates the tumor immune microenvironment, that MYC functionally regulates MPO, or that retrieved compounds have therapeutic relevance. The main contributions of this study are a reproducible computational workflow and a set of testable hypotheses that require external cohort validation and experimental confirmation.

Materials

List of materials used in this article

Name	Company	Catalog Number	Comments
CellChat	R package/Open source	https://github.com/sqjin/CellChat	Cell-cell communication analysis
ChIP-Atlas	Public database	https://chip-atlas.org/	TF target screening; 2021 update
clusterProfiler	Bioconductor	https://bioconductor.org/packages/clusterProfiler/	GO/KEGG enrichment analysis; v4.4.4
Cytoscape	Cytoscape Consortium	https://cytoscape.org/	Network visualization and topology analysis
DGIdb	Washington University/Public database	https://www.dgidb.org/	Drug-gene interaction retrieval
GDC/TCGA-BRCA	National Cancer Institute	https://portal.gdc.cancer.gov/	Bulk transcriptomic and clinical data source
Gene Expression Omnibus: GSE161529	NCBI	https://www.ncbi.nlm.nih.gov/geo/	Single-cell dataset source
GSEA/MSigDB	Broad Institute	https://www.gsea-msigdb.org/gsea/msigdb	Gene set enrichment analysis and gene-set reference; Version 3.0
GSVA	Bioconductor	https://bioconductor.org/packages/GSVA/	Gene set variation/ssGSEA-related scoring; Version 1.46.0
GTRD	Public database	http://gtrd.biouml.org/	TF target screening; 2021
KnockTF	Public database	http://www.licpathway.net/KnockTF/index.html	TF perturbation resource; Version 2.0
R	R Foundation for Statistical Computing	https://www.r-project.org/	Statistical computing environment
scTenifoldKnk	R package/Open source	https://github.com/cailab-tamu/scTenifoldKnk	Virtual knockdown analysis
Seurat	R package/Open source	https://satijalab.org/seurat/	Single-cell preprocessing and clustering
STRING	ELIXIR/Public database	https://string-db.org/	Protein-protein interaction analysis; v11
SwissADME	SIB Swiss Institute of Bioinformatics	http://www.swissadme.ch/	Drug-likeness assessment; 2017 release/web tool
TIMER	Public web resource	https://timer.cistrome.org/	Immune infiltration analysis; TIMER2.0
UCSC Xena or linked TCGA portal	UCSC	https://xenabrowser.net/	Exploratory data access/validation