Identification of novel neutrophil-extracellular-traps-related genes as biomarkers for breast cancer prognosis and immunotherapy

Introduction

Background

Breast cancer (BC) is a prevalent malignancy affecting women worldwide (1). While early-stage BC can often be cured with surgery and adjuvant therapy, the prognosis for patients with metastatic BC is extremely poor (2). Current prognostic indicators, such as the Breast Cancer Index (BCI) in correlation with traditional clinicopathological factors [e.g., tumor size, lymph node status, tumor grade, estrogen receptor expression, and human epidermal growth factor receptor 2 (HER2) expression], are useful for identifying patients at risk of relapse and predicting the efficacy of adjuvant chemotherapy and immunotherapy. However, these methods may be inaccurate in predicting early recurrence and metastasis. Once metastasis occurs, particularly in multiple organs, treatment becomes palliative rather than curative (3). Consequently, understanding the mechanisms underlying BC progression is crucial for reducing metastasis and improving patient survival.

Rationale and knowledge gap

Recent studies have highlighted the strong correlation between tumor onset and progression with chronic inflammation (4,5). Inflammation promotes cancer cell proliferation, adhesion, and migration by producing various growth factors and chemoattractants (6). Neutrophils, key players in inflammatory diseases, contribute to the pathological process of both acute and chronic inflammation through abnormal increases or decreases in their numbers. Initially considered as sensitive indicators of infection or early tissue damage in the peripheral circulation, neutrophils also play complex roles in the tumor microenvironment (TME). They can inhibit tumor progression by releasing antitumor factors (7). Furthermore, neutrophils release reticulated fiber structures known as neutrophil extracellular traps (NETs) (8). NETs are implicated in cancer progression, immune escape, metastasis, and cancer-associated thrombosis (9-11). In BC, recent research has unveiled the role of NETs in metastasis, with tumor-secreted proteinase C promoting breast-to-lung metastasis by regulating neutrophil recruitment and NET formation (12). Yang et al. demonstrated that NAMPT, secreted by tumors, induces SIRT1 activation in senescent neutrophils, facilitating lung metastasis in BC (13). Despite the emerging evidence of NETs’ significant role in BC metastasis, NET-related genes (NETRGs) have been infrequently reported as biomarkers for tumorigenesis and progression in BC.

Objective

This study aims to address this gap by downloading BC transcriptomes, clinical data, and NETRGs from public databases. Bioinformatics methods were used to identify prognostic genes associated with BC and NETs, and their prognostic value was evaluated. Additionally, the relationship between prognostic models, immune scores, and drug sensitivity was explored, offering insights relevant to BC prognosis and treatment. This manuscript adheres to the TRIPOD reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-1826/rc).

Methods

Datasets

The Cancer Genome Atlas-Breast cancer (TCGA-BRCA) dataset was retrieved from the TCGA repository (https://www.cancer.gov/ccg/research/genome-sequencing/tcga), comprising transcriptome data from 1,097 BC samples, with overall survival (OS) and clinical information available for 1,050 samples, along with 114 adjacent normal samples. The GSE42568 dataset, obtained from the GEO repository (https://www.ncbi.nlm.nih.gov/geo/), included transcriptome data, clinical, and survival data from 104 BC biopsy samples (GPL570), serving as the external validation set. Additionally, 137 NETRGs were sourced from a prior study (14). The flowchart of bioinformatics analysis is shown in Figure 1.

Figure 1 Flowchart of bioinformatics analysis. BRCA, breast cancer; DCA, decision curve analysis; DEGs, differential expressed genes; DE-NETRGs, differentially expressed NETRGs; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; LASSO, least absolute shrinkage and selection operator; NETRGs, neutrophil extracellular trap-related genes; qRT-PCR, quantitative reverse transcription polymerase chain reaction; ROC, receiver operating characteristic; TCGA, the cancer genome atlas; TIDE, tumor immune dysfunction and exclusion.

Identification and functional enrichment of differentially expressed NETRGs (DE-NETRGs)

Differentially expressed genes (DEGs) between BC and control groups were identified in the TCGA-BRCA dataset with a false discovery rate (FDR) <0.05 and |log₂fold change (FC)| ≥1. Differential analysis was conducted using the edgeR package (v 3.36.0) (15). Volcano plots and heatmaps were generated to visualize the results using the ggplot2 package (v 3.3.5) and the pheatmap package (v 1.0.12). DE-NETRGs were selected by overlapping DEGs and NETRGs. Enrichment analysis of DE-NETRGs was performed using the clusterProfiler package, with a focus on Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases (v 4.0.2) (16), applying a count threshold of ≥1 and an adjusted P<0.05.

Screening for biomarkers and construction of risk model

The 1,050 BC samples were split into training and internal validation sets using a 7:3 ratio. In order to evaluate whether the candidate genes were related to the patients’ survival, the expression levels of the candidate genes were taken as continuous variables and compared with the survival time. Univariate Cox analysis was performed using the survival package (v 3.2-13) to screen for genes related to survival (P<0.1) (17), Subsequently, the glmnet package (v 4.1-3) was used to conduct the least absolute shrinkage and selection operator (LASSO) regression analysis on the expression levels of the survival-related genes according to the patients’ survival information, and biomarkers were obtained from DE-NETRG. Based on risk scores calculated from the biomarkers, patients in the training set were divided into high- and low-risk groups according to median values.

Riskscore=∑1ncoef[genei×expression(genei)][1]

Kaplan-Meier (K-M) survival curves were generated using the survmine package (v 0.4.9), and the survivalROC package (v 1.0.3) (18) was used to calculate the area under the curve (AUC) values for the receiver operating characteristic (ROC) curves, evaluating the model’s predictive accuracy. The validity of the risk model was further confirmed through both internal and external validation using the GSE42568 dataset.

The relationship between biomarker expression and prognosis in BC

Biomarker expression levels were compared between BC and control samples in the TCGA-BRCA dataset. Based on the optimal expression threshold, BC samples were stratified into high- and low-expression groups, and K-M curves were generated to assess the relationship between biomarker expression and survival outcomes.

Clinical features and survival analysis

The differences in risk scores across various clinical subgroups were first analyzed. Survival analysis, stratified by clinical characteristics such as age, stage, race, and pathologic T, M, and N stages, was conducted using the survminer package (v 0.4.9) (https://cran.r-project.org/package=survminer), with K-M curves illustrating the results.

Independent prognostic analysis

Independent prognostic analyses based on clinical characteristics (e.g., age, gender, pathologic M) and risk scores were performed using the TCGA-BRCA dataset. Variables with statistical significance (P<0.05) were identified through univariate Cox analysis, followed by multivariate Cox analysis to determine independent prognostic factors. A nomogram was constructed and visualized using the rms package (v 6.1-0) (19). The model’s performance was evaluated through ROC curves, calibration curves, and decision curve analysis (DCA) curves.

Immune score, checkpoint, and immunotherapy sensitivity analyses

Using the estimate package (v 1.0.13) (https://R-Forge.R-project.org/projects/estimate/), the ESTIMATE algorithm was employed to calculate stromal, immune, and ESTIMATE scores as well as tumor purity in the two risk subgroups (P<0.05), with violin plots visualizing the results. Correlations between risk scores and the four scores were analyzed using the Spearman method, and scatter plots were used for visualization. Differences in the expression of 48 immune checkpoints between the two risk subgroups were assessed using the Wilcoxon test (P<0.05), with box plots presenting the findings using RColorBrewer package (v 1.1.3) (https://CRAN.R-project.org/package=RColorBrewer). Additionally, differences in immunotherapy responses were explored by comparing immunophenoscores (IPS) and immune checkpoint [immunomodulator (CP)] scores between the subgroups. Responder prediction, TIDE scores, dysfunction scores, and exclusion scores were also compared using the TIDE database using the ggplot2 package (version 3.3.5).

Immunoinfiltration analysis

The abundance of immune-infiltrating cells across all TCGA-BRCA samples was evaluated with the CIBERSORT algorithm (20). Variations in immune cell abundance between the risk subgroups were assessed via the Wilcoxon rank sum test (P<0.05), and the relationships between risk scores and differential immune cell populations were depicted using scatter plots. R>0.7 indicated strong correlation, 0.3<R<0.7 indicated moderate correlation, and R<0.3 indicated weak correlation.

Drug sensitivity analysis

The chemotherapy response determined by the half-maximal inhibitory concentration (IC50) for each patient sample was analyzed using the Genomics of Drug Sensitivity in Cancer (GDSC) database (https://www.cancerrxgene.org/) and the oncoPredict package (v 0.2), and the responses of the high- and low-risk groups to anticancer drugs were evaluated (21). Variations in anticancer drug sensitivity between the two risk subgroups were further evaluated to identify potential therapeutic implications.

Quantitative reverse transcription polymerase chain reaction (qRT-PCR)

A total of 10 BC and 10 control samples were obtained from the Ethics Committee of Hospital of Guangxi Zhuang Autonomous Region. Total RNA was isolated from the samples using TRIzol reagent (Ambion, Austin, TX, USA) according to the manufacturer’s protocol. RNA concentration was measured with the NanoPhotometer N50. cDNA was synthesized using the SureScript First Strand cDNA Synthesis Kit (Servicebio, Wuhan, China). The qRT-PCR reaction mix consisted of 3 µL of cDNA, 5 µL of 2× Universal Blue SYBR Green qPCR Master Mix, and 1 µL each of forward and reverse primers. Primer sequences are provided in Table 1. Gene expression levels were analyzed using the 2^-ΔΔCT method, with GAPDH as an internal control (22). Graphs were generated, and P values were calculated using GraphPad Prism 5.

Ethical statement

The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013). The study was approved by People’s Hospital of Guangxi Zhuang Autonomous Region (July 22, 2023/KY-KJT-2023-160) and informed consent was obtained from all individual participants.

Statistical analysis

Bioinformatics analyses were performed in R, and variance analysis was conducted using the Wilcoxon test.

Results

Identifying DE-NETRGs and exploring their enrichment functions and pathways

A total of 5,590 DEGs were identified between the BC and control groups, including 3,478 upregulated and 2,112 downregulated genes (Figure S1A,S1B; table available at https://cdn.amegroups.cn/static/public/tcr-24-1826-1.xlsx). By intersecting DEGs with NETRGs, 37 DE-NETRGs were identified (Figure S1C; Table S1). Functional enrichment analysis revealed that the DE-NETRGs were involved in 393 GO biological process (BP) terms, 18 GO cellular component (CC) terms, 39 GO molecular function (MF) terms, and 28 KEGG pathways. The GO enrichment analysis indicated that DE-NETRGs were mainly associated with leukocyte migration and the humoral immune response, etc. (Figure S1D; table available at https://cdn.amegroups.cn/static/public/tcr-24-1826-2.xlsx). KEGG pathway analysis revealed involvement in the NET formation and interleukin (IL)-17 signaling pathway, etc. (Figure S1E; Table S2). The enrichment analysis results suggest that these genes are significantly correlated with functions such as NETs, IL-17 signaling pathway, TNF signaling pathway, etc.

F2RL2, AZU1, IL33, and ELANE were biomarkers that might influence prognosis in patients with BC

Survival-associated genes were identified through univariate Cox analysis, yielding 4 significant genes (Figure 2A). Further, LASSO regression analysis identified F2RL2, AZU1, IL33, and ELANE as biomarkers (Figure 2B,2C). Based on the median risk score, patients were divided into high- (n=365) and low-risk (n=370) groups (Figure 2D). Survival analysis in the training set showed that the high-risk group exhibited lower survival rates (P<0.001), as indicated by the K-M curves (Figure 2E). The ROC curves demonstrated strong predictive performance, with AUC values greater than 0.60 for 1-, 2-, 3-, 4-, and 5-year survival (Figure 2F). Both internal and external validation sets were employed to assess the model’s predictive accuracy. The risk profile plots and survival curves from the validation sets aligned with those of the training set (Figure S2A-S2D). AUC values exceeded 0.60 for both validation sets (Figure S2E,S2F).

Figure 2 Development of a predictive model using DE-NETRGs in the TCGA-BRCA dataset. (A) Forest plot depicting the HRs for DE-NETRGs in BC. (B) Plot of the relationship between partial likelihood deviance and the natural logarithm of the tuning parameter (λ) derived from the LASSO Cox regression model. (C) Coefficients of selected features in the model shown in relation to the lambda parameter. (D) Distribution of the risk score curve, survival states, and gene expression in BC cohorts from the TCGA-BRCA dataset. (E) Kaplan-Meier survival curves comparing overall survival between patients in the high-risk and low-risk groups. (F) Time-dependent ROC curves assessing the predictive performance of the risk score. AUC, area under the curve; BC, breast cancer; CI, confidence interval; DE-NETRGs, differentially expressed NETRGs; LASSO, least absolute shrinkage and selection operator; NETRGs, neutrophil extracellular trap-related genes; OS, overall survival; ROC, receiver operating characteristic.

Prognostic analysis of biomarkers and survival analysis of risk score

Expression analysis revealed that ELANE and IL33 were significantly downregulated, while F2RL2 and AZU1 were upregulated in BC samples (Figure 3A). Survival analysis showed that higher expression levels of these biomarkers were associated with improved survival outcomes (Figure 3B-3E).

Figure 3 Prognostic analysis of biomarkers. (A) Expression analysis of model genes based on the TCGA-BRCA dataset, with yellow-brown and blue colors representing BC and control samples, respectively. **, P<0.01; ****, P<0.0001. (B-E) Kaplan-Meier survival curves demonstrating overall survival differences between high- and low-expression groups of the biomarkers. BC, breast cancer; BRCA, breast cancer; TCGA, the cancer genome atlas.

Clinical characteristic analysis revealed significant differences in risk scores across subgroups, including stage, race, radiation therapy, pathologic T stage, and surgical procedure (Figure S3A). Survival analysis further revealed that the pathologic N and stage III + IV subgroups exhibited notable differences in survival rates between the risk groups (Figure S3B).

The prognostic nomogram exhibited superior predictive accuracy

An independent prognostic analysis identified five key factors influencing prognosis: radiation therapy, age, tumor stage, pathologic N, and risk score (Figure 4A,4B). Using these factors, a nomogram was constructed to predict OS (Figure 4C). The calibration curve demonstrated good accuracy for the nomogram (Figure 4D), while DCA confirmed its superior clinical utility (Figure 4E). ROC curve analysis further validated the nomogram, with all AUC values exceeding 0.80 (Figure 4F-4H).

Figure 4 Creation and assessment of a nomogram model. (A,B) Univariate (A) and multivariate (B) Cox regression analyses of the TCGA-BRCA dataset, identifying significant factors (radiation therapy, age, tumor stage, pathologic N, and risk score) related to overall survival (P<0.05). (C) Prognostic nomogram incorporating the risk score and clinical factors to predict 1-, 3-, and 5-year survival outcomes. **, P<0.01. (D,E) Calibration (D) and DCA curves (E) for the prediction of 1-, 3-, and 5-year overall survival in patients with BC using the nomogram. (F-H) ROC curves for the prediction of 1-, 3-, and 5-year overall survival in patients with BCs using the nomogram. AUC, area under the curve; BC, breast cancer; BRCA, breast cancer; CI, confidence interval; DCA, decision curve analysis; HR, hazard ratio; OS, overall survival; Pr, probability; ROC, receiver operating characteristic; TCGA, The Cancer Genome Atlas.

Immune-related analysis between two risk subgroups

Following the identification of independent prognostic factors, immune-related analyses were conducted to assess the immune microenvironment between the two risk subgroups. The analysis revealed that the low-risk group exhibited elevated stromal and ESTIMATE scores, along with lower tumor purity (Figure S4A). Stromal cells play a more active role in nutrient supply, signal transduction, structural support, and preventing tumor cells from excessive and disorderly growth. A higher ESTIMATE score likely represented more immune cell infiltration or the presence of immune-related factors. The immune components inhibited tumor cell proliferation, migration, and invasion, which contributed to a lower tumor risk. In contrast, the high-risk group had higher tumor purity, indicating stronger tumor cell proliferation ability with fewer limitations from other components. Correlation analysis demonstrated a strong negative relationship between the risk score and stromal score (R=−0.31, P<2.2e−16). In the high-risk group, stromal components of the TME were reduced (Figure S4B). Differential analysis identified 30 immune checkpoints (such as TLA, IDO1, LAG3, etc.) with significant differences between the risk groups. Notably, the expression levels of IDO1, LAG3, TNFRSF9, CD80, and TNFSF15 were significantly higher in the high-risk group (Figure S4C). The impact of immunotherapy was also evaluated, revealing that the low-risk group exhibited higher IPS and CP scores (Figure S4D), stronger immune cell-mediated tumor cell killing, and higher immunogenicity. Moreover, using the TIDE algorithm, significant differences were found between the immune response (Respondent) and non-immune response (false), and the low-risk group was found to have higher TIDE and dysfunction scores (Figure S4E). A high TIDE score indicated an advantage in resisting immune evasion, while a high dysfunction score suggested that although immune cells had some dysfunction, tumor cells were still under control by the immune system and had not fully escaped.

F2RL2, AZU1, IL33, and ELANE were associated with multiple immune cells

Immune cell distribution across the two risk subgroups was analyzed, and the 22 immune cell types in each sample are displayed in Figure 5A. A total of 14 differential immune cell populations were identified between the risk subgroups (Figure 5B). Correlation analysis between the risk score and immune cells indicated that all risk scores were weakly correlated with differential immune cells (R<0.3), with resting mast cells showed the strongest negative correlation with the risk score (R=−0.27, P=1.3e−10), while M0 macrophages demonstrated the highest positive correlation (R=0.26, P=4.7e−10) (Figure 5C).

Figure 5 Immune infiltration analysis. (A) Distribution of 22 immune cell types in individual BC samples. (B) Boxplots showing the infiltration levels of 22 immune cell types between the high- and low-risk subgroups. ns, no significance; *, P<0.05; **, P<0.01; ***, P<0.001. (C) Correlation heatmap of the relationship between differentially expressed immune cells and risk score. BC, breast cancer.

Drug sensitivity of the risk groups

Drug sensitivity analysis indicated significant differences (P<0.05) in the sensitivities to 124 anticancer drugs between the two risk subgroups. The smaller the IC50 value of a drug was, the stronger its ability to inhibit cell growth, and the better its effect in treating cancer, with the high-risk group being more sensitive to 94 drugs and the low-risk group being more sensitive to 30 drugs. Drugs such as BMS-754807, sepantronium bromide, and UMI-77 were among those with notable differential responses in the high-risk group (Table S3).

Validation of the expression of four biomarkers

To further validate the findings, qRT-PCR was conducted on 10 BC and 10 control tissue samples to evaluate the expression of the biomarkers (Figure 6). The results indicated that, compared to the BC samples, ELANE expression was higher in the control samples (P=0.22), and the expression of IL33 (P=0.02) and F2RL2 (P=0.01) was significantly higher in the control group than in the BC group. AZU1 expression exhibited a significant difference between BC and control samples (P=0.04), but the trend was opposite to that observed in the public database, potentially due to tissue heterogeneity. In conclusion, the expression patterns of ELANE and IL33 were consistent with the bioinformatics analysis results.

Figure 6 Real-time quantitative polymerase chain reaction (qRT-PCR) analysis of biomarker expression. (A) AZU1 expression in BC and adjacent non-cancerous tissues. (B) IL33 expression in BC and adjacent non-cancerous tissues. (C) ELANE expression in BC and adjacent non-cancerous tissues. (D) F2RL2 expression in BC and adjacent non-cancerous tissues. ns, no significance; *, P<0.05. BC, breast cancer.

Discussion

The continued increase in the incidence and mortality rates of BC worldwide reflects the strongly aggressive and metastatic nature of BC. A recent study has reported that the patients with tumors may exhibit different prognoses and responses to treatment because of the molecular characteristics, despite having similar clinicopathological risk variables. Hence, it is necessary to explore novel molecular prognostic markers to supplement clinicopathological risk factors. The TME is significantly influenced by NETs, which promote the development, immune escape, and metastasis of tumor cells in different ways. In the TME, NETs can inhibit tumor cells from interactions with nearby antitumor immune cells, thereby altering the immune profile of the tumor and influencing its responsiveness to immunotherapy. However, NETRGs have been less frequently reported as biomarkers for predicting the occurrence and prognosis of BC. Therefore, it is essential to explore the relationship between prognostic genes based on NETs and the TME of patients with BC.

In this study, prognostic biomarkers for BC were explored using various analyses, and the results are of great significance. First, 137 NET biomarkers were initially screened and systematically applied using the TCGA-BRCA dataset to obtain 37 intersecting differential NET genes. To strengthen the clinical utility of our study, a risk model was constructed using four biomarkers (F2RL2, AZU1, IL33, and ELANE) that were found to be most highly associated with patient survival via univariate Cox regression, LASSO, and multivariate Cox regression analyses. The model was validated using GSE42568 data, which indicated that the risk model was able to predict the prognosis of BC. It was observed that individuals in the low-risk category consistently experienced a better prognosis compared to those in the high-risk category. Additionally, ROC curves proved the reliability of biomarkers, and the nomogram exhibited the superiority of the risk model in prognosis prediction. It was useful to assess the prognostic value of the risk model by means of ROC curves, with AUC values exceeding 0.7 at 1, 3, and 5 years, suggesting nomogram had a better ability to predict the prognosis of BC patients (11). These results suggested that the risk model had general applicability to BC and that NETs-associated biomarkers might be associated with BC development.

DE-NETRGs show great promise as potential predictive biomarkers for human cancer. Among these, F2RL2, which encodes protease-activated receptor-3 (PAR3), plays a pivotal role as a G-protein-coupled receptor and tumor suppressor. PAR3 is involved in homeostatic redox control and limits invasiveness in various cancers, including BC, glioblastoma, and prostate cancer (23,24). IL33 presents a dual role in tumour immunity (25-27). Numerous studies have shown that it is associated with tumourigenesis, proliferation and metastasis on the one hand, but on the other hand, it is able to activate the body’s immune effector mechanisms to inhibit tumour growth. In melanoma and lung cancer models, both overexpression and exogenous addition of IL33 significantly inhibited tumour lung metastasis (28). In BC, tumour cell-derived IL-3 induced interferon-γ production by natural killer (NK) and CD8⁺ T cells, which in turn inhibited the growth and metastasis of mouse BC. This phenomenon fits with our findings that we found IL33 to have a protective role in BC, as evidenced by the fact that higher protein levels of IL33 are associated with a better prognosis. Moreover, IL33 drives the anti-tumour immune response through activation of immune cells such as NK cells, T cells and macrophages, again in keeping with our findings. Together, these findings suggest that the role of IL33 in tumour immunity is complex and varied and that its manifestation in different tumour types has a certain regularity, which warrants further in-depth studies to explore its potential value in tumour therapy. The roles of ELANE in tumours are complex and diverse, ranging from the direct killing of cancer cells to the modulation of the tumour microenvironment and immune response. These findings provide new perspectives and potential targets for the development of new anti-cancer therapies (29). A recent study has indicated that in vitro exposure to azurocidin 1 (AZU1) promotes cellular sepsis by modulating the NF-κB/NLRP1/caspase-3/GSDMD axis in TNBC (30). Furthermore, AZU1 has been reported to inhibit aberrant proliferation in TNBC by regulating cellular pyroptosis (30).

These findings support the reliability of the risk models and offer new insights into the early diagnosis and prevention of BC. With the increasing use of cancer immunotherapies, such as immune checkpoint inhibitors, in BC treatment, the immune profile of the TME (TIME) has emerged as a critical factor influencing therapeutic response and clinical outcomes (31). A prior study has reported that TIME characteristics correlate with prognosis in BC (32). The present study demonstrated the significant roles of the four biomarkers (F2RL2, AZU1, IL33, and ELANE) in immunotherapy. Substantial differences in immune cell infiltration, immune functions, and immune checkpoint expression were identified between high- and low-risk BC patient groups. Specifically, the low-risk group exhibited increased levels of dendritic cells, activated NK cells, and CD8 T-cells (33-35), all of which are known to inhibit tumor growth and enhance immune responses and immunotherapy efficacy. Conversely, T follicular helper cells and regulatory T cells were more abundant in the high-risk group, promoting immune evasion and cancer metastasis (35). Similar to our expected results, patients in the low-risk group had higher IPS and CP scores, suggesting that patients in the low-risk group may respond better to immunotherapy. Surprisingly, patients in the high-risk group had lower TIDE scores, This result suggests that although the lower immune prognostic score in the high-risk group may imply a relatively weaker immune system, the low value of the TIDE score suggests that there may be a stronger mechanism of immune escape from the tumour in these patients. The low TIDE score may reflect the fact that the tumour escapes from immune surveillance through suppression of immune responses, which in turn leads to the diminished efficacy of immunotherapy. This phenomenon suggests that high-risk patients with poor immune function, despite showing poorer immune function in traditional immune prognostic scores, may actually be more resilient in terms of immune escape. The results indicate that biomarkers associated with TIME could guide targeted immunotherapy strategies, suggesting that even high-risk patients can significantly benefit from chemotherapy, immunotherapy, and targeted therapies (36).

NETs are a complex physiological and pathological process involving multiple molecular and cellular interactions, and have been shown to be closely associated with the development of a variety of diseases. Recent studies have shown that PAR2, a G protein-coupled receptor, plays an important role in the gene signature of neutrophils/NETs, and its gene expression is significantly and positively correlated with the formation of NETs. In particular, in BC cells, the TF/PAR2 signalling axis was found to promote the formation of NETs and thus tumour growth (37). Furthermore, in hepatocellular carcinoma (HCC), lenvatinib activated the Akt/mTOR signalling pathway by increasing the expression of NDUFA4L2, which in turn promoted the production of IL33, ultimately contributing to the formation of NETs in neutrophils (38). IL33 levels are elevated in rheumatoid arthritis (RA) patients and positively correlate with NETs formation.IL33 may further contribute to the release of NETs by activating signalling pathways such as MAPK and PI3K/Akt). In sepsis, proteomic analysis of the circulating pentraxin 3 (PTX3) complex revealed that its interaction with AZU1 and other components of NETs may play an important role in the formation of NETs (39). Notably, NETs consist of DNA backbone as well as particulate proteins such as ELANE, which not only play a role in the immune response, but also selectively kill a wide range of cancer cells and are less toxic to non-cancer cells (40,41). Therefore, the study of biomarkers such as PAR2, IL33, PTX3 and ELANE will not only help to gain insights into the mechanism of action of NETs in different diseases, but may also provide new targets for the early diagnosis and treatment of diseases such as BC and improve the prognosis of patients.

Although the four biomarkers (F2RL2, AZU1, IL33 and ELANE) in this study had the ability to predict prognosis, only ELANE and IL33 were expressed as expected in the control group in the TCGA database. Although validated by multiple databases and qRT-PCR, the link between prognostic-related biomarkers and their biological functions has not been deeply explored. This discrepancy is due to a number of factors, on the one hand, there are biases between databases and clinical samples, such as differences in collection methods and characteristics of the sample groups; on the other hand, the regulation of gene expression is complex and affected by unknown genes or signalling pathways, and the protective mechanisms of the four genes in BC may be more complex, and the study may not cover all the influencing factors. Therefore, in the future, larger datasets and larger cohort studies are needed to validate the accuracy of the NET risk model, collect more clinical samples, and conduct more ex vivo and in vivo experiments to verify the feasibility of the study.

Conclusions

In this study, four BC-prognosis-related biomarkers—F2RL2, AZU1, IL33, and ELANE—were identified, and a novel prognostic risk model was developed based on their data. This risk model demonstrated potential predictive utility for patients with BC, particularly in the context of immunotherapy. The model’s effectiveness was further supported by analyses of immune function, immune checkpoint genes, tumor-infiltrating immune cells, and drug sensitivity in patients with BC. Overall, this study identified new biomarkers associated with BC prognosis and provided potential therapeutic targets, contributing to improved treatment strategies and more accurate predictions of patient outcomes.

Acknowledgments

We extend our thanks to all colleagues who have provided assistance and support throughout the research process. Their collaboration and encouragement have been invaluable to our work.

Footnote

Reporting Checklist: The authors have completed the TRIPOD reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-1826/rc

Data Sharing Statement: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-1826/dss

Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-1826/prf

Funding: This work was supported by Guangxi Natural Science Foundation Project (No. 2024GXNSFBA010024 to Z.H.); Guangxi Zhuang Autonomous Region Health Commission self-funded research project (No. Z-A20230008 to Z.H.) and the Major Project of Science and Technology of Guangxi Zhuang Autonomous Region (No. Guike-AA22096018 to J.Y.).

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-24-1826/coif). Z.H. reports the funding from Guangxi Natural Science Foundation Project (No. 2024GXNSFBA010024) and Guangxi Zhuang Autonomous Region Health Commission (No. Z-A20230008). J.Y. reports the funding from the Major Project of Science and Technology of Guangxi Zhuang Autonomous Region (No. Guike-AA22096018). The other authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. The study was conducted in accordance with the Declaration of Helsinki (as revised in 2013). The study was approved by People’s Hospital of Guangxi Zhuang Autonomous Region (July 22, 2023/KY-KJT-2023-160) and informed consent was obtained from all individual participants.

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Wu B, Li Y, Shi B, et al. Temporal trends of breast cancer burden in the Western Pacific Region from 1990 to 2044: Implications from the Global Burden of Disease Study 2019. J Adv Res 2024;59:189-99. [Crossref] [PubMed]
2. Heer E, Harper A, Escandor N, et al. Global burden and trends in premenopausal and postmenopausal breast cancer: a population-based study. Lancet Glob Health 2020;8:e1027-37. [Crossref] [PubMed]
3. Wang T, McCullough LE, White AJ, et al. Prediagnosis aspirin use, DNA methylation, and mortality after breast cancer: A population-based study. Cancer 2019;125:3836-44. [Crossref] [PubMed]
4. Grivennikov SI, Greten FR, Karin M. Immunity, inflammation, and cancer. Cell 2010;140:883-99. [Crossref] [PubMed]
5. Metzemaekers M, Gouwy M, Proost P. Neutrophil chemoattractant receptors in health and disease: double-edged swords. Cell Mol Immunol 2020;17:433-50. [Crossref] [PubMed]
6. Kuang DM, Zhao Q, Wu Y, et al. Peritumoral neutrophils link inflammatory response to disease progression by fostering angiogenesis in hepatocellular carcinoma. J Hepatol 2011;54:948-55. [Crossref] [PubMed]
7. Hedrick CC, Malanchi I. Neutrophils in cancer: heterogeneous and multifaceted. Nat Rev Immunol 2022;22:173-87. [Crossref] [PubMed]
8. Brinkmann V, Reichard U, Goosmann C, et al. Neutrophil extracellular traps kill bacteria. Science 2004;303:1532-5. [Crossref] [PubMed]
9. Demers M, Wagner DD. NETosis: a new factor in tumor progression and cancer-associated thrombosis. Semin Thromb Hemost 2014;40:277-83. [Crossref] [PubMed]
10. Jung HS, Gu J, Kim JE, et al. Cancer cell-induced neutrophil extracellular traps promote both hypercoagulability and cancer progression. PLoS One 2019;14:e0216055. [Crossref] [PubMed]
11. Masucci MT, Minopoli M, Del Vecchio S, et al. The Emerging Role of Neutrophil Extracellular Traps (NETs) in Tumor Progression and Metastasis. Front Immunol 2020;11:1749. [Crossref] [PubMed]
12. Yu W, Wang Z, Dai P, et al. The activation of SIRT1 by resveratrol reduces breast cancer metastasis to lung through inhibiting neutrophil extracellular traps. J Drug Target 2023;31:962-75. [Crossref] [PubMed]
13. Yang C, Wang Z, Li L, et al. Aged neutrophils form mitochondria-dependent vital NETs to promote breast cancer lung metastasis. J Immunother Cancer 2021;9:e002875. [Crossref] [PubMed]
14. Wu J, Zhang F, Zheng X, et al. Identification of renal ischemia reperfusion injury subtypes and predictive strategies for delayed graft function and graft survival based on neutrophil extracellular trap-related genes. Front Immunol 2022;13:1047367. [Crossref] [PubMed]
15. Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 2010;26:139-40. [Crossref] [PubMed]
16. Wu T, Hu E, Xu S, et al. clusterProfiler 4.0: A universal enrichment tool for interpreting omics data. Innovation (Camb) 2021;2:100141. [Crossref] [PubMed]
17. Ye Y, Dai Q, Qi H. A novel defined pyroptosis-related gene signature for predicting the prognosis of ovarian cancer. Cell Death Discov 2021;7:71. [Crossref] [PubMed]
18. In J, Lee DK. Survival analysis: part II - applied clinical data analysis. Korean J Anesthesiol 2019;72:441-57. [Crossref] [PubMed]
19. Yang CC, Tsai CG, Yeh JA. Miniaturization of dielectric liquid microlens in package. Biomicrofluidics 2010;4:43006. [Crossref] [PubMed]
20. Newman AM, Liu CL, Green MR, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods 2015;12:453-7. [Crossref] [PubMed]
21. Maeser D, Gruener RF, Huang RS. oncoPredict: an R package for predicting in vivo or cancer patient drug response and biomarkers from cell line screening data. Brief Bioinform 2021;22:bbab260. [Crossref] [PubMed]
22. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. Methods 2001;25:402-8. [Crossref] [PubMed]
23. Liu S, Song A, Wu Y, et al. Analysis of genomics and immune infiltration patterns of epithelial-mesenchymal transition related to metastatic breast cancer to bone. Transl Oncol 2021;14:100993. [Crossref] [PubMed]
24. Lvu W, Fei X, Chen C, et al. In silico identification of the prognostic biomarkers and therapeutic targets associated with cancer stem cell characteristics of glioma. Biosci Rep 2020;40:BSR20201037. [Crossref] [PubMed]
25. Liu Q, Lu D, Hu J, et al. Dual role of interleukin-33 in tumors. Zhong Nan Da Xue Xue Bao Yi Xue Ban 2021;46:169-75. [Crossref] [PubMed]
26. Lu B, Yang M, Wang Q. Interleukin-33 in tumorigenesis, tumor immune evasion, and cancer immunotherapy. J Mol Med (Berl) 2016;94:535-43. [Crossref] [PubMed]
27. Akimoto M, Takenaga K. Role of the IL-33/ST2L axis in colorectal cancer progression. Cell Immunol 2019;343:103740. [Crossref] [PubMed]
28. Gao K. The effect and mechanism of transgenic over-expression Interleukin-33 in tumor metastastic mouse model. Peking Union Medical College, Chinese Academy of Medical Sciences; 2013.
29. Huang L, Rong Y, Tang X, et al. Engineered exosomes as an in situ DC-primed vaccine to boost antitumor immunity in breast cancer. Mol Cancer 2022;21:45. [Crossref] [PubMed]
30. Lei S, Li S, Xiao W, et al. Azurocidin 1 inhibits the aberrant proliferation of triple-negative breast cancer through the regulation of pyroptosis. Oncol Rep 2023;50:188. [Crossref] [PubMed]
31. Tower H, Ruppert M, Britt K. The Immune Microenvironment of Breast Cancer Progression. Cancers (Basel) 2019;11:1375. [Crossref] [PubMed]
32. Bareche Y, Buisseret L, Gruosso T, et al. Unraveling Triple-Negative Breast Cancer Tumor Microenvironment Heterogeneity: Towards an Optimized Treatment Approach. J Natl Cancer Inst 2020;112:708-19. [Crossref] [PubMed]
33. Emens LA. Breast cancer immunobiology driving immunotherapy: vaccines and immune checkpoint blockade. Expert Rev Anticancer Ther 2012;12:1597-611. [Crossref] [PubMed]
34. Gu-Trantien C, Loi S, Garaud S, et al. CD4⁺ follicular helper T cell infiltration predicts breast cancer survival. J Clin Invest 2013;123:2873-92. [Crossref] [PubMed]
35. Gutiérrez-Melo N, Baumjohann D. T follicular helper cells in cancer. Trends Cancer 2023;9:309-25. [Crossref] [PubMed]
36. Emens LA. Breast Cancer Immunotherapy: Facts and Hopes. Clin Cancer Res 2018;24:511-20. [Crossref] [PubMed]
37. Martins-Cardoso K, Maçao A, Souza JL, et al. TF/PAR2 Signaling Axis Supports the Protumor Effect of Neutrophil Extracellular Traps (NETs) on Human Breast Cancer Cells. Cancers (Basel) 2023;16:5. [Crossref] [PubMed]
38. Yi N, Zhang L, Huang X, et al. Lenvatinib-activated NDUFA4L2/IL33/PADI4 pathway induces neutrophil extracellular traps that inhibit cuproptosis in hepatocellular carcinoma. Cell Oncol (Dordr) 2024; Epub ahead of print. [Crossref]
39. Daigo K, Yamaguchi N, Kawamura T, et al. The proteomic profile of circulating pentraxin 3 (PTX3) complex in sepsis demonstrates the interaction with azurocidin 1 and other components of neutrophil extracellular traps. Mol Cell Proteomics 2012;11:M111.015073.
40. Sprung RW, Zhang Q, Kramer MH, et al. Stabilizing the Proteomes of Acute Myeloid Leukemia Cells: Implications for Cancer Proteomics. Mol Cell Proteomics 2024;23:100716. [Crossref] [PubMed]
41. Hao Y, Hu P, Zhang J. Genomic analysis of the prognostic effect of tumor-associated neutrophil-related genes across 15 solid cancer types: an immune perspective. Ann Transl Med 2020;8:1507. [Crossref] [PubMed]