
==== Front
Heliyon
Heliyon
Heliyon
2405-8440
Elsevier

S2405-8440(24)12476-0
10.1016/j.heliyon.2024.e36445
e36445
Research Article
Copper metabolism–related signature for prognosis prediction and MMP13 served as malignant factor for breast cancer
Han Chaojie HanCJ163@163.com
abc
Feng Zhangyang youngzhangfeng@126.com
a
Wang Yingjian 3194075@zju.edu.cn
b
Hu Mengsi 3415110@zju.edu.cn
b
Xu Shoufang xusf@zju.edu.cn
b
Jiang Feiyu 3320096@zju.edu.cn
b
Han Yetao 3414299@zju.edu.cn
b
Liu Zhiwei 3191038@zju.edu.cn
b⁎
Li Yunsen yunsenli@suda.edu.cn
a⁎⁎
a Institutes of Biology and Medical Sciences, Soochow University, 333 East Ganjiang Road, Suzhou, Jiangsu, 215127, China
b Department of Blood Transfusion, Sir Run Run Shaw Hospital, Zhejiang University School of Medicine, 3 East Qingchun Road, Hangzhou, Zhejiang, 310016, China
c Zhejiang Zhenyuan Biotech Co., LTD, 61 Yuedongbei Road, Shaoxing, Zhejiang, 312000, China
⁎ Corresponding author. Department of Blood Transfusion Sir Run Run Shaw Hospital, Zhejiang University School of Medicine, 3 East Qingchun Road, Hangzhou, Zhejiang, 310016, China. 3191038@zju.edu.cn
⁎⁎ Corresponding author. Institutes of Biology and Medical Sciences Soochow University, 333 East Ganjiang Road, Suzhou, Jiangsu, 215127, China. yunsenli@suda.edu.cn
17 8 2024
30 9 2024
17 8 2024
10 18 e3644518 6 2024
13 8 2024
15 8 2024
© 2024 The Authors
2024
https://creativecommons.org/licenses/by-nc-nd/4.0/ This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
Objectives

To comprehensively analyze the copper metabolism in Breast cancer, we established a prognostic signature for breast cancer (BC) related to copper metabolism.

Methods

Copper metabolism-related genes were sourced from previous literatures and were selected by the Univariate Cox regression. Cu-enrichment scores were calculated via ssGSEA. Differentially expressed genes were identified with limma between high and low Cu-enrichment scores group, then we used the Random Survival Forest and LASSO to build the CuScore for BC. Kaplan-Meier analysis, ROC curves, and Cox regression were used to evaluate CuScore. Genomic mutations were analyzed with GISTIC. Immune cells were examined using ESTIMATE, ssGSEA and TIMER. Enrichment analysis used clusterProfiler and GSVA. The GDSC database and oncoPredict package analyzed chemotherapeutic sensitivity. MMP13 was selected for in vitro assays.

Results

Four copper metabolism-related genes (UBE2D2, SLC31A1, ATP7A, and MAPK1) with prognostic value were identified. Higher expression levels of these genes were associated with higher Cu-enrichment scores, a factor of malignancy in breast cancer. Among 115 differentially expressed genes, 19 prognostic genes were identified, with three (CEACAM5, MMP13, and CRISP3) highlighted by Random Survival Forest and LASSO. Higher CuScores correlated with worse prognoses and were effective in predicting breast cancer outcomes. CuScore and metastasis were independent prognostic factors. Tumor-infiltrating immune cells were associated with lower CuScores. GO-GSEA analysis indicated six immune-related pathways might be regulated by CuScore. Patients with higher CuScores had lower TMB and were more sensitive to Sapitinib and LCL161, while those with lower CuScores might respond better to anti-PD1 therapy. High MMP13 expression in breast cancer was linked to malignancy, affecting cell proliferation and migration.

Conclusion

The identified copper metabolism-related gene signature has the potential to predict prognosis and guide clinical treatment for BC. Among these genes, MMP13 may act as a malignant factor in BC.

Keywords

Copper metabolism
Prognosis
Breast cancer
Treatment
MMP13
==== Body
pmc1 Introduction

Breast cancer (BC), the most prevalent malignancy in women globally and now the most common cancer worldwide, accounted for 24.5 % of new cancer cases and 15.5 % of cancer-related deaths among women [1], with 2.3 million new cases recorded and projections suggesting it will exceed 3 million cases by 2040 [2].The overall death rate from BC relative to its incidence is declining, primarily due to advancements in prevention, diagnosis, and treatment. Continued efforts in these areas could further reduce BC mortality [3]. Therefore, analyzing and evaluating prognostic factors is a crucial research direction, significantly contributing to the development of precision medicine for BC, and holds high clinical and social value.

Clinically, BC is categorized into four primary subtypes, determined by its histopathological features and the presence of receptors, such as estrogen receptor (ER), progesterone receptor (PR), and human epidermal growth factor receptor 2 (HER2) in the tumor cells [4].In recent years, with the gradual discovery of genomic markers (BRCA1, BRCA2, PIK3CA) and immune markers (tumor-infiltrating lymphocytes and PDL1), the molecular profile of BC has been further expanded and used to assess BC survival and prognosis [3]. In addition to BRCA1 or BRCA2, several tumor driver genes (including TP53, PALB2, PIK3CA, ESR1, CDH1, PTEN, ATM, CHEK2, RAD51C, and BARS1) are mutated in BC [[5], [6], [7], [8]]. In addition, a variety of non-genetic factors including age, breast lesions, high mammography, chest radiation, body mass index (BMI), exogenous hormone use, alcohol, menarche, and menopause also influence BC progression.

While traditional factors like tumor grade, size, and lymph node status offer valuable prognostic insights for newly diagnosed BC patients, medical experts generally agree that these prognostic factors no longer meet optimal patient management needs in the era of personalized cancer treatment [9]. As BC treatment advances, the need to pinpoint reliable biomarkers that can forecast patient outcomes and guide treatment decisions grows increasingly critical.

Bioinformatics is an interdisciplinary subject arising from the development of the human genome project, which mainly involves biology, computer, and mathematics [10]. The development of high-throughput data sequencing facilitates deeper genomic analysis [10]. Large amounts of medical data help us to acquire a more comprehensive knowledge of biology. Obtaining useful information in a large number of biological data has become a problem that needs to be solved, and data mining has become an important tool in bioinformatics research. Data mining is based on the existing data to obtain useful information and discover and reveal the biological law [10,11].

Copper is a transition metal with redox activity, capable of cycling between its reduced Cu+ and oxidized Cu2+ forms under standard chemical and physiological conditions [12]. Copper ions are involved in various biochemical reactions, acting as electron donors or acceptors [13]. They are capable of binding to various proteins and enzymes, functioning as cofactors or structural elements, and are essential in controlling energy metabolism, mitochondrial respiration, and antioxidant protection, among various other physiological functions [14,15]. Maintaining a dynamic balance of copper ions is essential, as imbalances can lead to oxidative stress and autophagy dysfunction [16], potentially triggering copper-related diseases.

In this study, we combined Kaplan-Meier survival analysis, Cox regression, Random Survival Forest(SVM), and Lasso algorithm to identify copper metabolism-related genes that affect BC prognosis from public cohorts. A copper-related prognostic assessment model was established and validated across multiple cohorts, demonstrating predictive power. Lastly, we evaluated the therapeutic properties of the prognostic model.

2 Materials and methods

2.1 Data collecting and data preprocessing for BC

Genomic information and complete clinicopathological annotations for BC were available from three database, including The Cancer Genome Atlas (TCGA, n = 1084), Molecular Taxonomy of Breast Cancer International Consortium (METABRIC, n = 1904) and Gene Expression Omnibus (GEO, GSE41119 n = 156, GSE103091 n = 107) were collected. Finally, a total of four BC transcriptional profile data cohorts were collected in this study, and patients with insufficient OS information were excluded. A total of 36 copper metabolism–related genes were obtained from previous literature [17].

2.2 Establishment of Cu-enrichment score in TCGA-BRCA

The Cox regression was used to facilitate dimensionality reduction and identify the significant prognostic value of 36 copper metabolism genes. Subsequently, copper enrichment scores for the prognostic genes were calculated using single sample gene set enrichment analysis (ssGSEA). According to the optimal cut-off, BC patients were then classified into high and low score groups.

2.3 Establishment of copper metabolism–related signature in TCGA-BRCA

To establish the signature, we firstly identified the genes that were differentially expressed between two high and low Cu-enrichment scores groups, employing criteria of logFC >1 and P-value< 0.05 by using limma package. Following this, the cox regression analysis was applied to pinpoint prognostic genes from the DEGs related to copper metabolism. Subsequently, the Random Survival Forest model, with a VIP threshold of greater than 0.3, was utilized to refine the selection of valuable prognostic genes. Finally, the LASSO method was employed to obtain the coefficients of the prognostic genes, thus establishing a signature related to copper metabolism [18], then we obtained the CuScore in the cohorts.

2.4 Efficacy evaluation for the CuScore of BC

For the TCGA cohort, we calculated the CuScore for each patients and classified them into two groups, according to the optimal prognostic cutoff of CuScore. To evaluate the correlation between CuScore and overall survival (OS), we utilized KM curves and Cox regression analyses. Additionally, Time-ROC analysis was performed to determine the prognostic accuracy and efficacy of CuScore. Additionally, to confirm the prognostic independence of CuScore in BC, Cox regression was performed using CuScore and various clinicopathological factors.

2.5 Analysis of genomic mutation for CuScore

Acquiring from TCGA database, the copy number variation and somatic mutations of patients were obtained for further analysis. We carried out evaluation of copy number of CNV landscape and amplification or delete peak increase or loss through GISTIC 2.0 [19].

2.6 Characteristic of tumor infiltrating immune cells for BC

For each BC patient, we used the ESTIMATE [20] and TIMER2.0 [21] web to calculate the ESTIMATE score and immune infiltration levels in BC. Additionally, we used the through Gene set variation analysis (GSVA) package to calculate the expression of 28 immune cells by the method of ssGSEA [22]. Additionally, seven categories of immunomodulators were identified from existing literature to explore the association between immune processes and CuScore.

2.7 Enrichment analysis and chemotherapeutics sensitivity analysis for CuScore

Gene Ontology (GO) functions analysis for CuScore were conducted by the clusterProfiler package. Kyoto Encyclopedia of Genes Genomes (KEGG) pathway analysis also performed. From the MSigDB, we also downloaded some target gene sets and calculate their abundance by using GSEA to show the function of CuScore [23]. For drug sensitivity of BC patients, we firstly extracted the pharmacogenomics from the database of GDSC [24]. Then, the drug responses were calculated using the drug sensitivity calculated by the oncoPredict package of R software [25].

2.8 Cell culture

Cell lines BT-474 and MCF-7, as well as normal human breast cells MCF-10A, were sourced from the Cell Bank of the Type Culture Collection of the Chinese Academy of Sciences. MCF-7 was transduced with pLKO.1-shRNA and subsequently selected using puromycin (4 μg/mL) over a 5-day period. All cell lines were cultured in DMEM (D5796) supplemented with 10 % fetal bovine serum (FBS, Gibco, 306.00301), 1 % amphotericin B (A2942), and 1 % penicillin-streptomycin (Sigma, P4333). Cells were incubated at 37 °C in a humidified atmosphere with 5 % CO2 and were passaged or provided with fresh medium every other day.

2.9 RNA extraction and RT-PCR

Total RNA was extracted from cultured BC cells using the RNeasy Plus kit (Qiagen, Germany) according to the manufacturer's guidelines. In brief, monolayer cultures were trypsinized, washed twice with PBS, and subjected to RNA extraction. The integrity and concentration of the extracted RNA were measured using a NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific, USA) prior to RNA sequencing. For RT-PCR analysis, total mRNA was isolated using the TRIzol reagent (Invitrogen). HiScript III All-in-one RT SuperMix was used to convert One microgram of mRNA into cDNA utilizing for qPCR, with PCR amplification quantified by SYBR Green (YEASEN Biotech). The primers used for RT-PCR are listed in Table S1.

2.10 Immunofluorescence staining

Cells (2.5 × 10^5) under optimal conditions were fixed with 4 % paraformaldehyde at 25 °C for 30 min and then washed with PBS. Permeabilization was carried out using 0.3 % Triton X-100 at 25 °C for 10 min, followed by a 30-min blocking step with 3 % bovine serum albumin at 37 °C. The cells were subsequently incubated with primary antibodies at 37 °C for 1 h. After washing, secondary antibodies were applied and incubated at 37 °C for another hour. Nuclei were stained with 1 μg/ml DAPI, and coverslips were finally mounted on glass slides for observation under a Carl Zeiss Axio ZI fluorescence microscope.

2.11 Western blot

Lysates were prepared from MCF-10A, BT-474, MCF7, MCF7-sh1, and MCF7-sh2 cells.Then we used the BCA protein assay kit (Beyotime, China) to measure the Protein concentrations. The primary antibodies used included MMP13 (ab39012) and GAPDH (ab8245).

2.12 Cell proliferation assay

At 37 °C, cells (3 × 103 per well) were seeded into 96-well plates, and 5 parallel plates were given to each group. Each well was injected with 10 μl from the Cell Counting Kit-8 (Vazyme, China) solution and incubated for 2 h at 37 °C. The absorbance of each cells plates was tested at 450 nm. Cell proliferation curves were constructed by absorbance at 24 h, 48 h, 72h and 96h.

2.13 Wound-healing migration and Matrigel invasion assays

MCF7-Con, MCF7-sh1, and MCF7-sh2 cells were seeded at a density of 10^6 cells well in six-well plates for wound-healing migration assays. After overnight serum starvation, a linear wound was created in the cell monolayer, which had reached 70%–80 % confluence, using a sterile 200 μl pipette tip. Wound closure was monitored at 0 and 36 h in serum-free DMEM. For the invasion assay, cells (1 × 10^5) in serum-free DMEM were seeded on a polycarbonate filter membrane without 8 μM polyvinylpyrrolidone, using filters precoated with Matrigel diluted 1:11 (BD Bioscience, San Jose, CA). The lower chamber contained DMEM with 10 % FBS, and the cells were incubated for 24 h. Invading cells were then quantified after crystal violet staining.

2.14 Statistical analysis

We used R software (version 4.3) to conduct the statistical analyses in the present study. Specifically, to compare data between two groups, student t-test was used for normally data and Wilcoxon test for non-normally data, and KM curves, estimating overall survival (OS) differences between groups, were created using the R survminer package. Cox regression analysis and Time-ROC curves, which assessed the predictive accuracy of the models, were generated using the R survival and timeROC package. Additionally, heatmaps were produced using the R pheatmap package. Unless otherwise specified, statistical significance was set at a p-value of 0.05.

3 Results

3.1 Characteristics of Cu-enrichment score in TCGA-BRCA cohort

A total of 36 copper metabolism-related genes were identified from the previous research, and through Univariate Cox regression analysis on these genes, four were found to have prognostic significance (Fig. 1A). All the four prognostic genes (UBE2D2, SLC31A1, ATP7A, and MAPK1) were pro-cancer factors because they all have an HR greater than 1 (Fig. 1A). The internal correlation of the four prognostic related copper metabolism genes in TCGA-BRCA was positive (Fig. 1B), and the correlation coefficients of ATP7A and MAPK1 were the highest (0.66). From the Kaplan-Meier curves (Fig. 1C), we could observe the significant difference in outcome between BC patients with low Cu-enrichment score or high Cu-enrichment score, which indicating the factor of malignancy for Cu-enrichment score. Heatmap showed the distribution of the Cu-enrichment score of TCGA-BRCA patients, along with the expression values of the four prognostic genes, stage and age (Fig. 1D). The higher expression value of SLC31A1, UBE2D2, ATP7A, and MAPK1 were associated with higher Cu-enrichment score (Fig. 1D).Fig. 1 Features of Cu-enrichment score in the TCGA cohort. (A) Forest plot of the prognostic significance of genes UBE2D2, SLC31A1, ATP7A, and MAPK1. (B) Shows the correlation matrix of these four genes linked to copper metabolism within the TCGA-BRCA dataset. (C) KM curves of Cu-enrichment scores. (D) Demonstrates the relationships between Cu-enrichment score, clinical parameters, and expression levels of the identified genes, with green depicting low expression and purple high expression.

Fig. 1

3.2 Establishment of copper metabolism–related signature in TCGA-BRCA

Between BC patients with low Cu-enrichment score and BC patients high Cu-enrichment score, we identified 115 DEGs, which also were reflected in volcano plot (Fig. 2A). Among the 115 DEGs, 19 genes with prognostic value were identified (Fig. 2B), including 13 malignant genes with HR > 1 (SEC24A, MTDH, DOP1B, AFF4, OPA1, DHFRP1, LINC01614, AL591806.4, MMP13, CRISP3, MMP1, CLEC3A, and CEACAM5) and six protective genes with HR < 1 (KRT14, SAA1, KRT15, AC020656.1, MIA, and AC005921.2). Therewith, we conducted Random Survival Forest algorithm on the 19 prognostic genes and showed the distribution of error rates generated in Fig. 2C, thus generating three genes (CEACAM5, MMP13, and CRISP3) with more prognostic potential (Fig. 2C). Based on the prognostic genes screened out by Random Survival Forest algorithm, We used the method of LASSO to calculate the gene coefficients with prognostic values of BC (Fig. 2D), thus calculating the CuScore of each patient with BC.Fig. 2 Establishment of copper metabolism–related signature in TCGA-BRCA. (A) Volcano plot of DEGs between high and low Cu-enrichment score groups. (B) Forest plot of 19 prognostic TRP-related genes. (C) Error rate distribution in a Random Survival Forest model (D) Lambda selection of Lasso analysis for highlighting optimal regularization.

Fig. 2

3.3 Efficacy evaluation for the CuScore of BC

After establishing the copper metabolism–related signature using CuScore, we assessed its prognostic value for BC. Prognosis was compared between patients with high and low CuScore, revealing that those with a lower CuScore had a more favorable prognosis, while a higher CuScore was associated with poorer outcomes, as shown by KM curves (Fig. 3A). The AUC values of CuScore for 3-year (0.668), 5-year (0.688), and 10-year (0.819) were all above 0.66, indicating the CuScore's moderate predictive accuracy for BC prognosis (Fig. 3B). Additionally, we identified CuScore (P = 0.001399) and metastasis (P = 0.006842) as the independent prognostic factors of BC (Fig. 3C). To further validate these findings, we analyzed four independent external datasets (METABRIC: n = 1904; GSE41119: n = 156; GSE103091: n = 107), all of which consistently showed worse clinical outcomes for patients with higher CuScore, as demonstrated in METABRIC (Fig. 3D), GSE41119 (Fig. 3E), and GSE103091 (Fig. 3F).Fig. 3 Efficacy evaluation of the CuScore in BC. (A) Kaplan-Meier (KM) survival curve comparing high and low CuScore subgroups in the TCGA-BRCA cohort. (B) Time-dependent ROC curves predicting 3-year, 5-year, and 10-year survival based on CuScore. (C) Forest plot assessing the prognostic value of CuScore alongside clinicopathologic features. (D) KM survival curve comparing high and low CuScore subgroups in the METABRIC cohort. (E) KM survival curve comparing high and low CuScore subgroups in the GSE41119 dataset. (F) KM survival curve comparing high and low CuScore subgroups in the GSE103091 dataset.

Fig. 3

3.4 Characteristic of tumor infiltrating immune cells for BC

Using several immune cell infiltration algorithms, a heatmap was generated to illustrate the infiltrating immune cell populations at two CuScore levels (Fig. 4A). Generally, lower CuScore was associated with higher levels of various immune cells, including CD8 T cells, T cells and several activated immune cell subsets. Conversely, higher CuScore was linked to increased levels of non-immune cells, such as neutrophils, fibroblasts and CD56 dim natural killer cells et al. Additionally, the Immune and ESTIMATE Score were higher in patients with lower CuScore, while Tumor Purity increased with higher CuScore (Fig. 4B). GO-GSEA analysis indicated that six immune-related pathways might be influenced by CuScore, including the receptor signaling pathway of T cell, B cell, interferon-gamma-mediated, chemokine-mediated, immune response, and NK cell activation (Fig. 4C).Fig. 4 Correlation of CuScore with immune cells for BC. (A) The heatmap of infiltrating immune cell at two levels of CuScore. (B) Violin diagrams of the expression differences between high- and low-CuScore on stromal score, immune score, ESTIMATE score and Tumor Purity. (C) GO-GSEA analysis for six pathways potentially regulated by CuScore.

Fig. 4

3.5 Analysis of genomic mutation for CuScore

We displayed top 30 gene mutation for BC patients in the groups with higher CuScore and lower CuScore (Fig. 5A). For patients with higher CuScore, the top 30 mutation genes were listed as below: PIK3CA, TP53, TIN, GATA3, MUC16, KMT2C, МАР3K1, CDH1, SYNE1, FLG, NCORT, USH2A, PTEN, SPTA1, HMCNI, NEB, MAP2K4, ZEHX4, ARID1A, RYR2, DMD, CSMD1, CSMD3, HUWE1, MUC4, RYR3, SPEN, SYNE2, AKT1, and RUNX1 (Fig. 5A). For patients with lower CuScore, the top 30 mutation genes were listed as below: TP53, РІKЗСА, TTN, CDH1, MUC16, GATA3, RYR2, HMCN1, KMT2C, MUC4, USH2A, DMD, CCDC168, SYNE1, АРОВ, МАР3K1, MUC17, SPTA1, ADGRV1, NEB, FLG, FAT3, NF1, PKHD1L1, PTEN, ZFHX4, UTBN, CSMD3, MDN1, and PCLO (Fig. 5A). The high mutation rate genes of the two subgroups were similar to some extent. GISTIC showed the Genomic landscape of High-CuScore group or Low-CuScore -group (Fig. 5B).Fig. 5 Genomic Mutation Analysis for CuScore. (A) Waterfall plot showing the top 30 mutation genes in the group with high and low CuScore.(B) Genomic landscape of patients with high and low CuScore.

Fig. 5

3.6 Analysis of immunotherapy and chemotherapy for CuScore

We evaluated the critical role of immunomodulators in BC. Overall, some immunomodulators were mainly associated with lower CuScore, while other immunomodulators (such as MICA, CD276, TNFSF4, VEGFA, and EDNRB) were associated with higher CuScore (Fig. 6A). We observed that patients with higher CuScore had lower levels of TMB (Fig. 6B). As we could see from the violin plot, the IC50 levels of Sapitinib and LCL161 were lower in BC patients with higher CuScore than that of BC patients with lower CuScore, suggesting BC patients with higher CuScore may be more sensitive to Sapitinib and LCL161 (Fig. 6C). We also performed the Submap analysis of CuScore to explore immunotherapy differences between the two subgroups related to CuScore. We observed that BC patients with lower CuScore may respond to anti-PD1 therapy (Fig. 6D).Fig. 6 Analysis of Immunotherapy and Chemotherapy in Relation to CuScore. (A) The relationship of CuScore with the immunomodulators. (B) Violin plot showing the IC50 and TMB in high- and low-CuScore groups. (C) Violin plot showing the IC50 for Sapitinib and LCL161 in high- and low-CuScore groups. (D) SubMap analysis for CuScore.

Fig. 6

3.7 Preliminary inquiry of MMP13 in TCGA-BRCA

On the basis of Random Survival Forest and LASSO algorithms, three genes with more prognostic potential, including CEACAM5, MMP13, and CRISP3. MMP13 was of great concern to us, and we selected MMP13 for preliminary inquiry in TCGA-BRCA cohort. MMP13 expressed higher in BC than normal (Fig. 7A) and had a significant poorer clinical outcomes (Fig. 7B). To further investigate the biological function of MMP13, we implemented KEGG analysis to find that MMP13 may affect the development and occurrence of BC through the following pathways: CELL CYCLE and JAK STAT SIGNALING PATHWAY. et al. (Fig. 7C).Fig. 7 Preliminary inquiry and In vitro cell assay for MMP13. (A) The difference of MMP13 between normal and tumor in TCGA. (B) KM curve of MMP13 in TCGA. (C) KEGG analysis of cancer-related signaling pathways positively modulated by MMP13. (D) The mRNA expression of MMP13 was detected in BT-474, MCF-10A, and MCF7 cells. (E) Immunofluorescence staining was carried out on BT-474, MCF-10A, and MCF7 cells. (F–G) Expression levels of MMP13 protein in BT-474, MCF-10A, and MCF7 cells. (H–I) Western blot examined knockdown of MMP13 in MCF7. (J) Effects of knockdown of MMP13 expression on viability of MCF-sh1 and MCF7-sh2 cells were assessed by CCK8. (K) Effects of MMP13 low expression on migration of MCF-sh1 and MCF7-sh2 cells were examined by Boyden chamber assay. (L) Effects of MMP13 low expression on migration of MCF-sh1 and MCF7-sh2 cells were examined by Wound healing assay.

Fig. 7

3.8 In vitro cell assay for MMP13

Initially, to validate the function of MMP3, compared to normal cells, we measured that MMP13 was elevated in both mRNA and protein expression in BC cells. Notably, MMP13 was predominantly localized in the nucleus (Fig. 7D–G). Two MCF7 cell lines with silenced MMP13 expression were successfully established (Fig. 7H and I). The CCK8 results indicated that MMP13 silencing could significantly inhibited the MCF7 cells proliferation (Fig. 7J). Additionally, Boyden chamber and wound-healing assays indicated that reducing MMP13 expression may impair the migratory ability of MCF7 cells (Fig. 7K and L).

4 Discussion

Breast cancer, the most common malignancy among females, continues to pose a huge threat to female's health and survival [26,27]. Prognostic evaluation has been a critical focus in research, typically relying on the TNM staging system and molecular markers such as Her-2 [28], ER, and PR [29] expression. However, the high heterogeneity of BC means that some patients do not respond well to treatments based on current evaluation methods [30]. This underscores the need for new biomarkers to improve BC stratification protocols. Copper metabolism has been linked to tumorigenesis [31], with cancer cells displaying a greater demand for copper compared to normal cells [[32], [33], [34], [35]]. Various cancers, including breast [33] and cervical cancers [35], show increased copper levels in tumor tissues and/or altered systemic copper distribution. In present study, we firstly identified four copper metabolism genes (UBE2D2, SLC31A1, ATP7A, and MAPK1) with prognostic relevance for BC. ATP7A plays a key role in regulating Cu2+ transport, excretion, and ATP hydrolysis, facilitating copper efflux from the cytoplasm across the basolateral membrane in non-hepatic tissues [36,37]. Abnormal activation of the MAPK signaling pathway Linc-RoR promotes [31,38], which can lead to altered metabolic pathways and gene expression, has been implicated in cancers such as lung [39] and pancreatic cancer [40]. Copper, an essential micronutrient with unique redox properties, is crucial for cell growth, making copper homeostasis vital. Disruption in copper metabolism not only affects tumor cells but also interacts with the tumor microenvironment, promoting tumor progression through modulation and remodeling of this environment [41,42].

Recently, previous researches had been devoted to developing and validating biomolecular markers that can not only provide prognostic information, but more importantly, predict patient response to treatment. Currently, several biomolecular markers based on multiple genes have been developed to beyond traditional clinical factors. For example, OncotypeDX and MammaPrint were shown to predict adjuvant chemotherapy responses in the early BC patients with positive of HR+/HER2- [43]. A prognostic model containing three genes (CEACAM5, MMP13, and CRISP3) was constructed by the Random Survival Forest algorithm and Lasso regression and the risk scoring formula was obtained. Then the reliability of the prediction ability of this prognosis model was verified by the external datasets.

Immunology therapy has improved the prognosis for BC patients, yet many still progress to metastatic disease, often because current anti-cancer drugs primarily target cancer cells. The complex interactions between these TME cells and cancer cells occur through direct contact, extracellular matrix components, and soluble factors, creating a microenvironment that can either promote or inhibit cancer progression. Tumor-infiltrating immune cells, for example, can eliminate immunogenic tumor cells, thereby hindering tumor growth, but they can also contribute to treatment resistance by fostering tumor immunogenicity and selecting for drug-resistant clones that evade the immune response. In our study, we examined the infiltrating abundance of immune cell populations across different CuScore groups and found that immunomodulators were predominantly associated with lower CuScore. Notably, BC patients with lower CuScore may be more responsive to anti-PD1 therapy. Our findings offer insights into the immune microenvironment and potential immunotherapy strategies for BC.

Among these genes, the MMPs are the main enzymes that degrade extracellular matrix, which are closely related to tumor invasion, metastasis, and prognosis [44]. MMP13 can reduce all extracellular matrix and basement membrane components except collagen, and promote tumor invasion and metastasis [44]. Huang et al. found that the positive rate of MMP13 in colorectal cancer tissues was 52.5 % (42/80), and its positive rate was significantly related to clinical stage, size, lymph node status and postoperative recurrence [45]. Stokes et al. showed that MMP13 was highly expressed in head and neck tumor tissues, and the larger the tumor size, the higher the expression intensity [46]. Chang et al. found that MMP13 was the specificity gene in the diagnosis of BC [47]. The degradation of extracellular matrix by MMP13 is related to the metastasize ability of cancer cells [48]. In this study, we found that MMP13 expression was increased in BC cell lines and MMP13 was mainly distributed in the nucleus. MMP13 was proved to affect the migration and proliferation of BC cell lines.

In Conclusion, the copper metabolism-related signature identified could predict prognosis and guide treatment for BC patients. MMP13 may serve as a malignant factor in BC.

Data availability statement

The RNA profiles are able to be gained from the public cohorts, including The Cancer Genome Atlas (TCGA), Molecular Taxonomy of Breast Cancer International Consortium and Gene Expression Omnibus (GEO). Further inquiries can be directed to the corresponding author.

Funding

This study was supported by Key Project of 10.13039/501100012175 Zhejiang Provincial Administration of Traditional Chinese Medicine (ZYJ23JS03 ) and 10.13039/501100004731 Natural Science Foundation of Zhejiang Province (LGF22H080009 ).

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

CRediT authorship contribution statement

Chaojie Han: Writing – review & editing, Writing – original draft, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Zhangyang Feng: Writing – original draft, Validation, Methodology, Investigation, Formal analysis, Data curation. Yingjian Wang: Writing – original draft, Methodology, Investigation, Formal analysis, Data curation. Mengsi Hu: Methodology, Formal analysis. Shoufang Xu: Software, Resources, Methodology. Feiyu Jiang: Software, Resources, Formal analysis, Data curation. Yetao Han: Methodology, Investigation, Formal analysis, Data curation. Zhiwei Liu: Writing – review & editing, Investigation, Funding acquisition, Conceptualization. Yunsen Li: Writing – review & editing, Writing – original draft, Data curation, Conceptualization.

Declaration of competing interest

The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:Zhiwei Liu reports financial support was provided by 10.13039/501100004835 Zhejiang University 10.13039/100007269 School of Medicine 10.13039/100017238 Sir Run Run Shaw Hospital . Yingjian Wang reports financial support was provided by Institutes of Biology and Medical Sciences, 10.13039/501100007824 Soochow University . If there are other authors, they declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix A Supplementary data

The following is the supplementary data to this article:Multimedia component 1

Multimedia component 1

Appendix A Supplementary data to this article can be found online at https://doi.org/10.1016/j.heliyon.2024.e36445.
==== Refs
References

1 Sung H. Global cancer statistics 2020: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries CA A Cancer J. Clin. 71 2021 209 249 10.3322/caac.21660
2 Lin H.J. Liu Y. Lofland D. Lin J. Breast cancer tumor microenvironment and molecular aberrations hijack tumoricidal immunity Cancers 14 2022 10.3390/cancers14020285
3 Loibl S. Poortmans P. Morrow M. Denkert C. Curigliano G. Breast cancer Lancet 397 2021 1750 1769 10.1016/s0140-6736(20)32381-3 33812473
4 Harbeck N. Gnant M. Breast cancer Lancet 389 2017 1134 1150 10.1016/s0140-6736(16)31891-8 27865536
5 Bertucci F. Genomic characterization of metastatic breast cancers Nature 569 2019 560 564 10.1038/s41586-019-1056-z 31118521
6 Antoniou A.C. Foulkes W.D. Tischkowitz M. Breast-cancer risk in families with mutations in PALB2 N. Engl. J. Med. 371 2014 1651 1652 10.1056/NEJMc1410673 25337756
7 André F. Alpelisib for PIK3CA-mutated, hormone receptor-positive advanced breast cancer N. Engl. J. Med. 380 2019 1929 1940 10.1056/NEJMoa1813904 31091374
8 Segal C.V. Dowsett M. Estrogen receptor mutations in breast cancer--new focus on an old target Clin. Cancer Res. 20 2014 1724 1726 10.1158/1078-0432.Ccr-14-0067 24583794
9 Duffy M.J. O'Donovan N. McDermott E. Crown J. Validated biomarkers: the key to precision treatment in patients with breast cancer Breast 29 2016 192 201 10.1016/j.breast.2016.07.009 27521224
10 Zaki M.J. Karypis G. Yang J. Data mining in bioinformatics (BIOKDD) Algorithm Mol. Biol. 2 4 2007 10.1186/1748-7188-2-4
11 Holzinger A. Dehmer M. Jurisica I. Knowledge Discovery and interactive Data Mining in Bioinformatics--State-of-the-Art, future challenges and research directions BMC Bioinf. 15 Suppl 6 2014 I1 10.1186/1471-2105-15-s6-i1
12 Linder M.C. Hazegh-Azam M. Copper biochemistry and molecular biology Am. J. Clin. Nutr. 63 1996 797s 811s 10.1093/ajcn/63.5.797 8615367
13 Li Y. Trush M.A. DNA damage resulting from the oxidation of hydroquinone by copper: role for a Cu(II)/Cu(I) redox cycle and reactive oxygen generation Carcinogenesis 14 1993 1303 1311 10.1093/carcin/14.7.1303 8392444
14 Huang W. Functional molecule-mediated assembled copper nanozymes for diabetic wound healing J. Nanobiotechnol. 21 2023 294 10.1186/s12951-023-02048-1
15 Wang X. Conversion of senescent cartilage into a pro-chondrogenic microenvironment with antibody-functionalized copper sulfate nanoparticles for efficient osteoarthritis therapy J. Nanobiotechnol. 21 2023 258 10.1186/s12951-023-02036-5
16 Guo H. Copper induces spleen damage through modulation of oxidative stress, apoptosis, DNA damage, and inflammation Biol. Trace Elem. Res. 200 2022 669 677 10.1007/s12011-021-02672-8 33740180
17 Ge E.J. Connecting copper and cancer: from transition metal signalling to metalloplasia Nat. Rev. Cancer 22 2022 102 113 10.1038/s41568-021-00417-2 34764459
18 Tibshirani R. The lasso method for variable selection in the Cox model Stat. Med. 16 1997 385 395 10.1002/(sici)1097-0258(19970228)16:4<385::aid-sim380>3.0.co;2-3 9044528
19 Mermel C.H. GISTIC2.0 facilitates sensitive and confident localization of the targets of focal somatic copy-number alteration in human cancers Genome Biol. 12 2011 R41 10.1186/gb-2011-12-4-r41 21527027
20 Yoshihara K. Inferring tumour purity and stromal and immune cell admixture from expression data Nat. Commun. 4 2013 2612 10.1038/ncomms3612 24113773
21 Li T. TIMER: a web server for comprehensive analysis of tumor-infiltrating immune cells Cancer Res. 77 2017 e108 e110 10.1158/0008-5472.Can-17-0307 29092952
22 Hänzelmann S. Castelo R. Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data BMC Bioinf. 14 7 2013 10.1186/1471-2105-14-7
23 Liberzon A. The Molecular Signatures Database (MSigDB) hallmark gene set collection Cell Syst 1 2015 417 425 10.1016/j.cels.2015.12.004 26771021
24 Yang W. Genomics of Drug Sensitivity in Cancer (GDSC): a resource for therapeutic biomarker discovery in cancer cells Nucleic Acids Res. 41 2013 D955 D961 10.1093/nar/gks1111 23180760
25 Maeser D. Gruener R.F. Huang R.S. oncoPredict: an R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data Briefings Bioinf. 22 2021 10.1093/bib/bbab260
26 Singh D. Assaraf Y.G. Gacche R.N. Long non-coding RNA mediated drug resistance in breast cancer Drug Resist. Updates : reviews and commentaries in antimicrobial and anticancer chemotherapy 63 2022 100851 10.1016/j.drup.2022.100851
27 Singh J. Sah B. Shen Y. Liu L. Histone methyltransferase inhibitor UNC0642 promotes breast cancer cell death by upregulating TXNIP-dependent oxidative stress Chem. Biol. Interact. 385 2023 110720 10.1016/j.cbi.2023.110720
28 Marchiò C. Evolving concepts in HER2 evaluation in breast cancer: heterogeneity, HER2-low carcinomas and beyond Semin. Cancer Biol. 72 2021 123 135 10.1016/j.semcancer.2020.02.016 32112814
29 Yip C.H. Rhodes A. Estrogen and progesterone receptors in breast cancer Future Oncol. 10 2014 2293 2301 10.2217/fon.14.110 25471040
30 Liang Y. Zhang H. Song X. Yang Q. Metastatic heterogeneity of breast cancer: molecular mechanism and potential therapeutic targets Semin. Cancer Biol. 60 2020 14 27 10.1016/j.semcancer.2019.08.012 31421262
31 Geng R. Copper deprivation enhances the chemosensitivity of pancreatic cancer to rapamycin by mTORC1/2 inhibition Chem. Biol. Interact. 382 2023 110546 10.1016/j.cbi.2023.110546
32 Atakul T. Altinkaya S.O. Abas B.I. Yenisey C. Serum copper and zinc levels in patients with endometrial cancer Biol. Trace Elem. Res. 195 2020 46 54 10.1007/s12011-019-01844-x 31399869
33 Feng Y. Serum copper and zinc levels in breast cancer: a meta-analysis J. Trace Elem. Med. Biol. 62 2020 126629 10.1016/j.jtemb.2020.126629
34 Pavithra V. Serum levels of metal ions in female patients with breast cancer J. Clin. Diagn. Res. 9 2015 BC25 c27 10.7860/jcdr/2015/11627.5476 25737978
35 Zhang M. Shi M. Zhao Y. Association between serum copper levels and cervical cancer risk: a meta-analysis Biosci. Rep. 38 2018 10.1042/bsr20180161
36 Li Y.Q. Yin J.Y. Liu Z.Q. Li X.P. Copper efflux transporters ATP7A and ATP7B: novel biomarkers for platinum drug resistance and targets for therapy IUBMB Life 70 2018 183 191 10.1002/iub.1722 29394468
37 Öhrvik H. Aaseth J. Horn N. Orchestration of dynamic copper navigation - new and missing pieces Metallomics 9 2017 1204 1229 10.1039/c7mt00010c 28685789
38 Icard P. The strategic roles of four enzymes in the interconnection between metabolism and oncogene activation in non-small cell lung cancer: therapeutic implications Drug Resist. Updates : reviews and commentaries in antimicrobial and anticancer chemotherapy 63 2022 100852 10.1016/j.drup.2022.100852
39 Leonetti A. Resistance mechanisms to osimertinib in EGFR-mutated non-small cell lung cancer Br. J. Cancer 121 2019 725 737 10.1038/s41416-019-0573-8 31564718
40 Meng L.D. Linc01232 promotes the metastasis of pancreatic cancer by suppressing the ubiquitin-mediated degradation of HNRNPA2B1 and activating the A-Raf-induced MAPK/ERK signaling pathway Cancer Lett. 494 2020 107 120 10.1016/j.canlet.2020.08.001 32814086
41 Wang Y. Cuproptosis: a novel therapeutic target for overcoming cancer drug resistance Drug Resist. Updates : reviews and commentaries in antimicrobial and anticancer chemotherapy 72 2024 101018 10.1016/j.drup.2023.101018
42 Chen L. A global meta-analysis of heavy metal(loid)s pollution in soils near copper mines: evaluation of pollution level and probabilistic health risks Sci. Total Environ. 835 2022 155441 10.1016/j.scitotenv.2022.155441
43 Kwon M.J. Emerging immune gene signatures as prognostic or predictive biomarkers in breast cancer Arch Pharm. Res. (Seoul) 42 2019 947 961 10.1007/s12272-019-01189-y
44 Li S. Pritchard D.M. Yu L.G. Regulation and function of matrix metalloproteinase-13 in cancer progression and metastasis Cancers 14 2022 10.3390/cancers14133263
45 Huang M.Y. MMP13 is a potential prognostic marker for colorectal cancer Oncol. Rep. 24 2010 1241 1247 20878116
46 Stokes A. Expression profiles and clinical correlations of degradome components in the tumor microenvironment of head and neck squamous cell carcinoma Clin. Cancer Res. 16 2010 2022 2035 10.1158/1078-0432.Ccr-09-2525 20305301
47 Chang H.J. MMP13 is potentially a new tumor marker for breast cancer diagnosis Oncol. Rep. 22 2009 1119 1127 10.3892/or_00000544 19787229
48 Ellsworth R.E. A gene expression signature that defines breast cancer metastases Clin. Exp. Metastasis 26 2009 205 213 10.1007/s10585-008-9232-9 19112599
