Mitophagy-Driven Immune Cell Infiltration Patterns Define Breast Cancer Subtypes with Differential Treatment Responses

Mitophagy-Driven Immune Cell Infiltration Patterns Define Breast Cancer Subtypes with Differential Treatment Responses | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Research Article Mitophagy-Driven Immune Cell Infiltration Patterns Define Breast Cancer Subtypes with Differential Treatment Responses Wenbin Guo, Zhiqiang Ye, Dan Wu This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-6960401/v1 This work is licensed under a CC BY 4.0 License Status: Posted Version 1 posted You are reading this latest preprint version Abstract Background The development of breast cancer (BC) entails intricate immunological and molecular interactions. Although mitophagy-related genes (MRGs) are essential for maintaining cellular homeostasis, little is known about their expression patterns, immunological effects, and therapeutic applications in BC. The purpose of this work was to describe immunological microenvironment changes, MRG dysregulation, and their predictive modeling potential in BC. Methods Using data from single-cell RNA sequencing (scRNA-seq) (GSE248288) and The Cancer Genome Atlas (TCGA), we examined immune cell infiltration MRG expression profiles, and functional pathways. Machine learning methods identified significant MRGs for predictive modeling, which were confirmed using qPCR in BC cell lines. Subtypes were identified using hierarchical clustering, and weighted gene co-expression network analysis (WGCNA) showed subtype-specific modules. Drug sensitivity and immunotherapy responses were evaluated using the oncoPredict and Tumor Immune Dysfunction and Exclusion (TIDE) algorithms. Results Tumor tissues showed elevated MRGs (e.g., SRC, PGAM5, FUNDC1) and altered immunological infiltration, with more activated T cells and macrophages but fewer naïve B cells. scRNA-seq demonstrated increased MRG activity in NK and B cells, which is connected to PI3K-AKT-mTOR signaling and allograft rejection. A nine-MRG predictive model (PGAM5, PRKN, TOMM40, FUNDC1, MAP1LC3B, PINK1, MTERF3, CSNK2A2, and SRC) demonstrated great diagnostic accuracy (AUC: 0.99). BC subgroups based on MRGs expression displayed unique molecular profiles: Cluster 1 (cell cycle dysregulation), Cluster 2 (cilia-related pathways), and Cluster 3 (small GTPase signaling). Cluster 2 showed potential immunotherapy responsiveness (low TIDE scores), whereas Cluster 1 was sensitive to chemotherapy (paclitaxel, gemcitabine). Conclusions MRGs are important regulators of breast cancer progression, regulating immunological dynamics, metabolic pathways, and treatment responses. The established nine-MRG model accurately predicts BC occurrence, whereas subtype-specific biochemical and immunological aspects provide insights for tailored therapy. These findings emphasize MRGs as possible biomarkers for diagnosis and personalized therapy options in BC. Mitophagy breast cancer Single-Cell Sequencing Immune Cells Immunotherapy Response Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Introduction Breast cancer remains one of the most prevalent malignancies affecting women worldwide, contributing significantly to cancer-related morbidity and mortality. Despite advancements in early detection and treatment, the high incidence of therapeutic resistance and metastasis poses substantial challenges in breast cancer management. A critical aspect of breast cancer pathology is the involvement of mitochondrial dysfunction and the subsequent process of mitophagy, which has garnered increasing attention for its role in cancer progression and treatment resistance[ 1 – 5 ]. Mitochondria are essential organelles involved in various cellular processes, including energy production, apoptosis, and redox homeostasis. Mitochondrial dysfunction is a hallmark of cancer cells, contributing to tumorigenesis and cancer progression by altering cellular metabolism and promoting resistance to apoptosis[ 6 – 10 ]. In breast cancer, oxidative stress-induced damage to mitochondrial DNA (mtDNA) and the resultant mitochondrial dysfunction play crucial roles in cancer cell survival and proliferation[ 1 ]. The regulation of mitochondrial dynamics and the maintenance of mitochondrial integrity through processes such as mitophagy are essential for cancer cells to adapt to their rapidly changing environment. Mitophagy, an evolutionarily conserved cellular process, plays a dual role in cancer by acting as both a tumor suppressor and promoter depending on the context and stage of cancer development [ 6 , 11 ]. Dysregulated mitophagy can lead to the accumulation of damaged mitochondria, thereby contributing to carcinogenesis and tumor progression. In breast cancer, mitophagy has been shown to facilitate the metabolic reprogramming of cancer cells, promoting a switch between aerobic glycolysis and oxidative phosphorylation to meet the high energy demands of proliferating tumor cells [ 11 ]. Key regulators of mitophagy, including PINK1, Parkin, BNIP3, and NIX , have been implicated in various cancers [ 11 – 14 ]. These proteins help maintain mitochondrial quality by targeting damaged mitochondria for autophagic degradation. The interplay between these mitophagy regulators and the broader mitochondrial dynamics is critical for cancer cell survival, particularly under stress conditions such as hypoxia and nutrient deprivation commonly found in the tumor microenvironment. The role of mitophagy in breast cancer is increasingly recognized. The field of mitoepigenetics, which focuses on the epigenetic regulation of mtDNA, has emerged as a significant area of research[ 1 ]. Mutational damage to mtDNA and defects in mitochondrial biogenesis and dynamics contribute to the oncogenic phenotype of breast cancer cells. Studies have shown that cancer cells often exhibit defects in mtDNA repair mechanisms, leading to increased reactive oxygen species (ROS) production and oxidative stress, which in turn promotes further mitochondrial damage and dysfunction. This creates a vicious cycle that fosters cancer cell survival and proliferation. Therapeutic targeting of these mitochondrial vulnerabilities, including the enhancement of mitophagy, represents a promising strategy for cancer treatment. Several mitochondrial-targeting agents, such as biguanides and oxidative phosphorylation (OXPHOS) inhibitors, are being explored for their potential to disrupt this cycle and inhibit breast cancer progression[ 1 , 15 ]. Targeting mitophagy in cancer therapy holds significant promise. By modulating mitophagy, it may be possible to eliminate cancer cells more effectively and overcome resistance to conventional therapies. For instance, the use of mitophagy inducers or inhibitors could help maintain mitochondrial homeostasis and enhance the efficacy of existing treatments. Furthermore, the identification of biomarkers related to mitophagy and mitochondrial dysfunction could aid in the early diagnosis and personalized treatment of breast cancer[ 11 ]. In this study, We undertook a comprehensive investigation of mitophagy-related gene expression in breast cancer. We then created a predictive model for reliable breast cancer diagnosis based on mitophagy-associated signature genes. This model was validated using a subset of the TCGA dataset. Breast cancer subgroups were identified using hierarchical clustering. In addition, we conducted extensive research into the prognosis, immunological infiltration, immunotherapy response, and drug sensitivity of each discovered breast cancer subtype. Materials and methods 2.1 Data acquisition and processing The Cancer Genome Atlas (TCGA) database provided RNA sequencing (RNA-seq) data and related clinical details for 1,119 breast cancer (BC) samples. The R tool "TCGAbiolinks"(version = 2.32.0) was used to preprocess and standardize the collected data. A training dataset and a test dataset were randomly selected from the BRCA samples at a 7:3 ratio. The remaining 336 samples were used for internal validation, leaving 783 samples in the training cohort. We used single-cell RNA-seq data from the Gene Expression Omnibus (GEO) database, accession number GSE248288, to further explore the function of mitophagy in breast cancer[ 16 ]. Additionally, a comprehensive list of mitophagy-related genes was retrieved from the Reactome database and cross-referenced with the GeneCards database ( https://www.genecards.org/ ) to ensure the inclusion of all relevant genes in our analysis. 2.2 Differential expression analysis of mitophagy-related genes To investigate the expression patterns of mitophagy-related genes (MRGs) in breast cancer, we first classified the samples from TCGA datasets into normal and tumor groups based on the available clinical information. The MRGs were obtained from the Reactome database, and their expression levels were visualized using a heatmap generated with the R package "pheatmap"(version = 1.0.12). Next, we performed differential expression analysis to identify differentially expressed genes (DEGs) between the normal and tumor samples (the threshold was set as log2 fold change > 0 or < 0 and adjusted p value < 0.05). The differential expression analysis was conducted using the R package "DESeq2"(version = 1.44.0), which employs a negative binomial distribution to model the read counts and provides a robust statistical framework for identifying DEGs. To visualize the results of the differential expression analysis, we constructed boxplots and volcano plots using the R packages "ggplot2"(version = 3.5.1) and "ggpubr"(version = 0.6.0). The boxplots allowed for a comparison of the expression levels of MRGs between normal and tumor samples, while the volcano plots provided an overview of the distribution of DEGs based on their log2fold change and adjusted p-values. Furthermore, to explore the relationships among the identified DEGs, we generated a correlation heatmap using the R package "corrplot"(version = 0.94). This heatmap visualized the pairwise correlations between the DEGs, enabling us to identify potential co-expression patterns and functional relationships among the MRGs in breast cancer. 2.3 Analysis of scRNA-seq data in BC To analyze the single-cell RNA sequencing (scRNA-seq) data in breast cancer (BC), we utilized the R package "Seurat"[ 17 , 18 ](version = 5.1.0). To assure high-quality data, we performed a series of filtering operations on the raw matrix. We specifically kept genes expressed in at least three cells, deleted cells expressing fewer than 200 genes, and excluded cells with mitochondrial gene expression greater than 25%. The scRNA-seq data was normalized using the "NormalizeData" function, which was then used to identify the top 2,000 highly variable genes using the "FindVariableFeatures" function. Principal component analysis (PCA) was then done with the "RunPCA" function from the "Seurat" package to reduce the dimensionality of the scRNA-seq data based on these highly variable genes. To integrate the two samples, we used the canonical correlation analysis (CCA) approach. ElbowPlot analysis was used to identify significant principal components, and cell clustering was done with the top ten principal components using PCA. The "FindNeighbors" function was used to create a k-nearest neighbor graph in the PCA space using the Euclidean distance, and the Uniform Manifold Approximation and Projection (UMAP) approach was used to reduce data dimensionality and visualize it. Cluster annotation was done with reference data from a previously published study[ 19 ]. To evaluate the gene expression scores of mitophagy-related genes (MRGs) for each cell in the scRNA-seq dataset, we used the "AddModuleScore" function to generate individual cell scores. In addition, we performed Gene Set Variation Analysis (GSVA) to examine the enrichment of signature pathways for each cell type using the Molecular Signatures Database (MSigDB)[ 20 , 21 ]. Finally, we conducted Gene Set Enrichment Analysis (GSEA) on the differentially expressed genes between the high and low groups based on mitophagy scores. 2.4 Construction and validation of a prediction model Machine learning models have developed as effective methods for dealing with enormous datasets and aiding rapid analysis and decision-making in disease diagnosis. In this study, we used a combinatorial technique leveraging support vector machines (SVM)[ 22 ] and random forest (RF) [ 23 ] algorithms to pick significant features from a list of 29 mitophagy-associated genes. Following that, a prediction model was created using logistic regression[ 24 ]. The prediction model's performance was examined using receiver operating characteristic (ROC) curve analysis to determine its sensitivity and specificity. The area under the curve (AUC) of the ROC curve was used to assess the model's prediction performance. To display the relative importance of the feature genes within the model and their contribution to the prediction results, a nomogram model was created using "rms" package (version = 6.8-2). To further validate the model's prediction performance, calibration curves were generated to determine the agreement between anticipated probabilities and observed outcomes. In addition, decision curve analysis (DCA), which involves calculating the net benefit of a predictive model at different decision thresholds while taking into account the false positives, false negatives, and costs associated with different decision thresholds, was used to assess the model's clinical value by estimating the net benefits at various probability thresholds[ 25 ]. The net benefit is calculated using the formula: Where: N is the total number of patients; Ture positives and False Positives are derived from the model’s predictions at a given threshold probability Pt. 2.5 Analysis based on the Human Protein Atlas database To further validate the prediction model, we utilized the Human Protein Atlas (HPA) database ( https://www.proteinatlas.org/ )[ 26 , 27 ] to investigate whether the expression levels of the feature genes differ between breast cancer tissues and normal tissues at the protein level. These resources provide comprehensive, publicly available data on protein expression across various human tissues and cancer types. 2.6 Correlation between feature genes and immunity To investigate the potential relationship between the feature genes and the immune microenvironment, we used the "CIBERSORT" package (version = 1.03) to estimate the level of immune cell infiltration in both normal and tumor samples. This computer method enables the deconvolution of complex gene expression data and the measurement of different immune cell subpopulations within tissue samples. In addition to looking at immune cell infiltration patterns, we investigated the links between feature genes, immune infiltration, and important inflammatory variables. The findings of these investigations were visually depicted using the "ggplot2" software. 2.7 Hierarchical Clustering To investigate the potential subgroups within breast cancer (BC) samples, we employed hierarchical clustering analysis using the "hclust" (‘stats’ R package, version = 4.4.0) R function on a cohort of 783 BC patients. The analysis was based on the expression profiles of nine feature genes, which were selected through a rigorous feature selection process. The parametric h was set to 4.5 to ensure optimal clustering performance. Following the identification of distinct BC subtypes, we explored the relationship between the expression levels of the nine feature genes, hereafter referred to as the molecular risk genes (MRGs), and the clinical pathological characteristics within each subtype. To quantify the relative expression of MRGs in each sample, we utilized the gene set variation analysis (GSVA) package (version = 1.52.3) to calculate sample-specific scores. The Kruskal test was then applied to assess the statistical significance of the differential expression of MRGs among the three identified BC subtypes. To evaluate the prognostic value of the identified subtypes, Kaplan-Meier survival analysis was conducted using the "survival"(version = 3.5-8) and "survminer"(version = 0.4.9) R packages. This analysis aimed to determine the impact of the BC subtypes on patient outcomes, providing insights into the potential clinical utility of the identified subgroups. 2.8 Analysis of immune infiltration in BC subtypes We utilized the CIBERSORT program in R to calculate the proportions of 22 immune cell subtypes among the identified BC subtypes. The LM22 gene expression matrix, which includes 547 genes identifying 22 human immune cell types, used as a reference for immune cell characterization. The CIBERSORT algorithm generated immune cell infiltration profiles for each BC subtype. In addition, we evaluated the expression levels of immune and stromal cells in the samples. The distribution of immune and stromal cells among BC subtypes was compared and displayed using box plots. 2.9 Weighted Gene Co-expression Network Analysis To identify hub genes within each breast cancer (BC) subtype, we performed consensus module analysis using the "WGCNA" package in R[ 28 ](version = 1.73). First, genes were ranked based on their median absolute deviation (MAD) values, and the top 5,000 genes were selected for further analysis. Subsequently, samples and genes with excessive missing values were filtered out to ensure data quality. A similarity matrix was constructed using the selected genes, which was then transformed into an adjacency matrix using a soft-thresholding power (β) of 5. This transformation emphasizes strong correlations and penalizes weak correlations, resulting in a scale-free network topology. Next, the “blockwiseModules” function was applied to perform modular analysis on the adjacency matrix. This function employs a dynamic tree-cutting algorithm to detect modules (clusters) of highly interconnected genes. The resulting gene co-expression network assigns genes to different modules based on their co-expression patterns. For each identified module, module eigengenes (MEs) were calculated as the first principal component, representing the overall expression profile of the module. To investigate the relationship between the identified modules and the phenotype data, Pearson's correlation coefficients were calculated between the MEs and the phenotype data. The correlation analysis reveals the association of each module with specific phenotypic traits, such as clinical features or outcomes. Finally, the correlation between the module genes and the phenotype data was visualized using a heat map, providing a comprehensive overview of the module-trait relationships. 2.10 Analysis of enrichment To obtain insight into the biological importance of the discovered gene modules, we chose genes from the module that had the strongest connection with each BC subtype based on WGCNA results. These gene sets were then subjected to enrichment analysis utilizing the Metascape database ( https://metascape.org/ )[ 29 ]. In addition, we used the "msigdbr" package in R to import gene sets from the biological processes (BP) category, and the "fgsea" package to perform biological process enrichment analysis for different BC subtypes. 2.11 Immune Checkpoint Gene Expression and Response Prediction We evaluated the expression of immune checkpoint genes across the identified BC subtypes and used the Tumor Immune Dysfunction and Exclusion (TIDE) database ( http://tide.dfci.harvard.edu/ )[ 30 , 31 ] to predict the response to immune checkpoint blockade (ICB) therapy. TIDE uses gene expression fingerprints from tumors and immune cells to assess two major strategies of tumor immune evasion: malfunction and cytotoxic T lymphocyte (CTL) exclusion. By taking into account these parameters, TIDE gives a thorough assessment of the probable efficacy of ICB therapy in each BC subtype, allowing for more accurate immunotherapy response prediction. 2.12 Drug sensitivity analysis Half-maximal inhibitory concentration (IC50) is a widely used concentration measure for evaluating the therapeutic efficacy or toxicity of drugs. To assess chemotherapy sensitivity in various BC subtypes, we utilized the R package "oncoPredict" in conjunction with expression matrices from the Genomics of Drug Sensitivity in Cancer (GDSC) database and drug treatment data as a training set[ 32 ]. 2.13 Cell culture, RNA extraction and quantitative Real-Time PCR (qPCR) Fujian Cancer Hospital supplied the human breast epithelial cell line MCF-10A, as well as breast cancer cell lines (MCF-7, MDA-MB-231, MDA-MB-468, Hs-578T, HCC1937, and BT-549). MCF-10A cells were cultured in DMEM/F12 medium (Gibco) with 5% horse serum (Gibco), 20 ng/mL epidermal growth factor (EGF, PeproTech), 0.5 µg/mL hydrocortisone (Sigma-Aldrich), 100 ng/mL cholera toxin (Sigma-Aldrich), and 10 µg/mL insulin (Sigma-Aldrich). Breast cancer cell lines (MCF-7, MDA-MB-231, MDA-MB-468, Hs-578T, HCC1937, and BT-549) were grown in Dulbecco's Modified Eagle Medium (DMEM, Gibco) with 10% fetal bovine serum (FBS, Gibco) and 1% penicillin-streptomycin. The cells were cultivated at 37°C in a humid atmosphere with 5% CO₂. Cell viability and confluence were assessed daily, and subcultures were carried out at 80–90% confluence using 0.25% trypsin-EDTA (Gibco). Total RNA was extracted from cells using the TRIzol reagent (Invitrogen) according to the manufacturer's instructions. A NanoDrop spectrophotometer (Thermo Fisher Scientific) was used to determine RNA content and purity, with A260/A280 ratios of 1.8 to 2.0 considered acceptable. cDNA was synthesized from 1 µg of total RNA using the PrimeScript RT Reagent Kit (Takara Bio) and oligo(dT) primers. SYBR Green Master Mix (Applied Biosystems) was utilized for quantitative PCR with a StepOnePlus Real-Time PCR System (Applied Biosystems). The PCR conditions were as follows: Initial denaturation at 95°C for 10 minutes, then 40 cycles at 95°C for 15 seconds and 60°C for 1 minute. Gene expression levels were normalized to GAPDH, and relative quantification was performed using the 2^(-ΔΔCt) approach. The primer sequences used in this study are shown in Supplementary file Table S1 . 2.14 Statistical analysis All statistical analyses were performed using the R program (version 4.3.2) or GraphPad Prism (version 9.5.0). The Wilcox rank sum test was used to identify significant differences between two independent groups, while the Kruskal test was used to identify significant differences between three groups. Results from qPCR were accessed using a one-way ANOVA. P-values < 0.05 were considered statistically significant. Results 3.1 Breast cancer (BC) is associated with upregulated Mitophagy-Related Genes (MRGs) and immune system activation. Using The Cancer Genome Atlas (TCGA) training dataset, we examined the 'gsva' scores of MRGs to evaluate the expression profiles of mitophagy genes in both tumor and normal tissues. Our findings showed that, in comparison to tumor tissues, normal tissues had substantially reduced MRG expression levels (Figure 1A). Seven of the 28 MRGs (CSNK2B, FUNDC1, MTERF3, PGAM5, SRC, ULK1, and TOMM40) were significantly elevated in tumor tissues, whereas eighteen were significantly downregulated. The remaining genes did not exhibit significant differences between normal and tumor tissues (Figure 1B-C, Table 1, supplement Fig.1A). Furthermore, correlation analysis was used to further understand the links between these MRGs (supplement Fig.1B). To examine the immune landscape in BC, we used the 'CIBERSORT' method to calculate the infiltration of different immune cell types in tumor and normal tissues. Our findings demonstrated a considerable increase in the presence of activated T cells, helper T cells, regulatory T cells (Tregs), M0 macrophages, and M1 macrophages in tumor tissues. Tumors had considerably lower levels of naïve B cells, plasma cells, resting T cells, active NK cells, monocytes, and M2 macrophages than normal tissues (Figure 1D). 3.2 Single-Cell RNA Sequencing Analysis Reveals the Role of MRGs in Different Cell Types We examined the single-cell RNA sequencing (scRNA-seq) data from two breast cancer tumor samples (GSM7910990 and GSM7910991) acquired from the GSE248288 dataset in order to look into the expression and functional implications of mitophagy-related genes (MRGs) at the single-cell level. Quality control techniques were used to guarantee that high-quality cells were included in the study (Figure 2A, Supplementary Figures 2A-C). Principal component analysis (PCA) was used to minimize dimensionality by selecting the top 2,000 highly variable genes (Supplementary Figure 2D). The cells were then clustered into 17 separate clusters (Supplementary Figure 2E). Based on biomarkers from the prior study[19], we classified the cells into major cell types, including T cells, natural killer (NK) cells, epithelial cells, endothelial cells, myeloid cells, plasma cells, B cells, and fibroblasts (Figures 2B-C). To analyze MRG expression in various cell types, we used the "AddModuleScore" tool to compute the MRG score for each cell type. Remarkably, NK cells trailed B cells in having the highest MRG score (Figure 2D). Gene set variation analysis (GSVA) was carried out to look more closely at the biological mechanisms connected to MRGs in various cell types. The findings showed that B cells had an increased allograft rejection pathway, fibroblasts primarily had an increased epithelial mesenchymal transition, and NK cells had increased PI3K-AKT-mTOR signaling, cholesterol homeostasis, and Myc-targets-V1 signaling (Supplementary Figure 3A). Gen set enrichment analysis (GSEA) was performed to obtain insights into the differential gene expression between samples with high and low mitophagy-related scores (Hscore and Lscore) (Figure 2E). Pathways such the ribosome, proteasome, MAPK signaling pathway, and taurine and hypotaurine metabolism had a substantial enrichment of the differentially expressed genes across the high and low groups (Figure 2F, Supplementary Figure 3B). These results suggest that MRGs may influence protein synthesis, degradation, signaling cascades, and metabolic processes in breast cancer cells. 3.3 Construction and Validation of a Predictive Model for Breast Cancer Occurrence To identify key mitophagy-related genes (MRGs) associated with breast cancer (BC) occurrence and develop a predictive model, we employed random forest (RF) algorithms and support vector machine (SVM). The SVM-recursive feature elimination (SVM-RFE) analysis was performed on the 28 MRGs, the top 10 ideal feature genes were chosen as a result. In parallel, a random forest analysis was carried out, and gene importance scores from the model were used to pick the top 10 genes based on the error rate curve (Figures 3A, B). Nine genes in total were obtained by intersecting the genes chosen by SVM-RFE and random forest: PGAM5, PRKN, TOMM40, FUNDC1, MAP1LC3B, PINK1, MTERF3, CSNK2A2, and SRC. qPCR analysis revealed distinct expression patterns of these genes across the panel of breast cancer cell lines or (MCF-7, MDA-MB-231, MDA-MB-468, Hs-578T, HCC1937, and BT-549) compared to the non-tumorigenic MCF-10A control (Supplement file Figure 7) . Gene expression levels of PGAM5, PRKN, TOMM40, FUNDC1, MAP1LC3B, PINK1, MTERF3, CSNK2A2, and SRC were assessed across six breast cancer cell lines (MCF-7, MDA-MB-231, MDA-MB-468, Hs-578T, HCC1937, and BT-549) and compared to the normal mammary epithelial cell line MCF-10A. We found that PGAM5, MAP1LC3B, PRKN, PINK1, SRC, FUNDC1, CSNK2A2, and MTERF3 showed significantly higher expression in most of the breast cancer cell lines compared to MCF-10A, particularly in MDA-MB-468, BT-549, and HCC1937 cells. TOMM40 expression was stable across both cancer and normal cell lines, showing no significant difference from MCF-10A. To assess the association between the selected genes and BC occurrence, a multivariable logistic regression analysis was performed, confirming the significance of these nine genes (Figure 3C). An area under the curve (AUC) of 0.99 was shown by the receiver operating characteristic (ROC) curve, which was used to assess the prediction ability of the model. An independent dataset used for internal validation of the model produced an AUC of 0.988(Figures 3D, E). Calibration plots and decision curve analysis (DCA) were employed to assess the calibration and clinical utility of the predictive model (Figures 3F-G). A nomogram model was created to offer a user-friendly tool for displaying the feature genes' predictive performance (Figure 3H). When all nine feature genes were combined into one model, the predictive value of the integrated model was higher than that of the feature genes alone. Then, we gathered the IHC staining of nine proteins linked to feature genes from the HPA database. These proteins were taken from both normal and breast cancer tissue. Our findings are supported by the observation that three of the feature genes (SRC, TOMM40, and PGAM5) had significantly greater levels of protein expression in tumor samples compared to normal samples (Supplementary Figure 4). There was no discernible variation in the levels of protein expression of the other feature genes between the normal and breast cancer samples. 3.4 Identification and Differential Analysis of Breast Cancer Subtypes All tumor samples were categorized into three different subtypes using hierarchical clustering based on the nine mitophagy-related characteristic genes expression matirx (Figure 4A). Significant variations in the expression levels of these characteristic genes were found among the three groups based on differential expression analysis (Figure 4B). The heatmap (Figure 4C) shows the distribution of clinical pathological features across the subtypes that have been discovered, as well as the feature genes' differential expression patterns. A noteworthy finding was that the predictive model performed well in differentiating between the three subtypes of breast cancer (Figure 4D). Gene set variation analysis (GSVA) was used to better describe the molecular variations between the subtypes. The findings showed that, in comparison to the other two clusters, Cluster 1 had a substantially higher overall expression level of mitophagy-related genes (MRGs) (Figure 4E). To evaluate the predictive value of the discovered subgroups, survival analysis was performed. Nevertheless, there were no appreciable variations in survival rates between the three groups (Figure 4F). This data suggests that patient survival outcomes may not be directly impacted by the subgroups that were discovered, despite their different molecular profiles. 3.5 Constructing a Gene Co-expression Network We used a weighted gene co-expression network analysis (WGCNA) to find co-expression modules that correlated best with each subtype of breast cancer. The analysis encompassed all samples, and 5 was determined to be the ideal soft-thresholding value. Thirteen unique co-expression modules were produced as a result of the WGCNA (Figures 5A, B, C). Interestingly, the blue module had the highest association with Cluster 1, while the turquoise module had a correlation with Cluster 2. Furthermore, the black module was discovered to be connected with Cluster 3 (Figure 5D). Using the Metascape database, functional enrichment analysis was performed on the modules that showed the strongest correlation with each cluster in order to obtain insights into the biological significance of the co-expression modules. In Cluster 1, the genes in the blue module were considerably enriched in pathways related to cell cycle biological processes (GO:000278, GO:0010564, hsa04110) (Supplementary Figure 5, top). This data implies that disrupted cell cycle regulation, a hallmark of cancer, may be present in Cluster 1. Cluster 2 showed a high enrichment of genes within the turquoise module in pathways associated with the synthesis of cilia and axonemal dynein complexes (Supplementary Figure 5, middle). Crucial organelles, cilia are engaged in a number of biological functions, including as movement and cell signaling. The predominance of cilium-related pathways in Cluster 2 implies potential abnormalities in ciliary function, which may contribute to the unique biological properties of this subtype. Lastly, In Cluster 3, the black module's genes were considerably enriched in the regulation of small GTPase-mediated signal transduction (Supplementary Figure 5, bottom). A number of essential cellular functions, including cell division, migration, and survival, are regulated by small GTPases like Ras and Rho. The pathway's enrichment in Cluster 3 raises the possibility that faulty small GTPase signaling contributes to the pathophysiology of this particular subtype of breast cancer. 3.6 Estimating the Immune Response and Drug Sensitivity to Chemotherapy To study differences in the immune microenvironment between the identified breast cancer (BC) subtypes, we used the CIBERSORT program to estimate the infiltration proportions of various immune and stromal cells. The research found substantial differences in infiltration levels of naïve B cells, resting T cells, activated T cells, regulatory T cells, helper T cells, monocytes, macrophages, resting dendritic cells, and resting mast cells amongst clusters. (Supplementary Figure 6A). These findings suggest that each BC subtype may exhibit distinct immune cell compositions, potentially influencing their response to immunotherapy. Tumor Immune Dysfunction and Exclusion (TIDE) algorithm was utilized to forecast the efficacy of immune checkpoint inhibition treatment. Interestingly, compared to other clusters, cluster 2 showed the lowest TIDE score, dysfunction and exclusion, but the highest MSI score, suggesting a possible improved response to immunotherapy (Figures 6A-D). To determine the sensitivity of various BC subtypes to chemotherapeutic medications, we used the “oncoPredict” program to compute half-maximal inhibitory concentration (IC50) values for many routinely used medicines. The findings showed that, in comparison to the other two clusters, cluster 1 had lower IC50 values for methotrexate, paclitaxel, gemcitabine, gemcitabine with carboplatin, and doxorubicin (Figures 6E-J). This finding suggests that patients in cluster 1 might benefit more from these specific chemotherapy regimens but not benefit from cyclophosphamide than other clusters. Conversely, cluster 2 exhibited significantly higher IC50 values in doxorubicin, paclitaxel, gemcitabine, gemcitabine plus carboplatin, and methotrexate, indicating that patients in cluster 2 might not benefits from these regents than other cluster patients. Given the effectiveness of immune checkpoint inhibitors in a subset of BC patients, we looked into the differences in immunotherapy checkpoint gene expression among the identified BC subtypes. Significant variations in the expression levels of several important checkpoint genes, such as BTNL9, BTN2A1, TNFSF9, PDCD1, CD47, PVR, TNFRSF4, CD274, and CTLA4, were found by our research (Supplementary Figure 6B). These findings suggest that the BC subtypes may exhibit varying responses to immune checkpoint blockade therapies, highlighting the potential for personalized immunotherapy approaches based on the molecular characteristics of each subtype. Discussion The intricate relationship between mitochondrial dynamics and cancer progression has emerged as a pivotal area of research, particularly in breast cancer (BC), where metabolic reprogramming and immune evasion are critical determinants of tumor aggressiveness. Our comprehensive investigation into mitophagy-related genes (MRGs) reveals their profound influence on BC heterogeneity, offering novel insights into diagnostic stratification, immune modulation, and therapeutic vulnerabilities. By integrating bulk transcriptomics, single-cell resolution analyses, and machine learning-driven modeling, we delineate a molecular framework that positions MRGs as central regulators of BC biology. The dysregulation of MRGs in BC tissues underscores their dual role in tumor survival and metabolic adaptation. Our data demonstrated significant upregulation of seven MRGs, including FUNDC1, SRC, and TOMM40, alongside the suppression of 18 others such as PINK1 and PRKN. FUNDC1 plays a critical role in tumor biology through its regulation of mitochondrial dynamics, particularly in modulating pathways related to cancer cell survival and proliferation. Notably, elevated FUNDC1 expression has been closely associated with poor clinical outcomes and chemoresistance in endometrial carcinoma.[33, 34] Conversely, the downregulation of PINK1-PRKN-mediated mitophagy, a pathway essential for eliminating dysfunctional mitochondria, may lead to ROS accumulation and genomic instability, hallmarks of aggressive BC phenotypes[35-37]. The elevated expression of SRC and TOMM40 further highlights their dual functionality: SRC kinase not only drives oncogenic signaling but also modulates mitochondrial fission-fusion dynamics[38-41], while TOMM40, a mitochondrial translocase, supports the import of nuclear-encoded proteins necessary for oxidative phosphorylation[42]. These findings align with prior studies linking mitochondrial dysfunction to chemoresistance and metastasis, yet our work uniquely identifies MRGs as integrative biomarkers that bridge metabolic and genomic instability in BC. A striking observation was the discordance between mRNA and protein levels for several MRGs, such as MAP1LC3B and ATG5. This discrepancy likely reflects post-transcriptional regulatory mechanisms, including miRNA-mediated repression or proteasomal degradation. For instance, miR-30a, known to target MAP1LC3B in hepatocellular carcinoma and diabetic cataract[43, 44], may similarly suppress its translation in BC[45], underscoring the necessity of multi-omics approaches to capture the full spectrum of gene regulation. Such insights emphasize that transcriptomic signatures alone may insufficiently predict functional protein activity, necessitating complementary proteomic and spatial analyses in future studies. The tumor microenvironment (TME) of BC exhibited a complex immune landscape characterized by infiltrating activated T cells, macrophages, and regulatory T cells (Tregs)[46-49]. This immune profile correlated with MRG expression patterns, suggesting a bidirectional interplay between mitophagy and immune evasion. Mitochondrial-derived vesicles (MDVs), which transport mitochondrial antigens to MHC-I molecules, could theoretically activate CD8+ T cells[50]; however, excessive mitophagy in tumor cells may deplete these antigens, thereby evading immune recognition—a phenomenon previously documented in melanoma. At the single-cell level, natural killer (NK) cells displayed high MRG scores coupled with enriched PI3K-AKT-mTOR signaling. Given that mTOR signaling governs NK cell metabolism and cytotoxicity, sustained mitophagy in these cells might paradoxically impair their anti-tumor function through metabolic exhaustion. This hypothesis is supported by the low TIDE dysfunction scores in Cluster 2, which exhibited preserved NK cell activity and heightened responsiveness to immune checkpoint inhibitors (ICIs). These findings position mitophagy as a double-edged sword: while promoting tumor survival, it may also sculpt immune cell fitness, presenting opportunities for therapeutic intervention. The development of a nine-MRG diagnostic model achieving an AUC of 0.99 represents a paradigm shift in BC stratification. Traditional biomarkers such as ER, PR, Her2, and Ki67, while invaluable for guiding endocrine and targeted therapies, fail to capture the metabolic heterogeneity underlying therapeutic resistance[51-54]. For instance, MRG-high tumors (Cluster 1) demonstrated heightened sensitivity to paclitaxel and gemcitabine, likely due to their reliance on mitochondrial metabolism—a feature invisible to conventional IHC profiling. The model’s strength lies in its integration of multi-gene interactions, which capture pathway-level dysregulation undetectable by single-marker assays. This approach mirrors the success of multi-gene panels like Oncotype DX, which outperform individual biomarkers in predicting chemotherapy benefit. However, transitioning from transcriptomic data to clinical utility requires addressing challenges such as platform-specific batch effects and the integration of MRG signatures with established molecular classifiers like PAM50. Preliminary analyses suggest overlap between MRG-high clusters and basal-like subtypes, known for their metabolic dependencies, highlighting the potential for unified molecular taxonomies that enhance prognostic accuracy. Hierarchical clustering identified three BC subtypes with distinct molecular and functional profiles. Cluster 1, characterized by elevated MRG expression and cell cycle activation, mirrors the “basal-like” subtype, renowned for genomic instability and chemosensitivity. The low IC50 values for taxanes and antimetabolites in this cluster align with their mechanism of targeting rapidly dividing cells. Cluster 2, enriched in cilium assembly pathways and exhibiting low TIDE scores, may represent an immune-inflamed phenotype responsive to ICIs. The reappearance of cilia in carcinomas, often linked to Hedgehog pathway activation, could recruit immunosuppressive macrophages, yet the high PD-L1 expression and microsatellite instability (MSI) in this cluster suggest a unique susceptibility to immunotherapy. Cluster 3, dominated by dysregulated small GTPase signaling, likely drives metastatic dissemination through cytoskeletal remodeling and invadopodia formation. Therapeutic targeting of Rho/ROCK inhibitors in this subset could disrupt invasion-proliferation crosstalk, offering a novel strategy to mitigate metastasis. Despite these advancements, our study has limitations. The retrospective nature of TCGA data introduces selection bias, necessitating prospective validation in cohorts with annotated treatment responses. Functional experiments, such as CRISPR-based gene editing or patient-derived organoid models, are essential to establish causality between MRGs and observed phenotypes. Additionally, the limited scope of single-cell RNA sequencing—restricted to two tumor samples—calls for expansion to diverse BC subtypes and metastatic lesions. Spatial multi-omics technologies could further elucidate the spatial distribution of MRGs within tumor niches, revealing interactions between cancer cells and immune infiltrates. Conclusions Our study repositions mitophagy as a central axis in BC biology, influencing metabolic plasticity, immune evasion, and therapeutic response. The MRG-based diagnostic model and subtype classification system transcend conventional paradigms, offering a roadmap for personalized therapy. Future research must prioritize translational validation, integrating these insights into clinical trials to assess their impact on patient outcomes. By aligning mitophagy inhibition with chemotherapy or immunotherapy, clinicians could exploit metabolic vulnerabilities while counteracting immune suppression, ultimately advancing precision oncology in breast cancer. Abbreviations AUC Area under the curve BC Breast cancer CCA Canonical correlation analysis CTL Cytotoxic T lymphocyte DCA Decision curve analysis DEGs Differentially expressed genes GDSC Genomics of Drug Sensitivity in Cancer GEO Gene Expression Omnibus database GSEA Gene Set Enrichment Analysis GSVA Gene Set Variation Analysis HPA Human Protein Atlas IC50 Half-maximal inhibitory concentration ICB Immune checkpoint blockade MAD Median absolute deviation Mes Module eigengenes MRGs Mitophagy-related genes MSI Microsatellite instability mtDNA Mitochondrial DNA OXPHOS Oxidative phosphorylation PCA Principal component analysis qPRC Quantitative Real-Time PCR ROC Receiver operating characteristic scRNA-seq Single-cell RNA sequencing SVM Support vector machines SVM-RFE SVM-recursive feature elimination TCGA The Cancer Genome Atlas TIDE Tumor Immune Dysfunction and Exclusion UMAP Uniform Manifold Approximation and Projection WGCNA Weighted gene co-expression network analysis Declarations Ethics approval and consent to participate ： Not applicable Consent for publication ： Not applicable Availability of data and materials ： Not applicable Competing interests ：The authors declare that they have no competing interests. Funding ： Not applicable Authors' contributions ： Wenbin Guo conceived, analyzed, and wrote the original manuscript; Zhiqiang Ye carried out the qPCR and cell culture; and Dan Wu re-edited it. All authors read and approved the final manuscript. Acknowledgments ： Not applicable. References Chen K, Lu P, Beeraka NM, Sukocheva OA, Madhunapantula SV, Liu J, Sinelnikov MY, Nikolenko VN, Bulygin KV, Mikhaleva LM et al (2022) Mitochondrial mutations and mitoepigenetics: Focus on regulation of oxidative stress-induced responses in breast cancers. Semin Cancer Biol 83:556–569 Gautam S, Marwaha D, Singh N, Rai N, Sharma M, Tiwari P, Urandur S, Shukla RP, Banala VT, Mishra PR (2023) Self-Assembled Redox-Sensitive Polymeric Nanostructures Facilitate the Intracellular Delivery of Paclitaxel for Improved Breast Cancer Therapy. Mol Pharm 20(4):1914–1932 Yu H, Yang C, Jian L, Guo S, Chen R, Li K, Qu F, Tao K, Fu Y, Luo F et al (2019) Sulfasalazine–induced ferroptosis in breast cancer cells is reduced by the inhibitory effect of estrogen receptor on the transferrin receptor. Oncol Rep 42(2):826–838 Cui L, Gouw AM, LaGory EL, Guo S, Attarwala N, Tang Y, Qi J, Chen YS, Gao Z, Casey KM et al (2021) Mitochondrial copper depletion suppresses triple-negative breast cancer in mice. Nat Biotechnol 39(3):357–367 Ramchandani D, Berisa M, Tavarez DA, Li Z, Miele M, Bai Y, Lee SB, Ban Y, Dephoure N, Hendrickson RC et al (2021) Copper depletion modulates mitochondrial oxidative phosphorylation to impair triple negative breast cancer metastasis. Nat Commun 12(1):7311 Sun Y, Shen W, Hu S, Lyu Q, Wang Q, Wei T, Zhu W, Zhang J (2023) METTL3 promotes chemoresistance in small cell lung cancer by inducing mitophagy. J Exp Clin Cancer Res 42(1):65 Wang J, Zhu P, Li R, Ren J, Zhou H (2020) Fundc1-dependent mitophagy is obligatory to ischemic preconditioning-conferred renoprotection in ischemic AKI via suppression of Drp1-mediated mitochondrial fission. Redox Biol 30:101415 Hu Z, Ding J, Ma ZC, Sun RP, Seoane JA, Shaffer JS, Suarez CJ, Berghoff AS, Cremolini C, Falcone A et al (2019) Quantitative evidence for early metastatic seeding in colorectal cancer. Nat Genet 51(7):1113– Miao Z, Miao Z, Wang S, Wu H, Xu S (2022) Exposure to imidacloprid induce oxidative stress, mitochondrial dysfunction, inflammation, apoptosis and mitophagy via NF-kappaB/JNK pathway in grass carp hepatocytes. Fish Shellfish Immunol 120:674–685 Li R, Xin T, Li D, Wang C, Zhu H, Zhou H (2018) Therapeutic effect of Sirtuin 3 on ameliorating nonalcoholic fatty liver disease: The role of the ERK-CREB pathway and Bnip3-mediated mitophagy. Redox Biol 18:229–243 Panigrahi DP, Praharaj PP, Bhol CS, Mahapatra KK, Patra S, Behera BP, Mishra SR, Bhutia SK (2020) The emerging, multifaceted role of mitophagy in cancer and cancer therapeutics. Semin Cancer Biol 66:45–58 Chen C, Liu Y, Liu L, Si CH, Xu YX, Wu XK, Wang CZ, Sun ZQ, Kang QZ (2023) Exosomal circTUBGCP4 promotes vascular endothelial cell tipping and colorectal cancer metastasis by activating Akt signaling pathway. J EXPERIMENTAL Clin CANCER Res 42(1) Lin Q, Li S, Jin H, Cai H, Zhu X, Yang Y, Wu J, Qi C, Shao X, Li J et al (2023) Mitophagy alleviates cisplatin-induced renal tubular epithelial cell ferroptosis through ROS/HO-1/GPX4 axis. Int J Biol Sci 19(4):1192–1210 Altea-Manzano P, Doglioni G, Liu YW, Cuadros AM, Nolan E, Fernández-García J, Wu Q, Planque M, Laue KJ, Cidre-Aranaz F et al (2023) A palmitate-rich metastatic niche enables metastasis growth via p65 acetylation resulting in pro-metastatic NF-κB signaling. Nat CANCER 4(3):344– Lu Y, Li Z, Zhang S, Zhang T, Liu Y, Zhang L (2023) Cellular mitophagy: Mechanism, roles in diseases and small molecule pharmacological regulation. Theranostics 13(2):736–766 Tan J, Egelston CA, Guo W, Stark JM, Lee PP (2024) STING signalling compensates for low tumour mutation burden to drive anti-tumour immunity. EBioMedicine 101:105035 Hao Y, Stuart T, Kowalski MH, Choudhary S, Hoffman P, Hartman A, Srivastava A, Molla G, Madad S, Fernandez-Granda C et al (2024) Dictionary learning for integrative, multimodal and scalable single-cell analysis. Nat Biotechnol 42(2):293–304 Hao Y, Hao S, Andersen-Nissen E, Mauck WM 3rd, Zheng S, Butler A, Lee MJ, Wilk AJ, Darby C, Zager M et al (2021) Integrated analysis of multimodal single-cell data. Cell 184(13):3573–3587 e3529 Liu YM, Ge JY, Chen YF, Liu T, Chen L, Liu CC, Ma D, Chen YY, Cai YW, Xu YY et al (2023) Combined Single-Cell and Spatial Transcriptomics Reveal the Metabolic Evolvement of Breast Cancer during Early Dissemination. Adv Sci (Weinh) 10(6):e2205395 Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, Paulovich A, Pomeroy SL, Golub TR, Lander ES et al (2005) Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A 102(43):15545–15550 Liberzon A, Birger C, Thorvaldsdottir H, Ghandi M, Mesirov JP, Tamayo P (2015) The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst 1(6):417–425 CORTES C (2005) Support-Vector Networks. Mach Learn 20:273–297 BREIMAN L: Random Forests. Mach Learn, (2001) 45:5–32 Marchevsky AM, Patel S, Wiley KJ, Stephenson MA, Gondo M, Brown RW, Yi ES, Benedict WF, Anton RC, Cagle PT (1998) Artificial neural networks and logistic regression as tools for prediction of survival in patients with Stages I and II non-small cell lung cancer. Mod Pathol 11(7):618–625 Van Calster B, Wynants L, Verbeek JFM, Verbakel JY, Christodoulou E, Vickers AJ, Roobol MJ, Steyerberg EW (2018) Reporting and Interpreting Decision Curve Analysis: A Guide for Investigators. Eur Urol 74(6):796–804 Uhlen M, Fagerberg L, Hallstrom BM, Lindskog C, Oksvold P, Mardinoglu A, Sivertsson A, Kampf C, Sjostedt E, Asplund A et al (2015) Proteomics. Tissue-based map of the human proteome. Science 347(6220):1260419 Sjostedt E, Zhong W, Fagerberg L, Karlsson M, Mitsios N, Adori C, Oksvold P, Edfors F, Limiszewska A, Hikmet F et al (2020) An atlas of the protein-coding genes in the human, pig, and mouse brain. Science 367(6482) Langfelder P, Horvath S (2008) WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics 9:559 Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, Benner C, Chanda SK (2019) Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun 10(1):1523 Fu J, Li K, Zhang W, Wan C, Zhang J, Jiang P, Liu XS (2020) Large-scale public data reuse to model immunotherapy response and resistance. Genome Med 12(1):21 Jiang P, Gu S, Pan D, Fu J, Sahu A, Hu X, Li Z, Traugh N, Bu X, Li B et al (2018) Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med 24(10):1550–1558 Maeser D, Gruener RF, Huang RS (2021) oncoPredict: an R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data. Brief Bioinform 22(6) Li J, Agarwal E, Bertolini I, Seo JH, Caino MC, Ghosh JC, Kossenkov AV, Liu Q, Tang HY, Goldman AR et al (2020) The mitophagy effector FUNDC1 controls mitochondrial reprogramming and cellular plasticity in cancer cells. Sci Signal 13(642) Hou H, Er P, Cheng J, Chen X, Ding X, Wang Y, Chen X, Yuan Z, Pang Q, Wang P et al (2017) High expression of FUNDC1 predicts poor prognostic outcomes and is a promising target to improve chemoradiotherapy effects in patients with cervical cancer. Cancer Med 6(8):1871–1881 Xie Y, Liu J, Kang R, Tang D (2020) Mitophagy Receptors in Tumor Biology. Front Cell Dev Biol 8:594203 Oshi M, Gandhi S, Yan L, Tokumaru Y, Wu R, Yamada A, Matsuyama R, Endo I, Takabe K (2022) Abundance of reactive oxygen species (ROS) is associated with tumor aggressiveness, immune response, and worse survival in breast cancer. Breast Cancer Res Treat 194(2):231–241 Guo M, Wang SM (2021) Genome Instability-Derived Genes Are Novel Prognostic Biomarkers for Triple-Negative Breast Cancer. Front Cell Dev Biol 9:701073 Djeungoue-Petga MA, Lurette O, Jean S, Hamel-Cote G, Martin-Jimenez R, Bou M, Cannich A, Roy P, Hebert-Chatelain E (2019) Intramitochondrial Src kinase links mitochondrial dysfunctions and aggressiveness of breast cancer cells. Cell Death Dis 10(12):940 Sirvent A, Benistant C, Roche S (2012) Oncogenic signaling by tyrosine kinases of the SRC family in advanced colorectal cancer. Am J Cancer Res 2(4):357–371 Pelaz SG, Tabernero A (2022) Src: coordinating metabolism in cancer. Oncogene 41(45):4917–4928 Lu KV, Zhu S, Cvrljevic A, Huang TT, Sarkaria S, Ahkavan D, Dang J, Dinca EB, Plaisier SB, Oderberg I et al (2009) Fyn and SRC are effectors of oncogenic epidermal growth factor receptor signaling in glioblastoma patients. Cancer Res 69(17):6889–6898 Yang W, Shin HY, Cho H, Chung JY, Lee EJ, Kim JH, Kang ES (2020) TOM40 Inhibits Ovarian Cancer Cell Growth by Modulating Mitochondrial Function Including Intracellular ATP and ROS Levels. Cancers (Basel) 12(5) Fu XT, Shi YH, Zhou J, Peng YF, Liu WR, Shi GM, Gao Q, Wang XY, Song K, Fan J et al (2018) MicroRNA-30a suppresses autophagy-mediated anoikis resistance and metastasis in hepatocellular carcinoma. Cancer Lett 412:108–117 Zhang L, Cheng R, Huang Y (2017) MiR-30a inhibits BECN1-mediated autophagy in diabetic cataract. Oncotarget 8(44):77360–77368 Xiao B, Shi X, Bai J (2019) miR-30a regulates the proliferation and invasion of breast cancer cells by targeting Snail. Oncol Lett 17(1):406–413 Liu J, Wang X, Deng Y, Yu X, Wang H, Li Z (2021) Research Progress on the Role of Regulatory T Cell in Tumor Microenvironment in the Treatment of Breast Cancer. Front Oncol 11:766248 Retecki K, Seweryn M, Graczyk-Jarzynka A, Bajor M (2021) The Immune Landscape of Breast Cancer: Strategies for Overcoming Immunotherapy Resistance. Cancers (Basel) 13(23) Nishikawa H, Koyama S (2021) Mechanisms of regulatory T cell infiltration in tumors: implications for innovative immune precision therapies. J Immunother Cancer 9(7) Huang P, Zhou X, Zheng M, Yu Y, Jin G, Zhang S (2023) Regulatory T cells are associated with the tumor immune microenvironment and immunotherapy response in triple-negative breast cancer. Front Immunol 14:1263537 Picca A, Guerra F, Calvani R, Coelho-Junior HJ, Landi F, Bucci C, Marzetti E (2023) Mitochondrial-Derived Vesicles: The Good, the Bad, and the Ugly. Int J Mol Sci 24(18) Hu X, Chen W, Li F, Ren P, Wu H, Zhang C, Gu K (2023) Expression changes of ER, PR, HER2, and Ki-67 in primary and metastatic breast cancer and its clinical significance. Front Oncol 13:1053125 Davey MG, Hynes SO, Kerin MJ, Miller N, Lowery AJ (2021) Ki-67 as a Prognostic Biomarker in Invasive Breast Cancer. Cancers (Basel) 13(17) Goncalves AC, Richiardone E, Jorge J, Polonia B, Xavier CPR, Salaroglio IC, Riganti C, Vasconcelos MH, Corbet C, Sarmento-Ribeiro AB (2021) Impact of cancer metabolism on therapy resistance - Clinical implications. Drug Resist Updat 59:100797 Kim J, DeBerardinis RJ (2019) Mechanisms and Implications of Metabolic Heterogeneity in Cancer. Cell Metab 30(3):434–446 Additional Declarations No competing interests reported. Supplementary Files supplementary.docx Cite Share Download PDF Status: Posted Version 1 posted You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-6960401","acceptedTermsAndConditions":true,"allowDirectSubmit":true,"archivedVersions":[],"articleType":"Research Article","associatedPublications":[],"authors":[{"id":475774483,"identity":"ae9c7d09-e35b-4655-992a-6747fc5a76a3","order_by":0,"name":"Wenbin Guo","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA60lEQVRIiWNgGAWjYBACNvnDxz98MJCQ42fvf/ggoaKGsBY+CbY0xhkFNsaSPWeYDR6cOUZYi5wEjxkzz4e0xA0zctgkH7YwE+Ew6ba0xzwGhxM38Jw9VpHYwMbA396dgF+LzOHjhnMMDhtvZ+9Lu5G4Q4ZB4szZDfi1MKQlSLwxOCy7s+eA2Y3EM2wMBhK5hLTkGEgAHca44UaCWUFiGzMRWiRyzCR5DNIUN9zIMWMgTgvPsWTDGQagQD6WLJFw5hgPQb/ItzcffPDhDygqmw9+/FFRI8ff3otfCwbgIU35KBgFo2AUjAKsAACYx07GvQM8ygAAAABJRU5ErkJggg==","orcid":"","institution":"Pingtan Comprehensive Experimental Area Hospital","correspondingAuthor":true,"prefix":"","firstName":"Wenbin","middleName":"","lastName":"Guo","suffix":""},{"id":475774484,"identity":"a09527b8-5f23-40d0-a142-760d0545b41c","order_by":1,"name":"Zhiqiang Ye","email":"","orcid":"","institution":"Fuzhou University Affiliated Provincial Hospital (Fujian Provincial Hospital)","correspondingAuthor":false,"prefix":"","firstName":"Zhiqiang","middleName":"","lastName":"Ye","suffix":""},{"id":475774485,"identity":"c26f14df-7e54-4377-a5de-470c9fe6307c","order_by":2,"name":"Dan Wu","email":"","orcid":"","institution":"Fujian Medical University Union Hospital","correspondingAuthor":false,"prefix":"","firstName":"Dan","middleName":"","lastName":"Wu","suffix":""}],"badges":[],"createdAt":"2025-06-24 01:23:15","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-6960401/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-6960401/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":85344646,"identity":"26e6d45e-8cc8-43d0-9ce7-3c3b22c76160","added_by":"auto","created_at":"2025-06-25 01:47:15","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":887232,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eMitophagy-related genes (MRGs) express differently in breast cancer (BC).\u003c/strong\u003e (A) Box plot showing the 'GSVA' scores for MRGs in each sample, emphasizing the variations in overall expression between the BC and normal samples. (B) A volcano plot showing the 28 MRGs that are differentially expressed in BC; upregulation is indicated by positive log2FC values and downregulation is indicated by negative log2FC values. (C) A box plot showing the differential expression of 28 MRGs in normal and BC samples. (D) Box plot showing how the expression of 22 immune cell types differs in BC samples compared to normal samples. The Wilcox test was utilized to compute the P-values. ***P\u0026lt;0.001, ****P\u0026lt;0.0001), **P\u0026lt;0.01, *P\u0026lt;0.05)\u003c/p\u003e","description":"","filename":"image1.png","url":"https://assets-eu.researchsquare.com/files/rs-6960401/v1/c7e9a166337468e1dec26b77.png"},{"id":85345079,"identity":"7c511edc-6b7a-4815-8ef5-b816899cc360","added_by":"auto","created_at":"2025-06-25 01:55:15","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":546410,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eExpression of mitophagy-related genes (MRGs) based on single-cell RNA sequencing (scRNA-seq) data.\u003c/strong\u003e(A) Quality control measures for the two samples included in the analysis. (B) Dot plot presenting the expression of biomarkers across different cell clusters. The cells were divided into 8 different types using (C) Uniform manifold approximation and projection (UMAP) dimensionality reduction techniques, where every color corresponded to an identified phenotype for each cluster. (D) MRG scoring in different cell types using AddModuleScore as a basis. (E) The high-scoring group's MRG scoring in various cell types is shown in the left panel, while the low-scoring group's scoring is shown in the right panel. (F) The KEGG pathway ribosome.\u003c/p\u003e","description":"","filename":"image2.png","url":"https://assets-eu.researchsquare.com/files/rs-6960401/v1/f6fbcdd9c3c64046633d34e7.png"},{"id":85344647,"identity":"b0d9855d-96e4-4f1f-91c6-b8d512a79db5","added_by":"auto","created_at":"2025-06-25 01:47:16","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":341375,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eDiagnostic model establishment and assessment. (\u003c/strong\u003eA) The Random Forest model's error rate curve, which shows how prediction error varies with the number of trees. (B) A ranking map showing the genes' significance scores as determined by the Random Forest model. Higher values in the MeanDecreaseGini index indicate a gene's larger significance in the model, which is used to quantify the importance of feature genes. A forest plot (C) displaying the nine genes that were chosen using logistic regression. (D) A Receiver Operating Characteristic (ROC) curve that assesses how well feature genes in the training data perform as diagnostics. The real positive rate is shown by the vertical axis, and the false positive rate by the horizontal axis. (E) ROC curve assessing feature genes' diagnostic performance based on validation data. (F) Calibration curve showing the predictive model's calibration performance. The observed value is represented by the vertical axis, and the expected value is shown by the horizontal axis. The model's calibration performance increases with the distance between the bias-corrected curve and the ideal dashed line. (G) The nomogram's clinical benefit is estimated using decision curve analysis (DCA). The net clinical benefits for several prediction algorithms at varying choice thresholds are plotted and compared. The DCA curves demonstrated that our model provided a higher net benefit compared to other strategies across a wide range of threshold probabilities (from 0.2 to 1). Specifically, the net benefit of our model was consistently superior to alternative strategies, indicating its potential to guide clinical decision-making by reducing unnecessary interventions while appropriately identifying patents who would benefit from treatment. (H) Nomogram for using feature genes to predict breast cancer (BC) risk.\u003c/p\u003e","description":"","filename":"image3.png","url":"https://assets-eu.researchsquare.com/files/rs-6960401/v1/b6ad3eb5f78e596dd3ebd591.png"},{"id":85344664,"identity":"e5ecd51c-affb-49ed-ae3d-5142a9019972","added_by":"auto","created_at":"2025-06-25 01:47:16","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":508445,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eIdentification and differentiation of mitophagy-related breast cancer (BC) subtypes. (\u003c/strong\u003eA) Hierarchical clustering with k (number of clusters) = 3 showing the agreement in clustering between samples. (B) Box plot showing how three BC subtypes differ in how nine mitophagy-related genes (MRGs) are expressed. (C) Heatmap showing the expression of genes associated with mitophagy as well as clinical data for various BC subtypes. (D) Box plot comparing three BC subtypes' diagnostic performance using prediction models. (E) Box plot showing the MRG gene set's gene set variation analysis (GSVA) scores across the three BC subtypes. (F) Survival curves showing the variations in survival between the three BC subtypes. P-values were calculated using the log-rank test. (P-values were estimated by the Kruskal-Wallis test. **P \u0026lt; 0.01, ***P \u0026lt; 0.001, ****P \u0026lt; 0.0001)\u003c/p\u003e","description":"","filename":"image4.png","url":"https://assets-eu.researchsquare.com/files/rs-6960401/v1/5d8eafc106bccce4bca643c9.png"},{"id":85344662,"identity":"282956e5-5b98-473f-95e9-e5a3b59ddad9","added_by":"auto","created_at":"2025-06-25 01:47:16","extension":"png","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":349268,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eBuilding a gene co-expression network with WGCNA in several breast cancer subtypes associated with mitophagy. \u003c/strong\u003e(A) A scale-free topology plot of WGCNA with a cutoff value of 0.9 is used to determine the ideal soft-thresholding power. (B) Eigengene adjacency heatmap showing how several modules are correlated. (C) A heatmap of clinical features and a dendrogram of samples. (D) Correlation matrix showing the correlation coefficient and matching p-value for each module in respect to various BC subtypes (clusters).\u003c/p\u003e","description":"","filename":"image5.png","url":"https://assets-eu.researchsquare.com/files/rs-6960401/v1/4ce0c4ca6db6480e4ed7a875.png"},{"id":85344648,"identity":"8f2cf135-5b17-4744-b1ff-cd2e11299781","added_by":"auto","created_at":"2025-06-25 01:47:16","extension":"png","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":520009,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003ePredicting immunotherapy outcomes and chemotherapeutic drug sensitivity in various subgroups.\u003c/strong\u003eBased on the TIDE database, box plots show the predicted values for the Tumor Immune Dysfunction and Exclusion (TIDE) score (A), Exclusion score (C), Dysfunction score (B), and Microsatellite Instability (MSI) score (D) in three BC subtypes. (E–J) Box plots with the half-maximal inhibitory concentration (IC50) values on the y-axis show the sensitivity of three mitophagy-related BC subtypes to the following drugs: cisplatin, doxorubicin, paclitaxel, gemcitabine, cyclophosphamide, decitabine+carboplatin, and methotrexate. (*P \u0026lt; 0.05, **P \u0026lt; 0.01, ***P \u0026lt; 0.001, ****P \u0026lt; 0.0001, P-values were estimated by the Kruskal-Wallis test)\u003c/p\u003e","description":"","filename":"image6.png","url":"https://assets-eu.researchsquare.com/files/rs-6960401/v1/a1330f6b7febc3500ad96908.png"},{"id":85403311,"identity":"7a010c1c-8498-4a18-84d0-73672a80fa93","added_by":"auto","created_at":"2025-06-25 12:32:09","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":4340031,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-6960401/v1/a4fcd76f-c593-4a41-ac73-6a56fef54329.pdf"},{"id":85344668,"identity":"062e5aee-ef8b-4bdf-bf69-d0cc1ae2a75c","added_by":"auto","created_at":"2025-06-25 01:47:17","extension":"docx","order_by":0,"title":"","display":"","copyAsset":false,"role":"supplement","size":32472119,"visible":true,"origin":"","legend":"","description":"","filename":"supplementary.docx","url":"https://assets-eu.researchsquare.com/files/rs-6960401/v1/c3b6394382d75aed8023a512.docx"}],"financialInterests":"No competing interests reported.","formattedTitle":"Mitophagy-Driven Immune Cell Infiltration Patterns Define Breast Cancer Subtypes with Differential Treatment Responses","fulltext":[{"header":"Introduction","content":"\u003cp\u003eBreast cancer remains one of the most prevalent malignancies affecting women worldwide, contributing significantly to cancer-related morbidity and mortality. Despite advancements in early detection and treatment, the high incidence of therapeutic resistance and metastasis poses substantial challenges in breast cancer management. A critical aspect of breast cancer pathology is the involvement of mitochondrial dysfunction and the subsequent process of mitophagy, which has garnered increasing attention for its role in cancer progression and treatment resistance[\u003cspan additionalcitationids=\"CR2 CR3 CR4\" citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eMitochondria are essential organelles involved in various cellular processes, including energy production, apoptosis, and redox homeostasis. Mitochondrial dysfunction is a hallmark of cancer cells, contributing to tumorigenesis and cancer progression by altering cellular metabolism and promoting resistance to apoptosis[\u003cspan additionalcitationids=\"CR7 CR8 CR9\" citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e]. In breast cancer, oxidative stress-induced damage to mitochondrial DNA (mtDNA) and the resultant mitochondrial dysfunction play crucial roles in cancer cell survival and proliferation[\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e]. The regulation of mitochondrial dynamics and the maintenance of mitochondrial integrity through processes such as mitophagy are essential for cancer cells to adapt to their rapidly changing environment.\u003c/p\u003e \u003cp\u003eMitophagy, an evolutionarily conserved cellular process, plays a dual role in cancer by acting as both a tumor suppressor and promoter depending on the context and stage of cancer development [\u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e, \u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e]. Dysregulated mitophagy can lead to the accumulation of damaged mitochondria, thereby contributing to carcinogenesis and tumor progression. In breast cancer, mitophagy has been shown to facilitate the metabolic reprogramming of cancer cells, promoting a switch between aerobic glycolysis and oxidative phosphorylation to meet the high energy demands of proliferating tumor cells [\u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e]. Key regulators of mitophagy, including \u003cem\u003ePINK1, Parkin, BNIP3, and NIX\u003c/em\u003e, have been implicated in various cancers [\u003cspan additionalcitationids=\"CR12 CR13\" citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e]. These proteins help maintain mitochondrial quality by targeting damaged mitochondria for autophagic degradation. The interplay between these mitophagy regulators and the broader mitochondrial dynamics is critical for cancer cell survival, particularly under stress conditions such as hypoxia and nutrient deprivation commonly found in the tumor microenvironment.\u003c/p\u003e \u003cp\u003eThe role of mitophagy in breast cancer is increasingly recognized. The field of mitoepigenetics, which focuses on the epigenetic regulation of mtDNA, has emerged as a significant area of research[\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e]. Mutational damage to mtDNA and defects in mitochondrial biogenesis and dynamics contribute to the oncogenic phenotype of breast cancer cells. Studies have shown that cancer cells often exhibit defects in mtDNA repair mechanisms, leading to increased reactive oxygen species (ROS) production and oxidative stress, which in turn promotes further mitochondrial damage and dysfunction. This creates a vicious cycle that fosters cancer cell survival and proliferation. Therapeutic targeting of these mitochondrial vulnerabilities, including the enhancement of mitophagy, represents a promising strategy for cancer treatment. Several mitochondrial-targeting agents, such as biguanides and oxidative phosphorylation (OXPHOS) inhibitors, are being explored for their potential to disrupt this cycle and inhibit breast cancer progression[\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e, \u003cspan citationid=\"CR15\" class=\"CitationRef\"\u003e15\u003c/span\u003e]. Targeting mitophagy in cancer therapy holds significant promise. By modulating mitophagy, it may be possible to eliminate cancer cells more effectively and overcome resistance to conventional therapies. For instance, the use of mitophagy inducers or inhibitors could help maintain mitochondrial homeostasis and enhance the efficacy of existing treatments. Furthermore, the identification of biomarkers related to mitophagy and mitochondrial dysfunction could aid in the early diagnosis and personalized treatment of breast cancer[\u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eIn this study, We undertook a comprehensive investigation of mitophagy-related gene expression in breast cancer. We then created a predictive model for reliable breast cancer diagnosis based on mitophagy-associated signature genes. This model was validated using a subset of the TCGA dataset. Breast cancer subgroups were identified using hierarchical clustering. In addition, we conducted extensive research into the prognosis, immunological infiltration, immunotherapy response, and drug sensitivity of each discovered breast cancer subtype.\u003c/p\u003e"},{"header":"Materials and methods","content":"\u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003e2.1 Data acquisition and processing\u003c/h2\u003e \u003cp\u003eThe Cancer Genome Atlas (TCGA) database provided RNA sequencing (RNA-seq) data and related clinical details for 1,119 breast cancer (BC) samples. The R tool \"TCGAbiolinks\"(version\u0026thinsp;=\u0026thinsp;2.32.0) was used to preprocess and standardize the collected data. A training dataset and a test dataset were randomly selected from the BRCA samples at a 7:3 ratio. The remaining 336 samples were used for internal validation, leaving 783 samples in the training cohort. We used single-cell RNA-seq data from the Gene Expression Omnibus (GEO) database, accession number GSE248288, to further explore the function of mitophagy in breast cancer[\u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e16\u003c/span\u003e]. Additionally, a comprehensive list of mitophagy-related genes was retrieved from the Reactome database and cross-referenced with the GeneCards database (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.genecards.org/\u003c/span\u003e\u003cspan address=\"https://www.genecards.org/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e) to ensure the inclusion of all relevant genes in our analysis.\u003c/p\u003e \u003c/div\u003e\n\u003ch3\u003e2.2 Differential expression analysis of mitophagy-related genes\u003c/h3\u003e\n\u003cp\u003eTo investigate the expression patterns of mitophagy-related genes (MRGs) in breast cancer, we first classified the samples from TCGA datasets into normal and tumor groups based on the available clinical information. The MRGs were obtained from the Reactome database, and their expression levels were visualized using a heatmap generated with the R package \"pheatmap\"(version\u0026thinsp;=\u0026thinsp;1.0.12). Next, we performed differential expression analysis to identify differentially expressed genes (DEGs) between the normal and tumor samples (the threshold was set as log2 fold change\u0026thinsp;\u0026gt;\u0026thinsp;0 or \u0026lt;\u0026thinsp;0 and adjusted p value\u0026thinsp;\u0026lt;\u0026thinsp;0.05). The differential expression analysis was conducted using the R package \"DESeq2\"(version\u0026thinsp;=\u0026thinsp;1.44.0), which employs a negative binomial distribution to model the read counts and provides a robust statistical framework for identifying DEGs. To visualize the results of the differential expression analysis, we constructed boxplots and volcano plots using the R packages \"ggplot2\"(version\u0026thinsp;=\u0026thinsp;3.5.1) and \"ggpubr\"(version\u0026thinsp;=\u0026thinsp;0.6.0). The boxplots allowed for a comparison of the expression levels of MRGs between normal and tumor samples, while the volcano plots provided an overview of the distribution of DEGs based on their log2fold change and adjusted p-values. Furthermore, to explore the relationships among the identified DEGs, we generated a correlation heatmap using the R package \"corrplot\"(version\u0026thinsp;=\u0026thinsp;0.94). This heatmap visualized the pairwise correlations between the DEGs, enabling us to identify potential co-expression patterns and functional relationships among the MRGs in breast cancer.\u003c/p\u003e\n\u003ch3\u003e2.3 Analysis of scRNA-seq data in BC\u003c/h3\u003e\n\u003cp\u003eTo analyze the single-cell RNA sequencing (scRNA-seq) data in breast cancer (BC), we utilized the R package \"Seurat\"[\u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e, \u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e18\u003c/span\u003e](version\u0026thinsp;=\u0026thinsp;5.1.0). To assure high-quality data, we performed a series of filtering operations on the raw matrix. We specifically kept genes expressed in at least three cells, deleted cells expressing fewer than 200 genes, and excluded cells with mitochondrial gene expression greater than 25%. The scRNA-seq data was normalized using the \"NormalizeData\" function, which was then used to identify the top 2,000 highly variable genes using the \"FindVariableFeatures\" function. Principal component analysis (PCA) was then done with the \"RunPCA\" function from the \"Seurat\" package to reduce the dimensionality of the scRNA-seq data based on these highly variable genes. To integrate the two samples, we used the canonical correlation analysis (CCA) approach. ElbowPlot analysis was used to identify significant principal components, and cell clustering was done with the top ten principal components using PCA. The \"FindNeighbors\" function was used to create a k-nearest neighbor graph in the PCA space using the Euclidean distance, and the Uniform Manifold Approximation and Projection (UMAP) approach was used to reduce data dimensionality and visualize it. Cluster annotation was done with reference data from a previously published study[\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e]. To evaluate the gene expression scores of mitophagy-related genes (MRGs) for each cell in the scRNA-seq dataset, we used the \"AddModuleScore\" function to generate individual cell scores. In addition, we performed Gene Set Variation Analysis (GSVA) to examine the enrichment of signature pathways for each cell type using the Molecular Signatures Database (MSigDB)[\u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e20\u003c/span\u003e, \u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e]. Finally, we conducted Gene Set Enrichment Analysis (GSEA) on the differentially expressed genes between the high and low groups based on mitophagy scores.\u003c/p\u003e\n\u003ch3\u003e2.4 Construction and validation of a prediction model\u003c/h3\u003e\n\u003cp\u003eMachine learning models have developed as effective methods for dealing with enormous datasets and aiding rapid analysis and decision-making in disease diagnosis. In this study, we used a combinatorial technique leveraging support vector machines (SVM)[\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e] and random forest (RF) [\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e] algorithms to pick significant features from a list of 29 mitophagy-associated genes. Following that, a prediction model was created using logistic regression[\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e]. The prediction model's performance was examined using receiver operating characteristic (ROC) curve analysis to determine its sensitivity and specificity. The area under the curve (AUC) of the ROC curve was used to assess the model's prediction performance. To display the relative importance of the feature genes within the model and their contribution to the prediction results, a nomogram model was created using \"rms\" package (version\u0026thinsp;=\u0026thinsp;6.8-2). To further validate the model's prediction performance, calibration curves were generated to determine the agreement between anticipated probabilities and observed outcomes. In addition, decision curve analysis (DCA), which involves calculating the net benefit of a predictive model at different decision thresholds while taking into account the false positives, false negatives, and costs associated with different decision thresholds, was used to assess the model's clinical value by estimating the net benefits at various probability thresholds[\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eThe net benefit is calculated using the formula:\u003c/p\u003e\n\u003cp\u003e\u003cimg src=\"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAacAAAAwCAYAAAC1+7XXAAAAAXNSR0IArs4c6QAAAARnQU1BAACxjwv8YQUAAAAJcEhZcwAADsMAAA7DAcdvqGQAAA3HSURBVHhe7Z3JixRJG4fT767S6klFBtuLIOjBDRFhBB1REEFwRzwIop4G3KY9K7Z3HQUFD4or6EHFBRTcaFcUHGYO6kWPOqMzf0B/9UTnr+ad6KysqraW7Kr3gajMjIyMeGN9IyMjosYMlkgcx3Ecp0D8Lz06juM4TmFw5eQ4juMUDldOjuM4TuFw5eQ4juMUDldOjuM4TuFw5eQ4juMUDldOjuM4TuFw5eQ4juMUDldOjuMM4+jRo8mYMWOSGTNmpDaO01pcOY1yaEDyzJs3b1KXrUcNnDU0dhcuXEhdNIabN28Gv798+ZLa5DNv3rzk4MGD6VVnsWLFimFpLvPo0aPUVXX27duX/PTTT8ny5ctTm8aTJSt50+gyS17jby1QhpCDMtWtUEez8qXVaeLKaRSjSvzw4cOEXajmzp2b7Ny5M5xjB7Nnzw7HdkADZ2X6/PlzaOw2btyYfPjwIXX1/bx9+zbp6elJJk6cmNokwX8qVRYvXrxIpk2bll51FufOnQvH8+fPhzS3ZWHx4sXhWA/NLD+xrO/fvw/Xa9euDcdG8fLly6S3tze9GoIOEsox5vfffw/HKVOmhGM3cubMmXBUu0K+ULdWrVpVbnPy6lfDKAXujFJKlXqwv78/nJcafvZIDHai1PNNz9pHLFOpwAc7js2kpBAHS4oxveoelL6UB8tI0gJ/Xr9+nV41HslaavxSm8FQnrFrJqQN6dHX15faOJasOio7tTc3btxoej65cuoQVFhsRRfcK71ZDJZ6j6nNECgvVVAUCG7UiGHPdazgKJz4Q1gcs8ITWQ2lGh9rRwO4fv36YI/hPL6PXLpPuKo4siOOQnbWgMK2DTV+4x9xVVyIO+5s44U8KDzccS9OF6WXwsNtOyCOpF81bD6qHNg0JH2JjyUvH6BaGsUgg80L4JnYDrkULn7bfAHKruKCwS2yqPxhVJ6snYzk5Mi19T+uF/gjdza94jJM3ATP2Hv4d+LEifRu8VAa2byVHemhNFF8MDa+jcKVU4dAhYorNdDg0mCocFm4poLFbvCLykaBsxWVSonBPRWOBiGvUMYyqVBbP/ELO8LHT/kbh6vGXvelSCQz9hbcEF4Mflm/kV9+4CfpwXPY2zAlN+ekDe55DuzzQMOjZ1uNGnFr4rQhHsRH8hPX2B12pJUlLx+qpVEWyKq84BnJYRt9zrFTYy5/rRuuldfcRy6Bn3G9QGaekewW6zf3VRbkJ+eEQVooTOJI3CUjR/wRkkHpi6LKKptFQeVZeUc6EF/iqDjITTNx5dQh2IqeBZXBVloVLhU2wA12VL4YPY97CqsqYFYFF3FDyfNUVAsVNW48qAi2YcQN11lhIVc9jQ8y2IYN1ABaSEtVThpkZADShmublmow293gkDfIofghf6wkJaviBuRlVh5k5VWlfKiWRjGS1RpkiPNGCs+CWyubVYoxyBE/H9cFoXIT+8PztjwC13KHX1JMpCvxsOnOPdzYNC8yyGnzhfQlb227IDfNxJVTB6BKFVdsC4XLVlLO4wbJNjAxVEYVVioabm1hjVHjU63Bxo0qtsB/2/jgF3LhNm4wsxofNQYxlRof/LSND+5sOvCMDGlGeLEfpD2VmPt56dJMpGRj2SzENU6buGwA/sSNaV4+KH0wldLIIlmzFJ1QA2jTU3loZeM+YZL+cR0grrEdcbDKQ2R1dICyYeOLf7qOG3LcxuUZcM99ws5LlyKgOMX5b5GbZuLKqQOgUlFQ8gq9rbhUcPU2LVTkrIoF+B83SHlIprzGh0YlrgQq9FkNvNzbxqYRjY91TxraXnEtFdWCX1lhtALiUC1s4mYVcV4eVCLOh3rTCJA1qwNhoSzGclAGKbtZyL3yToosrhc8n9VpQqa4ToB1T9zJYyGlUwvIQZyzwigSteRntTLSCHwqeQegqZ92KrWFaZ9//fVX+fz06dPJ/Pnzk/Hjx4cptax/wb5UmZOFCxcGdzGlRi959epVOGctCGuY8tYKSabp06eHYxZMUy5V/OTJkyfhGjlWr14dpp5rCjNrLjR9dezYseEoJPO4ceOSXbt2ldc5ff36NfjLNfZ6/sGDByGNuEZ+gfsffvghhL958+ak1OCU03LmzJnh+Ntvv4UjYeKn1mpx3LBhQziHCRMmpGethbgiS6UyIMiPd+/eBffEd2BgINhPnjw55Cf2iqvSD/LyoVoaxUjWagt8VRa1PuvkyZPJ/v37k+PHj4drwsEP/APKs+XZs2fpWVKOB1AXKDM8b6eT37lzp1wn7Joe3E+dOjXYE/6xY8fSO0kya9ascFTaICt+4jcQLnILyuVo4e+//07PKkN8iXtT1g2mSsoZpajnhrE94hju4YbeIT04em/0CPWmRM+Q60rQY6TXhx+4kz9Z1CoTEK785cizFsmtcO199USxt29PCj92r551PLRCXPLipLTJktHKj0Fe0qrV2HSK09CifLRpw9sWdpKbXjH3sZddXj5AXhrF1CorqJziFnlsPiMbdvKLcO39rHiAwo/dUy6yZMJdXpyQUTLEMlr5MXn1pigQV2RFbptuMUpHjs2I0xh+SgE4juM4TmHwYb0RwCsswzeskNaQg+M4jtM4RqScGEPVnksaawXOZV9pvHkkMG4rf61hbNeG3wr4VsF2KH/++Weyfv368L0EWVxJOY7jNI4RKacdO3aETSF7enqS7du3p7ZDH7j7+/vDB237kTgPPhzSuOehD6D4zSgk5vXr1+HD7o8//hjutYrLly8nS5cuDeco4EWLFoVzC8rUfnB3HMepBSZi2MkbFjriWfsBdirfNax34MCBsIlm3BAz86lW/vjjj/SsPlCEvLm0GuJrYTNNlKXdVPP58+fpmeM4Tu0sWLAgtB+xgkIx0RHftm1batP5fJdy4q2hr68vTK/MG16jN8CW67whMfVTQ34ct2zZEs41VFfrGwd+/vrrr+W3KtDUV30P4u1N00wJCzt6HjyLHFzHUyCJB89xD3/sfW27T3y5z1Ce/NGwHm5QYHKDyUL38oy/fTlOd8FSgFu3bv1HQUkx0dbVOiLVETBbbyQwfZCpmprOyxRKYLqlnXKphXBMNQWmWXKtxZm1LubCTWziqZ1M00QOZMIgI3aCa6ZHavq0phbb6ZLERc8gI+4lO8ThSn6OgnBi2RzHcWqF9ou2jOntcRvULXy3cgIpHCkm2zAzz1+KS6AApCDqUU7WX8Ik0/AfyEzc2DUGWnshkBkj9EwcD+wFiooCImI52qmcCNdNsU01sp6JTS1lKes5N8U3eahtsW1WLcRhFN1UoiFTyVeuXBm+/xw5ciT59u1bajsEM9sY5rLDVazqj93VC2Hyzevw4cPhWn8Sxh9iKRz+1E47I2QRr6bnT+tg0qRJZT8YOmQHgUYj//NMtWG9Uv65KbipRtYzseFPG6uR9Zyb4ptK8DmC70u0q+x0Uc/s56xwimwq0RDlBNrSg8bcwnYpJc0/TKBaKlw9sP0KlHobw8KqFf07avw8Y8CNJg4jyzQ6jRzHKT4oJr6N86/RKKX79++H70/1KKhO4LuU06dPn9KzobeQs2fPDntT2bRpU3L79u3y/lLad0t7T4l692hiUgNvaig+kBL8+eefy37jJs8/TZbQHlLLli0L0+ORj3sY5B5Jobh37144ZsXVcRwnCykmJlZpshczk7tSQZV66HWjiQR8z4n3Xsr63sK19mviaL8LAc9wj6P93iP47sN9awgbe+seWfg+JDf2Pt+fZK/vXZLJxgPZ+EYmexsXK4c+UFo/hL6HYbrxQ6bjOCODtoN2JgvaKNrIbsH31nMcx3EKR8O+OTlOvTCUq8kfdrjC2vtar+5Dec8wlkX23bRLQiPg80az04y8Uv7IsAb0u+pveH9ynDahvxuIhyuw76YhDOdftKSDYXE7bK/h8qyhf2c4pJM+c7SiLhGOPoMQtuq2PqMIhi1rWR7hb05OW2GpQamwhkkzduIIyw26aasW519Y0kGZYHLV3bt3U9sk+fjxY9ghodofKjpD7N69O/nll1/CnqSthjw6dOhQOL969Wo4ApPeGCXRnzTm4crJaSsopT179oRz2xCxNo5/63W6D2a6Msu3t7c3uXbtWmo71LAtWbIkvXKqgRLQP0q3E/0bNkN8c+bMCZ0OrUfNG/Zz5eS0DcbCWWRI4eV46tSpYE8jBCrUTvfAVGoWnbKR8t69e5OLFy+Wl3zQ2HqHZXRAnmkZDx0NYN3m+fPnQ6djsIa1nK6cnLbx+PHj8l+erFmzJrwtoZgGBgbCAkSn+3j69Gk571l3CFeuXAnlgo2YvcPSPOIJDVmm2gQHbXjNLjsM2d+4ceM//9jw4MGDmuu2KyenbdArXrhwYTjnW0JPT09y6dIlH77pYq5fv17usNg3arYn8w5Lc9HbTJ6ptmsN37fklp112GbOwttvrXXblZPTFjT5wY6J8xEchXXnzh0fvulSyHu9McHWrVvDG/WZM2e8wzLKodPJ96Za67YrJ6ctMFTDm5Jl3bp1YZYe+PBN98H6NvJ//Pjxqc3QBs+UEybOeIdlZLDJNv8a3gryNvT+559/wpH8Ja/jdWzDKL1+OU5L0V8BYOL1DmwdVWn7FqezUZnIWvPGNmFOfZBu2ooNQxo2a70TdVbhVNqyjf/H03ZvyBFvfRfj2xc5juM4hcOH9RzHcZzC4crJcRzHKRyunBzHcZzC4crJcRzHKRhJ8n9MMBQDsWTBbgAAAABJRU5ErkJggg==\" height=\"48\" width=\"423\"\u003e\u003c/p\u003e\n\u003cp\u003eWhere:\u003c/p\u003e \u003cp\u003eN is the total number of patients; Ture positives and False Positives are derived from the model\u0026rsquo;s predictions at a given threshold probability Pt.\u003c/p\u003e\n\u003ch3\u003e2.5 Analysis based on the Human Protein Atlas database\u003c/h3\u003e\n\u003cp\u003eTo further validate the prediction model, we utilized the Human Protein Atlas (HPA) database (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.proteinatlas.org/\u003c/span\u003e\u003cspan address=\"https://www.proteinatlas.org/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e)[\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e, \u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e] to investigate whether the expression levels of the feature genes differ between breast cancer tissues and normal tissues at the protein level. These resources provide comprehensive, publicly available data on protein expression across various human tissues and cancer types.\u003c/p\u003e \u003cdiv id=\"Sec8\" class=\"Section2\"\u003e \u003ch2\u003e2.6 Correlation between feature genes and immunity\u003c/h2\u003e \u003cp\u003eTo investigate the potential relationship between the feature genes and the immune microenvironment, we used the \"CIBERSORT\" package (version\u0026thinsp;=\u0026thinsp;1.03) to estimate the level of immune cell infiltration in both normal and tumor samples. This computer method enables the deconvolution of complex gene expression data and the measurement of different immune cell subpopulations within tissue samples. In addition to looking at immune cell infiltration patterns, we investigated the links between feature genes, immune infiltration, and important inflammatory variables. The findings of these investigations were visually depicted using the \"ggplot2\" software.\u003c/p\u003e \u003c/div\u003e\n\u003ch3\u003e2.7 Hierarchical Clustering\u003c/h3\u003e\n\u003cp\u003eTo investigate the potential subgroups within breast cancer (BC) samples, we employed hierarchical clustering analysis using the \"hclust\" (\u0026lsquo;stats\u0026rsquo; R package, version\u0026thinsp;=\u0026thinsp;4.4.0) R function on a cohort of 783 BC patients. The analysis was based on the expression profiles of nine feature genes, which were selected through a rigorous feature selection process. The parametric h was set to 4.5 to ensure optimal clustering performance. Following the identification of distinct BC subtypes, we explored the relationship between the expression levels of the nine feature genes, hereafter referred to as the molecular risk genes (MRGs), and the clinical pathological characteristics within each subtype. To quantify the relative expression of MRGs in each sample, we utilized the gene set variation analysis (GSVA) package (version\u0026thinsp;=\u0026thinsp;1.52.3) to calculate sample-specific scores. The Kruskal test was then applied to assess the statistical significance of the differential expression of MRGs among the three identified BC subtypes. To evaluate the prognostic value of the identified subtypes, Kaplan-Meier survival analysis was conducted using the \"survival\"(version\u0026thinsp;=\u0026thinsp;3.5-8) and \"survminer\"(version\u0026thinsp;=\u0026thinsp;0.4.9) R packages. This analysis aimed to determine the impact of the BC subtypes on patient outcomes, providing insights into the potential clinical utility of the identified subgroups.\u003c/p\u003e\n\u003ch3\u003e2.8 Analysis of immune infiltration in BC subtypes\u003c/h3\u003e\n\u003cp\u003eWe utilized the CIBERSORT program in R to calculate the proportions of 22 immune cell subtypes among the identified BC subtypes. The LM22 gene expression matrix, which includes 547 genes identifying 22 human immune cell types, used as a reference for immune cell characterization. The CIBERSORT algorithm generated immune cell infiltration profiles for each BC subtype. In addition, we evaluated the expression levels of immune and stromal cells in the samples. The distribution of immune and stromal cells among BC subtypes was compared and displayed using box plots.\u003c/p\u003e \u003cdiv id=\"Sec11\" class=\"Section2\"\u003e \u003ch2\u003e2.9 Weighted Gene Co-expression Network Analysis\u003c/h2\u003e \u003cp\u003eTo identify hub genes within each breast cancer (BC) subtype, we performed consensus module analysis using the \"WGCNA\" package in R[\u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e](version\u0026thinsp;=\u0026thinsp;1.73). First, genes were ranked based on their median absolute deviation (MAD) values, and the top 5,000 genes were selected for further analysis. Subsequently, samples and genes with excessive missing values were filtered out to ensure data quality. A similarity matrix was constructed using the selected genes, which was then transformed into an adjacency matrix using a soft-thresholding power (β) of 5. This transformation emphasizes strong correlations and penalizes weak correlations, resulting in a scale-free network topology. Next, the \u0026ldquo;blockwiseModules\u0026rdquo; function was applied to perform modular analysis on the adjacency matrix. This function employs a dynamic tree-cutting algorithm to detect modules (clusters) of highly interconnected genes. The resulting gene co-expression network assigns genes to different modules based on their co-expression patterns. For each identified module, module eigengenes (MEs) were calculated as the first principal component, representing the overall expression profile of the module. To investigate the relationship between the identified modules and the phenotype data, Pearson's correlation coefficients were calculated between the MEs and the phenotype data. The correlation analysis reveals the association of each module with specific phenotypic traits, such as clinical features or outcomes. Finally, the correlation between the module genes and the phenotype data was visualized using a heat map, providing a comprehensive overview of the module-trait relationships.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec12\" class=\"Section2\"\u003e \u003ch2\u003e2.10 Analysis of enrichment\u003c/h2\u003e \u003cp\u003eTo obtain insight into the biological importance of the discovered gene modules, we chose genes from the module that had the strongest connection with each BC subtype based on WGCNA results. These gene sets were then subjected to enrichment analysis utilizing the Metascape database (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://metascape.org/\u003c/span\u003e\u003cspan address=\"https://metascape.org/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e)[\u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e]. In addition, we used the \"msigdbr\" package in R to import gene sets from the biological processes (BP) category, and the \"fgsea\" package to perform biological process enrichment analysis for different BC subtypes.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec13\" class=\"Section2\"\u003e \u003ch2\u003e2.11 Immune Checkpoint Gene Expression and Response Prediction\u003c/h2\u003e \u003cp\u003eWe evaluated the expression of immune checkpoint genes across the identified BC subtypes and used the Tumor Immune Dysfunction and Exclusion (TIDE) database (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://tide.dfci.harvard.edu/\u003c/span\u003e\u003cspan address=\"http://tide.dfci.harvard.edu/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e)[\u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e30\u003c/span\u003e, \u003cspan citationid=\"CR31\" class=\"CitationRef\"\u003e31\u003c/span\u003e] to predict the response to immune checkpoint blockade (ICB) therapy. TIDE uses gene expression fingerprints from tumors and immune cells to assess two major strategies of tumor immune evasion: malfunction and cytotoxic T lymphocyte (CTL) exclusion. By taking into account these parameters, TIDE gives a thorough assessment of the probable efficacy of ICB therapy in each BC subtype, allowing for more accurate immunotherapy response prediction.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec14\" class=\"Section2\"\u003e \u003ch2\u003e2.12 Drug sensitivity analysis\u003c/h2\u003e \u003cp\u003eHalf-maximal inhibitory concentration (IC50) is a widely used concentration measure for evaluating the therapeutic efficacy or toxicity of drugs. To assess chemotherapy sensitivity in various BC subtypes, we utilized the R package \"oncoPredict\" in conjunction with expression matrices from the Genomics of Drug Sensitivity in Cancer (GDSC) database and drug treatment data as a training set[\u003cspan citationid=\"CR32\" class=\"CitationRef\"\u003e32\u003c/span\u003e].\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec15\" class=\"Section2\"\u003e \u003ch2\u003e2.13 Cell culture, RNA extraction and quantitative Real-Time PCR (qPCR)\u003c/h2\u003e \u003cp\u003eFujian Cancer Hospital supplied the human breast epithelial cell line MCF-10A, as well as breast cancer cell lines (MCF-7, MDA-MB-231, MDA-MB-468, Hs-578T, HCC1937, and BT-549). MCF-10A cells were cultured in DMEM/F12 medium (Gibco) with 5% horse serum (Gibco), 20 ng/mL epidermal growth factor (EGF, PeproTech), 0.5 \u0026micro;g/mL hydrocortisone (Sigma-Aldrich), 100 ng/mL cholera toxin (Sigma-Aldrich), and 10 \u0026micro;g/mL insulin (Sigma-Aldrich). Breast cancer cell lines (MCF-7, MDA-MB-231, MDA-MB-468, Hs-578T, HCC1937, and BT-549) were grown in Dulbecco's Modified Eagle Medium (DMEM, Gibco) with 10% fetal bovine serum (FBS, Gibco) and 1% penicillin-streptomycin. The cells were cultivated at 37\u0026deg;C in a humid atmosphere with 5% CO₂. Cell viability and confluence were assessed daily, and subcultures were carried out at 80\u0026ndash;90% confluence using 0.25% trypsin-EDTA (Gibco).\u003c/p\u003e \u003cp\u003eTotal RNA was extracted from cells using the TRIzol reagent (Invitrogen) according to the manufacturer's instructions. A NanoDrop spectrophotometer (Thermo Fisher Scientific) was used to determine RNA content and purity, with A260/A280 ratios of 1.8 to 2.0 considered acceptable. cDNA was synthesized from 1 \u0026micro;g of total RNA using the PrimeScript RT Reagent Kit (Takara Bio) and oligo(dT) primers. SYBR Green Master Mix (Applied Biosystems) was utilized for quantitative PCR with a StepOnePlus Real-Time PCR System (Applied Biosystems). The PCR conditions were as follows: Initial denaturation at 95\u0026deg;C for 10 minutes, then 40 cycles at 95\u0026deg;C for 15 seconds and 60\u0026deg;C for 1 minute. Gene expression levels were normalized to GAPDH, and relative quantification was performed using the 2^(-ΔΔCt) approach. The primer sequences used in this study are shown in Supplementary file Table \u003cspan refid=\"MOESM1\" class=\"InternalRef\"\u003eS1\u003c/span\u003e.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec16\" class=\"Section2\"\u003e \u003ch2\u003e2.14 Statistical analysis\u003c/h2\u003e \u003cp\u003eAll statistical analyses were performed using the R program (version 4.3.2) or GraphPad Prism (version 9.5.0). The Wilcox rank sum test was used to identify significant differences between two independent groups, while the Kruskal test was used to identify significant differences between three groups. Results from qPCR were accessed using a one-way ANOVA. P-values\u0026thinsp;\u0026lt;\u0026thinsp;0.05 were considered statistically significant.\u003c/p\u003e \u003c/div\u003e "},{"header":"Results","content":"\u003cp\u003e\u003cstrong\u003e3.1 Breast cancer (BC) is associated with upregulated Mitophagy-Related Genes (MRGs) and immune system activation.\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eUsing The Cancer Genome Atlas (TCGA) training dataset, we examined the \u0026apos;gsva\u0026apos; scores of MRGs to evaluate the expression profiles of mitophagy genes in both tumor and normal tissues. Our findings showed that, in comparison to tumor tissues, normal tissues had substantially reduced MRG expression levels (Figure 1A). Seven of the 28 MRGs (CSNK2B, FUNDC1, MTERF3, PGAM5, SRC, ULK1, and TOMM40) were significantly elevated in tumor tissues, whereas eighteen were significantly downregulated. The remaining genes did not exhibit significant differences between normal and tumor tissues (Figure 1B-C, Table 1, supplement Fig.1A). Furthermore, correlation analysis was used to further understand the links between these MRGs (supplement Fig.1B). To examine the immune landscape in BC, we used the \u0026apos;CIBERSORT\u0026apos; method to calculate the infiltration of different immune cell types in tumor and normal tissues. Our findings demonstrated a considerable increase in the presence of activated T cells, helper T cells, regulatory T cells (Tregs), M0 macrophages, and M1 macrophages in tumor tissues. Tumors had considerably lower levels of na\u0026iuml;ve B cells, plasma cells, resting T cells, active NK cells, monocytes, and M2 macrophages than normal tissues (Figure 1D).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003e3.2 Single-Cell RNA Sequencing Analysis Reveals the Role of MRGs in Different Cell Types\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eWe examined the single-cell RNA sequencing (scRNA-seq) data from two breast cancer tumor samples (GSM7910990 and GSM7910991) acquired from the GSE248288 dataset in order to look into the expression and functional implications of mitophagy-related genes (MRGs) at the single-cell level. Quality control techniques were used to guarantee that high-quality cells were included in the study (Figure 2A, Supplementary Figures 2A-C). Principal component analysis (PCA) was used to minimize dimensionality by selecting the top 2,000 highly variable genes (Supplementary Figure 2D). The cells were then clustered into 17 separate clusters (Supplementary Figure 2E). Based on biomarkers from the prior study[19], we classified the cells into major cell types, including T cells, natural killer (NK) cells, epithelial cells, endothelial cells, myeloid cells, plasma cells, B cells, and fibroblasts (Figures 2B-C). To analyze MRG expression in various cell types, we used the \u0026quot;AddModuleScore\u0026quot; tool to compute the MRG score for each cell type. Remarkably, NK cells trailed B cells in having the highest MRG score (Figure 2D). Gene set variation analysis (GSVA) was carried out to look more closely at the biological mechanisms connected to MRGs in various cell types. The findings showed that B cells had an increased allograft rejection pathway, fibroblasts primarily had an increased epithelial mesenchymal transition, and NK cells had increased PI3K-AKT-mTOR signaling, cholesterol homeostasis, and Myc-targets-V1 signaling (Supplementary Figure 3A). Gen set enrichment analysis (GSEA) was performed to obtain insights into the differential gene expression between samples with high and low mitophagy-related scores (Hscore and Lscore) (Figure 2E). Pathways such the ribosome, proteasome, MAPK signaling pathway, and taurine and hypotaurine metabolism had a substantial enrichment of the differentially expressed genes across the high and low groups (Figure 2F, Supplementary Figure 3B). These results suggest that MRGs may influence protein synthesis, degradation, signaling cascades, and metabolic processes in breast cancer cells.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003e3.3 Construction and Validation of a Predictive Model for Breast Cancer Occurrence\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTo identify key mitophagy-related genes (MRGs) associated with breast cancer (BC) occurrence and develop a predictive model, we employed random forest (RF) algorithms and support vector machine (SVM). The SVM-recursive feature elimination (SVM-RFE) analysis was performed on the 28 MRGs, the top 10 ideal feature genes were chosen as a result. In parallel, a random forest analysis was carried out, and gene importance scores from the model were used to pick the top 10 genes based on the error rate curve (Figures 3A, B). Nine genes in total were obtained by intersecting the genes chosen by SVM-RFE and random forest: PGAM5, PRKN, TOMM40, FUNDC1, MAP1LC3B, PINK1, MTERF3, CSNK2A2, and SRC. qPCR analysis revealed distinct expression patterns of these genes across the panel of breast cancer cell lines or (MCF-7, MDA-MB-231, MDA-MB-468, Hs-578T, HCC1937, and BT-549) compared to the non-tumorigenic MCF-10A control \u003cstrong\u003e(Supplement file Figure 7)\u003c/strong\u003e. Gene expression levels of PGAM5, PRKN, TOMM40, FUNDC1, MAP1LC3B, PINK1, MTERF3, CSNK2A2, and SRC were assessed across six breast cancer cell lines (MCF-7, MDA-MB-231, MDA-MB-468, Hs-578T, HCC1937, and BT-549) and compared to the normal mammary epithelial cell line MCF-10A. We found that PGAM5, MAP1LC3B, PRKN, PINK1, SRC, FUNDC1, CSNK2A2, and MTERF3 showed significantly higher expression in most of the breast cancer cell lines compared to MCF-10A, particularly in MDA-MB-468, BT-549, and HCC1937 cells. TOMM40 expression was stable across both cancer and normal cell lines, showing no significant difference from MCF-10A. To assess the association between the selected genes and BC occurrence, a multivariable logistic regression analysis was performed, confirming the significance of these nine genes (Figure 3C). An area under the curve (AUC) of 0.99 was shown by the receiver operating characteristic (ROC) curve, which was used to assess the prediction ability of the model. An independent dataset used for internal validation of the model produced an AUC of 0.988(Figures 3D, E). Calibration plots and decision curve analysis (DCA) were employed to assess the calibration and clinical utility of the predictive model (Figures 3F-G). A nomogram model was created to offer a user-friendly tool for displaying the feature genes\u0026apos; predictive performance (Figure 3H). When all nine feature genes were combined into one model, the predictive value of the integrated model was higher than that of the feature genes alone. Then, we gathered the IHC staining of nine proteins linked to feature genes from the HPA database. These proteins were taken from both normal and breast cancer tissue. Our findings are supported by the observation that three of the feature genes (SRC, TOMM40, and PGAM5) had significantly greater levels of protein expression in tumor samples compared to normal samples (Supplementary Figure 4). There was no discernible variation in the levels of protein expression of the other feature genes between the normal and breast cancer samples.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003e3.4 Identification and Differential Analysis of Breast Cancer Subtypes\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eAll tumor samples were categorized into three different subtypes using hierarchical clustering based on the nine mitophagy-related characteristic genes expression matirx (Figure 4A). Significant variations in the expression levels of these characteristic genes were found among the three groups based on differential expression analysis (Figure 4B). The heatmap (Figure 4C) shows the distribution of clinical pathological features across the subtypes that have been discovered, as well as the feature genes\u0026apos; differential expression patterns. A noteworthy finding was that the predictive model performed well in differentiating between the three subtypes of breast cancer (Figure 4D). Gene set variation analysis (GSVA) was used to better describe the molecular variations between the subtypes. The findings showed that, in comparison to the other two clusters, Cluster 1 had a substantially higher overall expression level of mitophagy-related genes (MRGs) (Figure 4E). To evaluate the predictive value of the discovered subgroups, survival analysis was performed. Nevertheless, there were no appreciable variations in survival rates between the three groups (Figure 4F). This data suggests that patient survival outcomes may not be directly impacted by the subgroups that were discovered, despite their different molecular profiles.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003e3.5 Constructing a Gene Co-expression Network\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eWe used a weighted gene co-expression network analysis (WGCNA) to find co-expression modules that correlated best with each subtype of breast cancer. The analysis encompassed all samples, and 5 was determined to be the ideal soft-thresholding value. Thirteen unique co-expression modules were produced as a result of the WGCNA (Figures 5A, B, C). Interestingly, the blue module had the highest association with Cluster 1, while the turquoise module had a correlation with Cluster 2. Furthermore, the black module was discovered to be connected with Cluster 3 (Figure 5D). Using the Metascape database, functional enrichment analysis was performed on the modules that showed the strongest correlation with each cluster in order to obtain insights into the biological significance of the co-expression modules. In Cluster 1, the genes in the blue module were considerably enriched in pathways related to cell cycle biological processes (GO:000278, GO:0010564, hsa04110) (Supplementary Figure 5, top). This data implies that disrupted cell cycle regulation, a hallmark of cancer, may be present in Cluster 1. Cluster 2 showed a high enrichment of genes within the turquoise module in pathways associated with the synthesis of cilia and axonemal dynein complexes (Supplementary Figure 5, middle). Crucial organelles, cilia are engaged in a number of biological functions, including as movement and cell signaling. The predominance of cilium-related pathways in Cluster 2 implies potential abnormalities in ciliary function, which may contribute to the unique biological properties of this subtype. Lastly, In Cluster 3, the black module\u0026apos;s genes were considerably enriched in the regulation of small GTPase-mediated signal transduction (Supplementary Figure 5, bottom). A number of essential cellular functions, including cell division, migration, and survival, are regulated by small GTPases like Ras and Rho. The pathway\u0026apos;s enrichment in Cluster 3 raises the possibility that faulty small GTPase signaling contributes to the pathophysiology of this particular subtype of breast cancer.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003e3.6 Estimating the Immune Response and Drug Sensitivity to Chemotherapy\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eTo study differences in the immune microenvironment between the identified breast cancer (BC) subtypes, we used the CIBERSORT program to estimate the infiltration proportions of various immune and stromal cells. The research found substantial differences in infiltration levels of na\u0026iuml;ve B cells, resting T cells, activated T cells, regulatory T cells, helper T cells, monocytes, macrophages, resting dendritic cells, and resting mast cells amongst clusters. (Supplementary Figure 6A). These findings suggest that each BC subtype may exhibit distinct immune cell compositions, potentially influencing their response to immunotherapy. Tumor Immune Dysfunction and Exclusion (TIDE) algorithm was utilized to forecast the efficacy of immune checkpoint inhibition treatment. Interestingly, compared to other clusters, cluster 2 showed the lowest TIDE score, dysfunction and exclusion, but the highest MSI score, suggesting a possible improved response to immunotherapy (Figures 6A-D). To determine the sensitivity of various BC subtypes to chemotherapeutic medications, we used the \u0026ldquo;oncoPredict\u0026rdquo; program to compute half-maximal inhibitory concentration (IC50) values for many routinely used medicines. The findings showed that, in comparison to the other two clusters, cluster 1 had lower IC50 values for methotrexate, paclitaxel, gemcitabine, gemcitabine with carboplatin, and doxorubicin (Figures 6E-J). This finding suggests that patients in cluster 1 might benefit more from these specific chemotherapy regimens but not benefit from cyclophosphamide than other clusters. Conversely, cluster 2 exhibited significantly higher IC50 values in doxorubicin, paclitaxel, gemcitabine, gemcitabine plus carboplatin, and methotrexate, indicating that patients in cluster 2 might not benefits from these regents than other cluster patients. Given the effectiveness of immune checkpoint inhibitors in a subset of BC patients, we looked into the differences in immunotherapy checkpoint gene expression among the identified BC subtypes. Significant variations in the expression levels of several important checkpoint genes, such as BTNL9, BTN2A1, TNFSF9, PDCD1, CD47, PVR, TNFRSF4, CD274, and CTLA4, were found by our research (Supplementary Figure 6B). These findings suggest that the BC subtypes may exhibit varying responses to immune checkpoint blockade therapies, highlighting the potential for personalized immunotherapy approaches based on the molecular characteristics of each subtype.\u003c/p\u003e"},{"header":"Discussion","content":"\u003cp\u003eThe intricate relationship between mitochondrial dynamics and cancer progression has emerged as a pivotal area of research, particularly in breast cancer (BC), where metabolic reprogramming and immune evasion are critical determinants of tumor aggressiveness. Our comprehensive investigation into mitophagy-related genes (MRGs) reveals their profound influence on BC heterogeneity, offering novel insights into diagnostic stratification, immune modulation, and therapeutic vulnerabilities. By integrating bulk transcriptomics, single-cell resolution analyses, and machine learning-driven modeling, we delineate a molecular framework that positions MRGs as central regulators of BC biology.\u0026nbsp;\u003c/p\u003e\n\u003cp\u003eThe dysregulation of MRGs in BC tissues underscores their dual role in tumor survival and metabolic adaptation. Our data demonstrated significant upregulation of seven MRGs, including FUNDC1, SRC, and TOMM40, alongside the suppression of 18 others such as PINK1 and PRKN. FUNDC1 plays a critical role in tumor biology through its regulation of mitochondrial dynamics, particularly in modulating pathways related to cancer cell survival and proliferation. Notably, elevated FUNDC1 expression has been closely associated with poor clinical outcomes and chemoresistance in endometrial carcinoma.[33, 34] Conversely, the downregulation of PINK1-PRKN-mediated mitophagy, a pathway essential for eliminating dysfunctional mitochondria, may lead to ROS accumulation and genomic instability, hallmarks of aggressive BC phenotypes[35-37]. The elevated expression of SRC and TOMM40 further highlights their dual functionality: SRC kinase not only drives oncogenic signaling but also modulates mitochondrial fission-fusion dynamics[38-41], while TOMM40, a mitochondrial translocase, supports the import of nuclear-encoded proteins necessary for oxidative phosphorylation[42]. These findings align with prior studies linking mitochondrial dysfunction to chemoresistance and metastasis, yet our work uniquely identifies MRGs as integrative biomarkers that bridge metabolic and genomic instability in BC.\u003c/p\u003e\n\u003cp\u003eA striking observation was the discordance between mRNA and protein levels for several MRGs, such as MAP1LC3B and ATG5. This discrepancy likely reflects post-transcriptional regulatory mechanisms, including miRNA-mediated repression or proteasomal degradation. For instance, miR-30a, known to target MAP1LC3B in hepatocellular carcinoma and diabetic cataract[43, 44], may similarly suppress its translation in BC[45], underscoring the necessity of multi-omics approaches to capture the full spectrum of gene regulation. Such insights emphasize that transcriptomic signatures alone may insufficiently predict functional protein activity, necessitating complementary proteomic and spatial analyses in future studies.\u003c/p\u003e\n\u003cp\u003eThe tumor microenvironment (TME) of BC exhibited a complex immune landscape characterized by infiltrating activated T cells, macrophages, and regulatory T cells (Tregs)[46-49]. This immune profile correlated with MRG expression patterns, suggesting a bidirectional interplay between mitophagy and immune evasion. Mitochondrial-derived vesicles (MDVs), which transport mitochondrial antigens to MHC-I molecules, could theoretically activate CD8+ T cells[50]; however, excessive mitophagy in tumor cells may deplete these antigens, thereby evading immune recognition—a phenomenon previously documented in melanoma. At the single-cell level, natural killer (NK) cells displayed high MRG scores coupled with enriched PI3K-AKT-mTOR signaling. Given that mTOR signaling governs NK cell metabolism and cytotoxicity, sustained mitophagy in these cells might paradoxically impair their anti-tumor function through metabolic exhaustion. This hypothesis is supported by the low TIDE dysfunction scores in Cluster 2, which exhibited preserved NK cell activity and heightened responsiveness to immune checkpoint inhibitors (ICIs). These findings position mitophagy as a double-edged sword: while promoting tumor survival, it may also sculpt immune cell fitness, presenting opportunities for therapeutic intervention.\u003c/p\u003e\n\u003cp\u003eThe development of a nine-MRG diagnostic model achieving an AUC of 0.99 represents a paradigm shift in BC stratification. Traditional biomarkers such as ER, PR, Her2, and Ki67, while invaluable for guiding endocrine and targeted therapies, fail to capture the metabolic heterogeneity underlying therapeutic resistance[51-54]. For instance, MRG-high tumors (Cluster 1) demonstrated heightened sensitivity to paclitaxel and gemcitabine, likely due to their reliance on mitochondrial metabolism—a feature invisible to conventional IHC profiling. The model’s strength lies in its integration of multi-gene interactions, which capture pathway-level dysregulation undetectable by single-marker assays. This approach mirrors the success of multi-gene panels like Oncotype DX, which outperform individual biomarkers in predicting chemotherapy benefit. However, transitioning from transcriptomic data to clinical utility requires addressing challenges such as platform-specific batch effects and the integration of MRG signatures with established molecular classifiers like PAM50. Preliminary analyses suggest overlap between MRG-high clusters and basal-like subtypes, known for their metabolic dependencies, highlighting the potential for unified molecular taxonomies that enhance prognostic accuracy.\u003c/p\u003e\n\u003cp\u003eHierarchical clustering identified three BC subtypes with distinct molecular and functional profiles. Cluster 1, characterized by elevated MRG expression and cell cycle activation, mirrors the\u0026nbsp;“basal-like”\u0026nbsp;subtype, renowned for genomic instability and chemosensitivity. The low IC50 values for taxanes and antimetabolites in this cluster align with their mechanism of targeting rapidly dividing cells. Cluster 2, enriched in cilium assembly pathways and exhibiting low TIDE scores, may represent an immune-inflamed phenotype responsive to ICIs. The reappearance of cilia in carcinomas, often linked to Hedgehog pathway activation, could recruit immunosuppressive macrophages, yet the high PD-L1 expression and microsatellite instability (MSI) in this cluster suggest a unique susceptibility to immunotherapy. Cluster 3, dominated by dysregulated small GTPase signaling, likely drives metastatic dissemination through cytoskeletal remodeling and invadopodia formation. Therapeutic targeting of Rho/ROCK inhibitors in this subset could disrupt invasion-proliferation crosstalk, offering a novel strategy to mitigate metastasis.\u003c/p\u003e\n\u003cp\u003eDespite these advancements, our study has limitations. The retrospective nature of TCGA data introduces selection bias, necessitating prospective validation in cohorts with annotated treatment responses. Functional experiments, such as CRISPR-based gene editing or patient-derived organoid models, are essential to establish causality between MRGs and observed phenotypes. Additionally, the limited scope of single-cell RNA sequencing—restricted to two tumor samples—calls for expansion to diverse BC subtypes and metastatic lesions. Spatial multi-omics technologies could further elucidate the spatial distribution of MRGs within tumor niches, revealing interactions between cancer cells and immune infiltrates.\u003c/p\u003e"},{"header":"Conclusions","content":"\u003cp\u003eOur study repositions mitophagy as a central axis in BC biology, influencing metabolic plasticity, immune evasion, and therapeutic response. The MRG-based diagnostic model and subtype classification system transcend conventional paradigms, offering a roadmap for personalized therapy. Future research must prioritize translational validation, integrating these insights into clinical trials to assess their impact on patient outcomes. By aligning mitophagy inhibition with chemotherapy or immunotherapy, clinicians could exploit metabolic vulnerabilities while counteracting immune suppression, ultimately advancing precision oncology in breast cancer.\u003c/p\u003e"},{"header":"Abbreviations","content":" \u003cdiv class=\"gridtable\"\u003e\u003ctable float=\"No\" id=\"Taba\" border=\"1\"\u003e \u003ccolgroup cols=\"2\"\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c1\" colnum=\"1\"\u003e\u003c/div\u003e \u003cdiv align=\"left\" class=\"colspec\" colname=\"c2\" colnum=\"2\"\u003e\u003c/div\u003e \u003cthead\u003e \u003ctr\u003e \u003cth align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eAUC\u003c/div\u003e \u003c/th\u003e \u003cth align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eArea under the curve\u003c/div\u003e \u003c/th\u003e \u003c/tr\u003e \u003c/thead\u003e \u003ctbody\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eBC\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eBreast cancer\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eCCA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eCanonical correlation analysis\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eCTL\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eCytotoxic T lymphocyte\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eDCA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eDecision curve analysis\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eDEGs\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eDifferentially expressed genes\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eGDSC\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eGenomics of Drug Sensitivity in Cancer\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eGEO\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eGene Expression Omnibus database\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eGSEA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eGene Set Enrichment Analysis\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eGSVA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eGene Set Variation Analysis\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eHPA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eHuman Protein Atlas\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eIC50\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eHalf-maximal inhibitory concentration\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eICB\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eImmune checkpoint blockade\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eMAD\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eMedian absolute deviation\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eMes\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eModule eigengenes\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eMRGs\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eMitophagy-related genes\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eMSI\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eMicrosatellite instability\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003emtDNA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eMitochondrial DNA\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eOXPHOS\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eOxidative phosphorylation\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003ePCA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003ePrincipal component analysis\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eqPRC\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eQuantitative Real-Time PCR\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eROC\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eReceiver operating characteristic\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003escRNA-seq\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eSingle-cell RNA sequencing\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eSVM\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eSupport vector machines\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eSVM-RFE\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eSVM-recursive feature elimination\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eTCGA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eThe Cancer Genome Atlas\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eTIDE\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eTumor Immune Dysfunction and Exclusion\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eUMAP\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eUniform Manifold Approximation and Projection\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003ctr\u003e \u003ctd align=\"left\" colname=\"c1\"\u003e \u003cdiv class=\"SimplePara\"\u003eWGCNA\u003c/div\u003e \u003c/td\u003e \u003ctd align=\"left\" colname=\"c2\"\u003e \u003cdiv class=\"SimplePara\"\u003eWeighted gene co-expression network analysis\u003c/div\u003e \u003c/td\u003e \u003c/tr\u003e \u003c/tbody\u003e \u003c/colgroup\u003e \u003c/table\u003e\u003c/div\u003e "},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eEthics approval and consent to participate\u003c/strong\u003e\u003cstrong\u003e：\u003c/strong\u003eNot applicable\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eConsent for publication\u003c/strong\u003e\u003cstrong\u003e：\u003c/strong\u003eNot applicable\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAvailability of data and materials\u003c/strong\u003e\u003cstrong\u003e：\u003c/strong\u003eNot applicable\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCompeting interests\u003c/strong\u003e：The authors declare that they have no competing interests.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eFunding\u003c/strong\u003e\u003cstrong\u003e：\u003c/strong\u003eNot applicable\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAuthors\u0026apos; contributions\u003c/strong\u003e\u003cstrong\u003e：\u003c/strong\u003eWenbin Guo conceived, analyzed, and wrote the original manuscript; Zhiqiang Ye carried out the qPCR and cell culture; and Dan Wu re-edited it. All authors read and approved the final manuscript.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAcknowledgments\u003c/strong\u003e\u003cstrong\u003e：\u003c/strong\u003eNot applicable.\u0026nbsp;\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eChen K, Lu P, Beeraka NM, Sukocheva OA, Madhunapantula SV, Liu J, Sinelnikov MY, Nikolenko VN, Bulygin KV, Mikhaleva LM et al (2022) Mitochondrial mutations and mitoepigenetics: Focus on regulation of oxidative stress-induced responses in breast cancers. Semin Cancer Biol 83:556\u0026ndash;569\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGautam S, Marwaha D, Singh N, Rai N, Sharma M, Tiwari P, Urandur S, Shukla RP, Banala VT, Mishra PR (2023) Self-Assembled Redox-Sensitive Polymeric Nanostructures Facilitate the Intracellular Delivery of Paclitaxel for Improved Breast Cancer Therapy. Mol Pharm 20(4):1914\u0026ndash;1932\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYu H, Yang C, Jian L, Guo S, Chen R, Li K, Qu F, Tao K, Fu Y, Luo F et al (2019) Sulfasalazine\u0026ndash;induced ferroptosis in breast cancer cells is reduced by the inhibitory effect of estrogen receptor on the transferrin receptor. Oncol Rep 42(2):826\u0026ndash;838\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCui L, Gouw AM, LaGory EL, Guo S, Attarwala N, Tang Y, Qi J, Chen YS, Gao Z, Casey KM et al (2021) Mitochondrial copper depletion suppresses triple-negative breast cancer in mice. Nat Biotechnol 39(3):357\u0026ndash;367\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRamchandani D, Berisa M, Tavarez DA, Li Z, Miele M, Bai Y, Lee SB, Ban Y, Dephoure N, Hendrickson RC et al (2021) Copper depletion modulates mitochondrial oxidative phosphorylation to impair triple negative breast cancer metastasis. Nat Commun 12(1):7311\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSun Y, Shen W, Hu S, Lyu Q, Wang Q, Wei T, Zhu W, Zhang J (2023) METTL3 promotes chemoresistance in small cell lung cancer by inducing mitophagy. J Exp Clin Cancer Res 42(1):65\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang J, Zhu P, Li R, Ren J, Zhou H (2020) Fundc1-dependent mitophagy is obligatory to ischemic preconditioning-conferred renoprotection in ischemic AKI via suppression of Drp1-mediated mitochondrial fission. Redox Biol 30:101415\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHu Z, Ding J, Ma ZC, Sun RP, Seoane JA, Shaffer JS, Suarez CJ, Berghoff AS, Cremolini C, Falcone A et al (2019) Quantitative evidence for early metastatic seeding in colorectal cancer. Nat Genet 51(7):1113\u0026ndash;\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMiao Z, Miao Z, Wang S, Wu H, Xu S (2022) Exposure to imidacloprid induce oxidative stress, mitochondrial dysfunction, inflammation, apoptosis and mitophagy via NF-kappaB/JNK pathway in grass carp hepatocytes. Fish Shellfish Immunol 120:674\u0026ndash;685\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi R, Xin T, Li D, Wang C, Zhu H, Zhou H (2018) Therapeutic effect of Sirtuin 3 on ameliorating nonalcoholic fatty liver disease: The role of the ERK-CREB pathway and Bnip3-mediated mitophagy. Redox Biol 18:229\u0026ndash;243\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePanigrahi DP, Praharaj PP, Bhol CS, Mahapatra KK, Patra S, Behera BP, Mishra SR, Bhutia SK (2020) The emerging, multifaceted role of mitophagy in cancer and cancer therapeutics. Semin Cancer Biol 66:45\u0026ndash;58\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen C, Liu Y, Liu L, Si CH, Xu YX, Wu XK, Wang CZ, Sun ZQ, Kang QZ (2023) Exosomal circTUBGCP4 promotes vascular endothelial cell tipping and colorectal cancer metastasis by activating Akt signaling pathway. J EXPERIMENTAL Clin CANCER Res 42(1)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLin Q, Li S, Jin H, Cai H, Zhu X, Yang Y, Wu J, Qi C, Shao X, Li J et al (2023) Mitophagy alleviates cisplatin-induced renal tubular epithelial cell ferroptosis through ROS/HO-1/GPX4 axis. Int J Biol Sci 19(4):1192\u0026ndash;1210\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAltea-Manzano P, Doglioni G, Liu YW, Cuadros AM, Nolan E, Fern\u0026aacute;ndez-Garc\u0026iacute;a J, Wu Q, Planque M, Laue KJ, Cidre-Aranaz F et al (2023) A palmitate-rich metastatic niche enables metastasis growth via p65 acetylation resulting in pro-metastatic NF-κB signaling. Nat CANCER 4(3):344\u0026ndash;\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLu Y, Li Z, Zhang S, Zhang T, Liu Y, Zhang L (2023) Cellular mitophagy: Mechanism, roles in diseases and small molecule pharmacological regulation. Theranostics 13(2):736\u0026ndash;766\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTan J, Egelston CA, Guo W, Stark JM, Lee PP (2024) STING signalling compensates for low tumour mutation burden to drive anti-tumour immunity. EBioMedicine 101:105035\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHao Y, Stuart T, Kowalski MH, Choudhary S, Hoffman P, Hartman A, Srivastava A, Molla G, Madad S, Fernandez-Granda C et al (2024) Dictionary learning for integrative, multimodal and scalable single-cell analysis. Nat Biotechnol 42(2):293\u0026ndash;304\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHao Y, Hao S, Andersen-Nissen E, Mauck WM 3rd, Zheng S, Butler A, Lee MJ, Wilk AJ, Darby C, Zager M et al (2021) Integrated analysis of multimodal single-cell data. Cell 184(13):3573\u0026ndash;3587 e3529\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiu YM, Ge JY, Chen YF, Liu T, Chen L, Liu CC, Ma D, Chen YY, Cai YW, Xu YY et al (2023) Combined Single-Cell and Spatial Transcriptomics Reveal the Metabolic Evolvement of Breast Cancer during Early Dissemination. Adv Sci (Weinh) 10(6):e2205395\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSubramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, Paulovich A, Pomeroy SL, Golub TR, Lander ES et al (2005) Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A 102(43):15545\u0026ndash;15550\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiberzon A, Birger C, Thorvaldsdottir H, Ghandi M, Mesirov JP, Tamayo P (2015) The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst 1(6):417\u0026ndash;425\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCORTES C (2005) Support-Vector Networks. Mach Learn 20:273\u0026ndash;297\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBREIMAN L: Random Forests. Mach Learn, (2001) 45:5\u0026ndash;32\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMarchevsky AM, Patel S, Wiley KJ, Stephenson MA, Gondo M, Brown RW, Yi ES, Benedict WF, Anton RC, Cagle PT (1998) Artificial neural networks and logistic regression as tools for prediction of survival in patients with Stages I and II non-small cell lung cancer. Mod Pathol 11(7):618\u0026ndash;625\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eVan Calster B, Wynants L, Verbeek JFM, Verbakel JY, Christodoulou E, Vickers AJ, Roobol MJ, Steyerberg EW (2018) Reporting and Interpreting Decision Curve Analysis: A Guide for Investigators. Eur Urol 74(6):796\u0026ndash;804\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eUhlen M, Fagerberg L, Hallstrom BM, Lindskog C, Oksvold P, Mardinoglu A, Sivertsson A, Kampf C, Sjostedt E, Asplund A et al (2015) Proteomics. Tissue-based map of the human proteome. \u003cem\u003eScience\u003c/em\u003e 347(6220):1260419\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSjostedt E, Zhong W, Fagerberg L, Karlsson M, Mitsios N, Adori C, Oksvold P, Edfors F, Limiszewska A, Hikmet F et al (2020) An atlas of the protein-coding genes in the human, pig, and mouse brain. \u003cem\u003eScience\u003c/em\u003e 367(6482)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLangfelder P, Horvath S (2008) WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics 9:559\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, Benner C, Chanda SK (2019) Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun 10(1):1523\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFu J, Li K, Zhang W, Wan C, Zhang J, Jiang P, Liu XS (2020) Large-scale public data reuse to model immunotherapy response and resistance. Genome Med 12(1):21\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eJiang P, Gu S, Pan D, Fu J, Sahu A, Hu X, Li Z, Traugh N, Bu X, Li B et al (2018) Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med 24(10):1550\u0026ndash;1558\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMaeser D, Gruener RF, Huang RS (2021) oncoPredict: an R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data. Brief Bioinform 22(6)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi J, Agarwal E, Bertolini I, Seo JH, Caino MC, Ghosh JC, Kossenkov AV, Liu Q, Tang HY, Goldman AR et al (2020) The mitophagy effector FUNDC1 controls mitochondrial reprogramming and cellular plasticity in cancer cells. \u003cem\u003eSci Signal\u003c/em\u003e 13(642)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHou H, Er P, Cheng J, Chen X, Ding X, Wang Y, Chen X, Yuan Z, Pang Q, Wang P et al (2017) High expression of FUNDC1 predicts poor prognostic outcomes and is a promising target to improve chemoradiotherapy effects in patients with cervical cancer. Cancer Med 6(8):1871\u0026ndash;1881\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eXie Y, Liu J, Kang R, Tang D (2020) Mitophagy Receptors in Tumor Biology. Front Cell Dev Biol 8:594203\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eOshi M, Gandhi S, Yan L, Tokumaru Y, Wu R, Yamada A, Matsuyama R, Endo I, Takabe K (2022) Abundance of reactive oxygen species (ROS) is associated with tumor aggressiveness, immune response, and worse survival in breast cancer. Breast Cancer Res Treat 194(2):231\u0026ndash;241\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGuo M, Wang SM (2021) Genome Instability-Derived Genes Are Novel Prognostic Biomarkers for Triple-Negative Breast Cancer. Front Cell Dev Biol 9:701073\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDjeungoue-Petga MA, Lurette O, Jean S, Hamel-Cote G, Martin-Jimenez R, Bou M, Cannich A, Roy P, Hebert-Chatelain E (2019) Intramitochondrial Src kinase links mitochondrial dysfunctions and aggressiveness of breast cancer cells. Cell Death Dis 10(12):940\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSirvent A, Benistant C, Roche S (2012) Oncogenic signaling by tyrosine kinases of the SRC family in advanced colorectal cancer. Am J Cancer Res 2(4):357\u0026ndash;371\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePelaz SG, Tabernero A (2022) Src: coordinating metabolism in cancer. Oncogene 41(45):4917\u0026ndash;4928\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLu KV, Zhu S, Cvrljevic A, Huang TT, Sarkaria S, Ahkavan D, Dang J, Dinca EB, Plaisier SB, Oderberg I et al (2009) Fyn and SRC are effectors of oncogenic epidermal growth factor receptor signaling in glioblastoma patients. Cancer Res 69(17):6889\u0026ndash;6898\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYang W, Shin HY, Cho H, Chung JY, Lee EJ, Kim JH, Kang ES (2020) TOM40 Inhibits Ovarian Cancer Cell Growth by Modulating Mitochondrial Function Including Intracellular ATP and ROS Levels. Cancers (Basel) 12(5)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFu XT, Shi YH, Zhou J, Peng YF, Liu WR, Shi GM, Gao Q, Wang XY, Song K, Fan J et al (2018) MicroRNA-30a suppresses autophagy-mediated anoikis resistance and metastasis in hepatocellular carcinoma. Cancer Lett 412:108\u0026ndash;117\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang L, Cheng R, Huang Y (2017) MiR-30a inhibits BECN1-mediated autophagy in diabetic cataract. Oncotarget 8(44):77360\u0026ndash;77368\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eXiao B, Shi X, Bai J (2019) miR-30a regulates the proliferation and invasion of breast cancer cells by targeting Snail. Oncol Lett 17(1):406\u0026ndash;413\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiu J, Wang X, Deng Y, Yu X, Wang H, Li Z (2021) Research Progress on the Role of Regulatory T Cell in Tumor Microenvironment in the Treatment of Breast Cancer. Front Oncol 11:766248\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRetecki K, Seweryn M, Graczyk-Jarzynka A, Bajor M (2021) The Immune Landscape of Breast Cancer: Strategies for Overcoming Immunotherapy Resistance. \u003cem\u003eCancers (Basel)\u003c/em\u003e 13(23)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNishikawa H, Koyama S (2021) Mechanisms of regulatory T cell infiltration in tumors: implications for innovative immune precision therapies. J Immunother Cancer 9(7)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHuang P, Zhou X, Zheng M, Yu Y, Jin G, Zhang S (2023) Regulatory T cells are associated with the tumor immune microenvironment and immunotherapy response in triple-negative breast cancer. Front Immunol 14:1263537\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePicca A, Guerra F, Calvani R, Coelho-Junior HJ, Landi F, Bucci C, Marzetti E (2023) Mitochondrial-Derived Vesicles: The Good, the Bad, and the Ugly. Int J Mol Sci 24(18)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHu X, Chen W, Li F, Ren P, Wu H, Zhang C, Gu K (2023) Expression changes of ER, PR, HER2, and Ki-67 in primary and metastatic breast cancer and its clinical significance. Front Oncol 13:1053125\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDavey MG, Hynes SO, Kerin MJ, Miller N, Lowery AJ (2021) Ki-67 as a Prognostic Biomarker in Invasive Breast Cancer. \u003cem\u003eCancers (Basel)\u003c/em\u003e 13(17)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGoncalves AC, Richiardone E, Jorge J, Polonia B, Xavier CPR, Salaroglio IC, Riganti C, Vasconcelos MH, Corbet C, Sarmento-Ribeiro AB (2021) Impact of cancer metabolism on therapy resistance - Clinical implications. Drug Resist Updat 59:100797\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKim J, DeBerardinis RJ (2019) Mechanisms and Implications of Metabolic Heterogeneity in Cancer. Cell Metab 30(3):434\u0026ndash;446\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":true,"hideJournal":true,"highlight":"","institution":"","isAcceptedByJournal":false,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true},"keywords":"Mitophagy, breast cancer, Single-Cell Sequencing, Immune Cells, Immunotherapy Response","lastPublishedDoi":"10.21203/rs.3.rs-6960401/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-6960401/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003ch2\u003eBackground\u003c/h2\u003e \u003cp\u003eThe development of breast cancer (BC) entails intricate immunological and molecular interactions. Although mitophagy-related genes (MRGs) are essential for maintaining cellular homeostasis, little is known about their expression patterns, immunological effects, and therapeutic applications in BC. The purpose of this work was to describe immunological microenvironment changes, MRG dysregulation, and their predictive modeling potential in BC.\u003c/p\u003e\u003ch2\u003eMethods\u003c/h2\u003e \u003cp\u003eUsing data from single-cell RNA sequencing (scRNA-seq) (GSE248288) and The Cancer Genome Atlas (TCGA), we examined immune cell infiltration MRG expression profiles, and functional pathways. Machine learning methods identified significant MRGs for predictive modeling, which were confirmed using qPCR in BC cell lines. Subtypes were identified using hierarchical clustering, and weighted gene co-expression network analysis (WGCNA) showed subtype-specific modules. Drug sensitivity and immunotherapy responses were evaluated using the oncoPredict and Tumor Immune Dysfunction and Exclusion (TIDE) algorithms.\u003c/p\u003e\u003ch2\u003eResults\u003c/h2\u003e \u003cp\u003eTumor tissues showed elevated MRGs (e.g., SRC, PGAM5, FUNDC1) and altered immunological infiltration, with more activated T cells and macrophages but fewer na\u0026iuml;ve B cells. scRNA-seq demonstrated increased MRG activity in NK and B cells, which is connected to PI3K-AKT-mTOR signaling and allograft rejection. A nine-MRG predictive model (PGAM5, PRKN, TOMM40, FUNDC1, MAP1LC3B, PINK1, MTERF3, CSNK2A2, and SRC) demonstrated great diagnostic accuracy (AUC: 0.99). BC subgroups based on MRGs expression displayed unique molecular profiles: Cluster 1 (cell cycle dysregulation), Cluster 2 (cilia-related pathways), and Cluster 3 (small GTPase signaling). Cluster 2 showed potential immunotherapy responsiveness (low TIDE scores), whereas Cluster 1 was sensitive to chemotherapy (paclitaxel, gemcitabine).\u003c/p\u003e\u003ch2\u003eConclusions\u003c/h2\u003e \u003cp\u003eMRGs are important regulators of breast cancer progression, regulating immunological dynamics, metabolic pathways, and treatment responses. The established nine-MRG model accurately predicts BC occurrence, whereas subtype-specific biochemical and immunological aspects provide insights for tailored therapy. These findings emphasize MRGs as possible biomarkers for diagnosis and personalized therapy options in BC.\u003c/p\u003e","manuscriptTitle":"Mitophagy-Driven Immune Cell Infiltration Patterns Define Breast Cancer Subtypes with Differential Treatment Responses","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2025-06-25 01:47:11","doi":"10.21203/rs.3.rs-6960401/v1","editorialEvents":[{"type":"communityComments","content":0}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"c1337ef2-53f3-4f42-8b4f-ab8ab0933dca","owner":[],"postedDate":"June 25th, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"posted","subjectAreas":[],"tags":[],"updatedAt":"2025-06-25T12:23:51+00:00","versionOfRecord":[],"versionCreatedAt":"2025-06-25 01:47:11","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-6960401","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-6960401","identity":"rs-6960401","version":["v1"]},"buildId":"8U1c8b4HqxoKbykW_rLl7","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}