Research Paper Volume 13, Issue 19 pp 23072—23095

Establishment of a novel ferroptosis-related lncRNA pair prognostic model in colon adenocarcinoma

Hong Li1, *, , Lili Liu2, *, , Tianyi Huang3, *, , Ming Jin3, , Zhen Zheng3, , Hui Zhang3, , Meng Ye4, , Kaitai Liu3, ,

  • 1 Department of General Surgery, The Affiliated Lihuili Hospital, Ningbo University, Ningbo, China
  • 2 Department of Medical Oncology, The Second Affiliated Hospital of Dalian Medical University, Dalian, China
  • 3 Department of Radiation Oncology, The Affiliated Lihuili Hospital, Ningbo University, Ningbo, China
  • 4 Department of Oncology and Hematology, The Affiliated Hospital of Medical School of Ningbo University, Ningbo, China
* Equal contribution

Received: June 14, 2021       Accepted: September 20, 2021       Published: October 5, 2021      

https://doi.org/10.18632/aging.203599
How to Cite

Copyright: © 2021 Li et al. This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Abstract

Long non-coding RNAs (lncRNAs) have been reported to be prognostic factors for cancer. Ferroptosis is an iron-dependent process of programmed cell death. Here, we established a ferroptosis-related lncRNA (frlncRNA) pair signature and revealed its prognostic value in colon adenocarcinoma (COAD) by analyzing the data from The Cancer Genome Atlas (TCGA). FrlncRNAs were identified based on co-expression analysis using the Pearson correlation. Differentially expressed frlncRNAs (DEfrlncRNAs) were recognized and paired, followed by prognostic assessment using univariate Cox regression analysis. The least absolute shrinkage and selection operator (LASSO) penalized Cox analysis was used to determine and construct a risk score prognostic model, by which the receiver operating characteristic (ROC) curves for predicting the overall survival (OS) were conducted. Following the evaluation of whether it was an independent prognostic factor, correlations between the risk score model and clinicopathological characteristics, hypoxia- and immune-related factors, and somatic variants were investigated. In total, 148 DEfrlncRNA pairs were identified, 25 of which were involved in a risk score prognostic signature. The area under ROC curves (AUCs) representing the predictive effect for 1-, 3-, and 5-year survival rates were 0.860, 0.885, and 0.934, respectively. The risk score model was confirmed as an independent prognostic factor and was significantly superior to the clinicopathological characteristics. Correlation analyses showed disparities in clinicopathological characteristics, hypoxia- and immune-related factors, and somatic variants, as well as specific signaling pathways between high- and low-risk groups. The novel risk score prognostic model constructed by pairing DEfrlncRNAs showed promising clinical prediction value in COAD.

Introduction

Colon cancer (CC) is one of the most frequently diagnosed malignant tumors worldwide and is a leading cause of cancer-related deaths in China [12]. Approximately 80%–90% of CC cases are colon adenocarcinoma (COAD) based on the pathologic classification [3]. In general, 20%–25% of CC patients are diagnosed with unresectable metastatic disease and have a poor prognosis [4]. In recent years, the 5-year over survival (OS) rate of CC patients with local stage disease was 90.3%, while that of patients with distant metastases was 12.5% [5]. Thus, developing promising prognostic-related risk assessment signatures is of great clinical significance in the management of colon cancer.

Long non-coding RNAs (lncRNAs), which are a subset of RNAs more than 200 nucleotides in length, account for approximately 80% of the human transcriptome [6]. LncRNAs do not code proteins but exert biological functions by regulating gene expression [7]. Studies have reported that lncRNAs are involved in the malignant phenotypes of various cancers, including proliferation, invasion, and treatment resistance [8]. Evidence has shown abnormal lncRNA expression in tumor samples compared to the corresponding adjacent normal tissues [9, 10]. In addition, studies have demonstrated that lncRNAs contribute to tumorigenesis not only through physiological and biochemical cellular processes, such as ferroptosis, antioxidant capacity, apoptosis, and autophagy, but also by altering the immune microenvironment in cancers [1113]. Ferroptosis is an iron-dependent cell death first proposed by Dixon in 2012 and is characterized by an intracellular iron-dependent accumulation of reactive oxygen species (ROS) and lipid peroxidation [14]. It is a new form of cell death that differs from apoptosis and autophagy [15, 16]. Recent studies have revealed that ferroptosis is a critical factor in metabolism, redox biology, and cell death and has been gradually confirmed as a novel feature for cancer therapy, particularly for those resistant to traditional therapies [17, 18].

Previous studies have reported lncRNA prognostic signatures that have promising predictive and prognostic values in cancer. Tang et al. developed a signature that included 25 differentially expressed ferroptosis-related lncRNAs to predict the prognosis of head and neck squamous cell carcinoma (HNSCC) [19]. The area under the receiver operator characteristic (ROC) curve (AUC) of the lncRNA signature was 0.782, showing a promising prediction value for HNSCC. Moreover, Shen et al. identified and validated an immune-related lncRNA prognostic signature for breast cancer. In this prior study, the prognostic signature comprised 11 lncRNAs and was associated with the infiltration of immune cell subtypes [20]. Furthermore, Soudeh Ghafouri-Fard et al. reviewed the lncRNA signature in gastric cancer based on the available literature and concluded that these transcripts deserve further evaluation as therapeutic targets in gastric cancer [21]. In addition, Wei et al. constructed an autophagy-related lncRNA signature comprising eight lncRNAs and predicted unfavorable prognosis in colorectal cancer. The AUC of this signature was 0.689, indicating that the risk model was effective [22]. Zeng et al. established a differentially expressed lncRNA signature to evaluate the outcome of patients with colorectal cancer. They found that 20 lncRNAs closely related to OS in patients with COAD, four of which were involved in a prognostic model, could serve as an independent factor for survival in COAD [23]. Although the lncRNA signature showed promising prediction performance in cancer prognosis, it was not perfect. The signature required specific expression levels of selected lncRNAs, which should be normalized to reduce the batch effects between different testing platforms before clinical application. Additionally, Hong et al. utilized a novel modeling algorithm, pairing, and iteration to construct an immune-related lncRNA pair prognosis signature in human hepatocellular carcinoma [24]. The AUCs of 1-, 3-, and 5-year survival rates were 0.865, 0.851, and 0.904, respectively. However, there are few reports on lncRNA pair models for cancer prognostic prediction.

In the present study, we developed a novel ferroptosis-related lncRNA pair (frlncRNA pair) prognostic model for COAD. Correlations between the model and clinicopathological characteristics, hypoxia-, immune-related factors, and somatic variants were investigated.

Results

Data characteristics

The expression data of 480 colon adenocarcinoma and 41 adjacent normal samples were included in the present study. The clinical information of patients (n=459) including age, gender, stage, T status, N status, and M status, are shown in Table 1. A total of 259 ferroptosis-related genes (frGenes) were downloaded from the FerrDb website (Supplementary Table 1). Overall, 896 lncRNAs were identified as ferroptosis-related lncRNAs (frlncRNAs) (Supplementary Table 2). Subsequently, 165 (including 158 upregulated and 7 downregulated) of these were identified as differentially expressed ferroptosis-related lncRNAs (DEfrlncRNAs) (Supplementary Table 3), which were visualized using a heatmap (Supplementary Figure 1) and a volcano plot (Supplementary Figure 2).

Table 1. The clinical characteristics of COAD patients in the TCGA dataset.

VariablesNumber of cases
Total459
Age (years)
<60 / ≥60126/333
Gender
Male/Female216/243
Stage
I/II/III/IV/NA76/177/129/65/11
T
T0/T1/T2/T3/T41/11/78/313/56
N
N0/N1/N2270/106/83
M
M0/M1/NA337/65/57

Establishment of a frlncRNA pair risk score prognostic model

To explore a more objective prognostic evaluation model that did not require the specific expression values to be normalized, a 0-or-1–matrix of 7543 frlncRNA pairs was constructed (Supplementary Table 4). According to univariate Cox proportional hazards regression analyses, 148 frlncRNA pairs were considered to be prognostic-associated lncRNA pairs (Supplementary Table 5). In training cohort, after the least absolute shrinkage and selection operator (LASSO) regression analysis, a prognostic signature including 44 frlncRNA pairs was established (Supplementary Table 6). The AUC representing the predictive value of the risk score model for 1-year survival rate in training and validation cohort were 0.904 and 0.728, respectively (Supplementary Figure 3A, 3B). Survival analyses showed significant difference between the high- and low-risk groups both in training and validation cohort (Supplementary Figure 3C, 3D). Subsequently, we constructed a risk score model based on the data of entire cohort (Figure 1). The list of 25 frlncRNA pairs and their corresponding calculation coefficients are shown in Table 2. The AUCs for 1-, 3-, and 5-year survival rates were 0.860, 0.885, and 0.934, respectively (Figure 2A). Furthermore, we identified the maximum inflection point of 1.307 as the optimal cut-off point on the 5-year receiver operator characteristic (ROC) curve (Figure 2B). Moreover, our results showed that the risk score model was significantly superior to the common clinicopathological characteristics, including age, gender, T status, N status, and M status in predicting the OS of COAD patients (Figure 2C).

Establishment of a prognostic model based on (A, B) LASSO regression analysis; (C) Univariate Cox regression analysis.

Figure 1. Establishment of a prognostic model based on (A, B) LASSO regression analysis; (C) Univariate Cox regression analysis.

Table 2. The list of lncRNA pairs and corresponding calculation coefficients.

LncRNA pairCoefficient
LMO7-AS1|LINC005130.493078
LINC01614|AC145423.20.746247
LINC01703|FENDRR0.375638
LINC02487|AC010973.2-0.47787
LINC02195|AC048344.4-0.89094
AC020907.4|AC010973.2-0.67771
MHENCR|AC025857.20.728441
AL031716.1|AL117379.1-0.76209
MIR17HG|AL161729.4-0.64218
AC127024.4|AL355802.30.774316
AC010973.2|AL031673.10.521472
AL021578.1|AL133243.20.426426
AC016831.4|AC011676.1-0.94739
AC090116.1|AL353804.20.49057
AP002336.2|AC093732.1-0.56376
AC011676.1|AC092168.20.835864
AP005899.1|GK-AS10.378791
LINC-PINT|LINC005130.605353
AF117829.1|SNHG22-0.57948
AC245100.7|LINC018110.77997
AL354836.1|SNHG40.551547
SCARNA9|AC104695.4-0.39001
AL137782.1|AP001469.3-0.95082
AP001469.3|ABALON0.598743
ABALON|CD44-AS1-0.64299
(A) The ROC curves for predicting the 1-, 3-, and 5-year OS; (B) identification of the maximum inflection point as the optimal cut-off value on the 5-year ROC curve; (C) comparison of the risk score model and clinicopathological characteristics in predicting the 5-year OS.

Figure 2. (A) The ROC curves for predicting the 1-, 3-, and 5-year OS; (B) identification of the maximum inflection point as the optimal cut-off value on the 5-year ROC curve; (C) comparison of the risk score model and clinicopathological characteristics in predicting the 5-year OS.

Predictive assessment and clinical correlation of the prognostic model

According to the cut-off point recognized previously, 226 patients were classified into the low-risk group and 192 into the high-risk group (Supplementary Table 7). The risk assessment model for prognosis prediction demonstrated that the number of deaths increased with an increase in the risk score (Figure 3A, 3B). Survival analysis revealed that the high-risk group had significantly worse OS than low-risk group (Figure 3C). Age, T status, N status, M status, and risk score model were identified as significant risk factors in the univariate analysis (all P<0.01) (Figure 3D). The risk score model, T status and M status were confirmed as independent prognostic factors by multivariate analysis (all P<0.01) (Figure 3E). Furthermore, our results showed that the risk score model was significantly related to T status, N status, M status, and stage (Figure 4). Moreover, an accurate prognostic nomogram incorporating the risk score model and common clinicopathological characteristics was established for predicting 1-, 3-, and 5-year OS probability, which might be well applied in the clinical evaluation of COAD patients (Figure 5).

Risk scores (A) and survival outcomes (B) of each patient; (C) survival curves of high-risk group and low-risk group patients; (D) univariate and (E) multivariate Cox regression analyses of the risk score model and clinicopathological characteristics.

Figure 3. Risk scores (A) and survival outcomes (B) of each patient; (C) survival curves of high-risk group and low-risk group patients; (D) univariate and (E) multivariate Cox regression analyses of the risk score model and clinicopathological characteristics.

Correlations between the risk score model and clinicopathological characteristics, represented by a heatmap (A), and box diagrams (B).

Figure 4. Correlations between the risk score model and clinicopathological characteristics, represented by a heatmap (A), and box diagrams (B).

Prognostic nomogram incorporating the risk score model and clinicopathological characteristics.

Figure 5. Prognostic nomogram incorporating the risk score model and clinicopathological characteristics.

Functional enrichment analyses of differently expressed frGenes (DEfrGenes)

We identified 142 DEfrgenes in COAD tissues (Supplementary Table 8). Based on these genes, we explored the underlying biological functions by GO annotation and KEGG pathway analyses using Metascape. In this study, GO pathway and process enrichment analysis included molecular functions (functional set), biological processes (pathway), and cellular components (structural complex). The top 20 clusters with their representative enriched terms are shown in Figure 6A. Our results showed that biological processes related to oxygen metabolism were significantly associated with differentially expressed ferroptosis-related genes, including GO:0006979 (response to oxidative stress), GO:0072593 (reactive oxygen species metabolic process), GO:0055114 (oxidation-reduction process), GO:0070482 (response to oxygen levels), GO:0016491 (oxidoreductase activity), and GO:1901615 (organic hydroxy compound metabolic process). In addition, autophagy-related biological processes, such as GO:0006914 (autophagy), GO:0000422 (autophagy of mitochondria), GO:0055072 (iron ion homeostasis), and GO:0000407 (phagophore assembly site) were prominently related to differentially expressed ferroptosis-related genes. Moreover, DEfrGenes in COAD also dramatically affected GO:0097190 (apoptotic signaling pathway), GO:0070997 (neuron death), and several metabolism-related biological processes.

(A) Top 20 enriched clusters of DEfrGenes in GO annotation analysis; (B) top 20 enriched KEGG pathways of DEfrGenes.

Figure 6. (A) Top 20 enriched clusters of DEfrGenes in GO annotation analysis; (B) top 20 enriched KEGG pathways of DEfrGenes.

The top 20 KEGG pathways (P<0.05) for the DEfrGenes were identified (Figure 6B). Classical cancer-related pathways such as “foxo signaling pathway,” “MicroRNAs in cancer,” “necroptosis,” “platinum drug resistance,” “AMPK signaling pathway,” “TGF-beta signaling pathway” and “PPAR signaling pathway” were correlated with the functions of the DEfrGenes in COAD. It is worth noting that “HIF-1 signaling pathway” was the fourth most significantly correlated signaling pathway.

Correlations between the risk score model and hypoxia-related factors

GO and KEGG analysis showed that frGenes were notably related to biological processes of oxygen metabolism. Therefore, we further investigated the relationship between the risk score model and hypoxia-related factors, including hypoxia-inducible factors, and hypoxia-related genes (Supplementary Table 9). The results revealed that the high-risk group was significantly correlated with high expression of ARNT, HIF3A, VEGFA, TUBB6, and low expression of TPI1, MRPS17, LDHA, ENO1, and CDKN3 (Figure 7).

Correlations between the risk score model and hypoxia-related factors.

Figure 7. Correlations between the risk score model and hypoxia-related factors.

Correlations between the risk score model and immune-related factors

In consideration of the increasing evidence on the correlation between immunological features and survival in malignant tumors, we discussed the correlations between the risk score model and immune-related factors. The immune-related factors included tumor-infiltrating immune cells (TIICs) and immune checkpoint genes (ICGs). The 16 ICGs were chosen based on a previous study by Danilova [25] with minor modifications by adding ligands (PVR and NECTIN2) and competing receptor (CD226) of TIGIT (Supplementary Table 10) [26, 27]. We discovered that the low-risk group was related to more TIICs such as CD8+ T cells, CD4+ T cells, B cells, and neutrophils, whereas the high-risk group was related to more tumor-infiltrating immune cells such as NK cells, macrophages, and T cell regulatory (Tregs) (Figure 8A). In addition, the high-risk group had significantly downregulated expression of CD47 and markedly upregulated expression of CD276 and NECTIN2 (Figure 8B).

Correlations between the risk score model and (A) tumor-infiltrating immune cells; (B) immune checkpoint genes.

Figure 8. Correlations between the risk score model and (A) tumor-infiltrating immune cells; (B) immune checkpoint genes.

Correlations between the risk score model and somatic variants

Studies have shown that cancers harboring more non-synonymous variants, which were defined as high tumor mutation burden (TMB), were associated with favorable survival outcomes with cancer immunotherapy. Given the prognostic value of TMB, we sought to investigate the relationship between the risk score model and TMB level. As shown in Figure 9A, there was no significant difference in TMB levels between the high- and low-risk groups. Furthermore, our study found no significant correlation between TMB and prognosis in patients with COAD (Figure 9B). Next, we evaluated the predictive power of the risk score model in the low and high TMB subgroups. The results revealed that the prognostic model presented consistent predictive value in both the low and high TMB subgroups, indicating that the TMB status did not interfere with the prediction of this model (Figure 9C). Moreover, our study explored the mutation rates of reported prognostic-related genes in low-and high-risk groups. Analysis results demonstrated that APC, SMAD4, DOCK2, TMEM132D, and VCAN genes had more frequent mutations in the low-risk group, whereas TP53 and BRAF mutations were more frequent in the high-risk group (Figure 9D, 9E).

Correlations between the risk score model and somatic variants. (A) TMB levels between the high- and low-risk groups; (B) correlation between TMB and prognosis in patients with COAD; (C) prognostic predictive value in different TMB subgroups; (D, E) the mutation rates of reported prognostic-related genes in low- and high-risk groups.

Figure 9. Correlations between the risk score model and somatic variants. (A) TMB levels between the high- and low-risk groups; (B) correlation between TMB and prognosis in patients with COAD; (C) prognostic predictive value in different TMB subgroups; (D, E) the mutation rates of reported prognostic-related genes in low- and high-risk groups.

Gene set enrichment analysis (GSEA)

The GSEA results indicated that the high-risk score group had markedly negative correlations with 26 enrichment pathways (Supplementary Table 11). As showed in Figure 10, the top 10 KEGG pathways included “peroxisome,” “citrate cycle,” “cell cycle,” “oxidative phosphorylation,” “nucleotide excision repair,” and metabolism-related signaling pathways such as “propanoate metabolism,” “valine leucine and isoleucine degradation,” “pyruvate metabolism,” “pyrimidine metabolism,” and so on.

Gene set enrichment analysis of enriched signaling pathways.

Figure 10. Gene set enrichment analysis of enriched signaling pathways.

Discussion

Ferroptosis is an iron-dependent cell death, which is characterized by an intracellular iron-dependent accumulation of ROS and lipid peroxidation. Shen et al. found resibufogenin inhibited colorectal cancer cell growth and tumorigenesis through triggering ferroptosis, suggesting that suppression of ferroptosis is closely related to the proliferation of colorectal cancer cells [28]. Sui and colleagues demonstrated that RSL3-induced cell ferroptosis was relevant to colorectal cancer progression [29]. They revealed that RSL3 can promote CRC cell ferroptosis by promoting ROS production. Sujeong Park et al. reported that bromelain effectively inhibited cell growth and proliferation by stimulating ferroptosis, especially in Kras mutant colorectal cancer cells [30]. The findings indicate that ferroptosis may be involved in carcinogenesis in mutant colorectal cancer cells, compared to Kras wild-type colorectal cancer cells. These studies showed that ferroptosis is crucial in regulating the growth and proliferation of colon cancer cells, and drug activation of ferroptosis related signaling pathways is a promising strategy for treating colon cancer. Moreover, studies have shown that ferroptosis related gene signatures can be effectively used for prognostic prediction, optimizing survival risk assessment and facilitating personalized management in colon cancer [3133]. These results indicate that the research of ferroptosis can provide ideas for the clinical diagnosis and treatment of colon cancer. LncRNA prognostic signatures have been reported to have promising predictive value in cancer. In recent years, prognostic prediction models based on specific classes of gene-related lncRNAs have attracted particular attention from researchers. However, most of these prognostic models have limited predictive efficacy. For example, the AUC was 0.782 for a signature that included 25 differentially expressed ferroptosis-related lncRNAs in predicting the prognosis of HNSCC [19]. An autophagy-related lncRNA signature had an AUC value of 0.689, which comprised eight lncRNAs and predicted an unfavorable prognosis in colorectal cancer [22]. A prognostic model including four differentially expressed lncRNAs had an AUC value of 0.706 for evaluating the outcome of patients with colorectal cancer [23]. In addition to unsatisfactory predictive performance, these models also have practical shortcomings. These prognostic models were established based on the specific expression values of the identified lncRNAs. The measured values must be normalized to reduce batch effects between different testing platforms before clinical application. In the present study, we developed a novel ferroptosis-related lncRNA pair prognostic signature for colon adenocarcinoma. The AUCs representing the excellent predictive value for 1-, 3-, and 5-year survival rates were 0.860, 0.885, and 0.934, respectively. Our prognostic signature included 25 frlncRNA pairs and was confirmed as an independent prognostic factor by Cox regression analysis. Notably, it was superior to common clinicopathological characteristics such as age, sex, T status, N status, and M status in predicting the OS for COAD. More importantly, the signature was established by a novel modeling algorithm, pairing, and iteration; therefore, it could be better applied in clinical practice.

To explore the underlying biological functions of our signature, we performed GO annotation and KEGG pathway analyses of the identified DEfrGenes. Both function and pathway enrichment analyses highlighted the hypoxia and oxygen metabolism processes, such as “HIF-1 signaling pathway,” “response to oxidative stress,” “reactive oxygen species metabolic process,” “oxidation-reduction process,” “response to oxygen levels,” “oxidoreductase activity” and “organic hydroxy compound metabolic process.” Studies have found that hypoxia is a poor prognostic factor that regulates the tumor microenvironment. Hypoxia activates a series of signaling pathways and a large panel of genes involved in apoptosis, autophagy, DNA damage, mitochondrial activity, p53, and drug efflux, which affects the survival of cancer cells and confers resistance to conventional anticancer treatments [34, 35]. Preclinical studies have shown that inhibition of HIF-1 signaling activity can significantly reduce cancer growth, and research is currently underway to identify HIF-1 inhibitors and validate their efficacy in antitumor therapy [36]. Several hypoxia-related signatures have been constructed for the clinical prediction of diagnosis, prognosis, and recurrence in hepatocellular carcinoma. These hypoxia-related signatures were revealed to be closely associated with the immune microenvironment, contributing to better clinical management in patients with hepatocellular carcinoma [3739].

For the hypoxia and oxygen metabolism processes, as indicated by both the GO and KEGG analyses, we further explored the relationship between the risk score model and hypoxia-inducible factors and hypoxia-related genes. The results revealed that the high-risk group exhibited a significant correlation with high expression of ARNT, HIF3A, VEGFA, TUBB6 and low expression of TPI1, MRPS17, LDHA, ENO1, and CDKN3. Studies have demonstrated that high levels of HIF-1/2α proteins are associated with increased radiotherapy and chemotherapy resistance, resulting in cancer progression [40, 41]. In addition to inducing hypoxia, the HIF signaling pathway also participates in various bioregulatory processes, including ROS, cytokines, growth factors, oncogene dysfunction, and tumor suppressors, providing insights into the targeting of the HIF signaling pathway as a potential therapeutic approach in cancer management [42]. Studies have reported that hypoxia can promote disease development and metastasis by activating the HIF1α/VEGFA pathway in breast and colorectal cancer [43, 44]. Liang and colleagues found that hypoxia activated the HIF1α/VEGFA axis in breast cancer angiogenesis by inducing miR-153 expression and decreased expression of HIF1α and VEGFA, resulting in suppression of tumor angiogenesis. Chen et al. demonstrated that hypoxia could induce angiogenesis in colorectal cancer by activating the HIF-1α/VEGF-A pathway. LDHA and ENO1 are crucial glycolytic enzymes, which have been revealed as hypoxia-related factors in previous studies. Wei et al. found that inhibition of LDHA-mediated aerobic glycolysis could markedly suppress the growth of bladder cancer cells [45]. Cui et al. reported that HIF1/2α could activate LDHA expression, and high expression of LDHA promotes the growth and migration of pancreatic cancer cells [46]. Nevertheless, Wang et al. reported that high LDHA expression was associated with better PFS (18.3 vs. 10.1 months, P=0.005) and overall response rate (72.2% vs. 15.4%, P=0.006) of metastatic colorectal cancer patients receiving first-line chemotherapy [47]. Wang et al. discovered the role of ENO1 in hypoxia-induced gemcitabine chemoresistance, decreased expression of which restored sensitivity to chemotherapy by modulating redox homeostasis in pancreatic cancer cells [48]. However, our results showed lower expression of LDHA and ENO1 in the high-risk score group, which suggests some unknown underlying signal regulation mechanisms, which are worth further exploration.

The immune microenvironment plays an important role in tumorigenesis. Infiltrating immune cells may act as tumor-antagonizing or tumor-promoting factors [49, 50]. Cancer cells eventually gain the ability to inhibit the tumor-antagonizing functions of immune cells and escape immunological surveillance, resulting in cancer development [51]. In recent years, immunotherapy targeting immune microenvironment regulation and immune checkpoint modulation has shown promising efficacy in cancer treatment [52, 53]. In addition, iron and immunity are closely linked [54]. Many of the genes and proteins involved in iron homoeostasis play a vital role in controlling iron fluxes. Cells of the innate immune system, monocytes, macrophages, microglia and lymphocytes, are able to combat bacterial insults by carefully controlling their iron fluxes, which are mediated by hepcidin and ferroportin. The cell involved in iron overload with the greatest effect on immunity is the macrophage. Intriguing evidence indicated that parenchymal iron overload is linked to genes classically associated with the immune system [55]. In the present study, we investigated the correlation between the risk score model and immune-related factors. The results showed that the low-risk group was related to more tumor-infiltrating immune cells such as CD8+ T cells, CD4+ T cells, B cells, and neutrophils, whereas the high-risk group was related to more tumor-infiltrating immune cells such as NK cells, macrophages, and Tregs. Moreover, our analyses showed that the high-risk group had significantly decreased expression of CD47 and markedly upregulated expression of CD276 and NECTIN2. Studies have shown that CD47 plays a key role in immune regulation and tumor development. Overexpression of CD47 could protect tumor cells from phagocytosis and is a promising therapeutic target in cancer therapy [56, 57]. CD276, also known as B7-H3, is an important immune checkpoint member of the B7 and CD28 families. CD276 was found to be overexpressed in various tumor cells, which acts as an inhibitor of T cell function and indicates poor prognosis in cancer patients [58, 59]. Studies have shown that NECTIN2 is an adhesion molecule that participates in tumor growth, metastasis, and tumor angiogenesis [60]. There were significant differences in tumor-infiltrating immune cells and ICGs between the high-risk and low-risk groups in our prognostic signature, which may serve as potential molecular markers for predicting the efficacy of immunotherapy.

TMB is not only a prognostic predictor, but also a predictor of the efficacy of immunotherapy. TMB has been demonstrated to be a useful biomarker for predicting the efficacy of immune checkpoint inhibitors in some cancer types [61]. Studies have shown a positive association between TMB and the response to immunotherapy in melanoma and NSCLC [62, 63]. Cao et al. found that patients with high TMB may benefit more from immunotherapy and experience better survival time in various cancer types [64]. However, the predictive role of TMB in the prognosis and efficacy of immunotherapy in colon cancer remains controversial. Lee et al. revealed that high TMB indicated better outcomes for patients treated with curative surgery followed by adjuvant fluoropyrimidine and oxaliplatin chemotherapy [65]. In contrast, Pai et al. found that the low TMB group had improved progression-free survival time (9.9 months) compared to the intermediate/high TMB group (5.8 months) in metastatic colorectal cancer patients treated with first-line chemotherapy [66]. Zhou et al. constructed a nomogram model for prognosis prediction in CC patients, which combined TMB profiles, immunocyte infiltration status, and clinicopathological data. The results showed that patients with low TMB had better outcomes [67]. Indeed, the prognostic value of TMB and its predictive effect on immunotherapy are controversial. Recently, a study found that TMB could only predict response to immunotherapy in cancer types where CD8 T cell levels positively correlated with neoantigen load, such as melanoma, lung, and bladder cancers [68]. The present study showed that TMB was not a significant prognostic factor for CC. There was no significant difference in TMB levels between the high- and low-risk groups. Furthermore, our study suggests that the prognostic role of the signature is independent of TMB levels. We further discovered that the high-risk group had more frequent mutations in TP53 and BRAF, as well as less frequent mutations in APC, SMAD4, DOCK2, TMEM132D, and VCAN genes. Studies have shown that TP53 mutations may accelerate the progression of CRC by activating oncogenic and inflammatory pathways [69]. BRAF mutation was an independent factor of recurrence in microsatellite-stabilized colon cancer, and the outcome of CC patients with BRAF gene mutations was significantly poor [70]. Studies have reported that APC mutation is the most common mutation, accounting for 60% of gene mutations in CC [71]. The current results confirm the role of these genes in predicting prognosis in patients with colon cancer, which is consistent with a previous study [72].

The present study establishes a novel ferroptosis-related lncRNA pair prognostic model that may provide better management in patients with COAD. However, some limitations still need to be considered. First, we only draw conclusions from bioinformatics analysis, and further experimental verification is required. Second, due to the small sample size, there are some statistical defects, for example, the hazard ratio of this the risk score model has a big error range in univariate and multivariate analysis.

Conclusions

The novel risk score model constructed by pairing DEfrlncRNAs showed promising clinical prediction value in COAD, which is worthy of further investigation.

Materials and Methods

Data collection

High-throughput sequencing (HTSeq) transcriptome profiling harmonized to fragments per kilobase million (FPKM) and simple nucleotide variation data of COAD as well as clinical characteristics of patients were obtained from The Cancer Genome Atlas (TCGA) database (https://portal.gdc.cancer.gov). Cases with survival data were randomly divided into the training cohort and validation cohort in 2:1 ratio. The GTF file of long non-coding RNA gene annotation (version GRCh38.p13) was downloaded from the GENCODE human website (https://www.gencodegenes.org/human/). The ferroptosis-related genes involved in the current study were obtained from the website of FerrDb, a public database of ferroptosis regulators and markers, and ferroptosis-disease associations (http://www.zhounan.org/ferrdb/) [73].

Identification of frlncRNAs and pairing DEfrlncRNAs

Identification of the frlncRNAs was performed using Pearson correlation to assess the relationship between the ferroptosis-related genes and long non-coding RNAs. The absolute value of correlation coefficients >0.4 and P<0.001 were considered statistically significant. The significant thresholds of DEfrlncRNAs were set as |log2FC| >1.5 and false discovery rate (FDR) <0.001. Subsequently, we established frlncRNA pairs based on these DEfrlncRNAs, as previously described [24]. All the DEfrlncRNAs were cyclically paired and assigned a value according to pairwise comparison in accordance with the following rules: assume that lncRNA A and lncRNA B were paired together as lncRNA pair C, which was assigned a value of 1 if the expression level of lncRNA A was higher than lncRNA B; otherwise, it was assigned a value of 0. If a lncRNA pair had a 0 or 1 ratio of less than 20% or greater than 80% in all samples, it was filtered. Cox regression analyses were performed to evaluate the prognostic value of frlncRNA pairs (P<0.01).

Construction and validation of a frlncRNA pair prognostic signature

Prognostically associated frlncRNA pairs were used to establish the LASSO regression model in the training cohort. Subsequently, a risk score model of these frlncRNA pairs was constructed and the risk score of each patient was calculated according to the following formula: Risk score = coefficient lncRNA pair1 × expression lncRNA pair1 + coefficient lncRNA pair2 × expression lncRNA pair2 + coefficient lncRNA pair3 × expression lncRNA pair3 +……+ coefficient lncRNA pairn × expression lncRNA pairn. The 1-year ROC curve for predicting the OS was constructed according to the risk score model. Using the median value of risk scores, patients in training- and validation cohort were divided into high- and low-risk groups. Survival curves were conducted using the Kaplan-Meier method. The validation cohort was applied for cross test to assess the stability of this model. In order to obtain a more accurate model with a larger sample size, we further constructed the final model based on the data of entire cohort. The 1-, 3-, and 5-year ROC curves were conducted. The maximum inflection point of the 5-year ROC curve was considered as the optimal cut-off point for the classification of different risk groups. A risk assessment model for prognosis prediction was generated to show the risk scores and survival outcomes of each patient. Moreover, univariate and multivariate regression analyses were used to evaluate whether the model was an independent prognostic factor for OS in COAD patients.

Clinical correlation analysis of the risk score model and development of a nomogram

The correlations between the risk score model and traditional clinicopathological features were assessed using the chi-squared test and presented as a heatmap. The risk score differences among different groups of these clinicopathological features were calculated using the Wilcoxon signed-rank test and visualized using box diagrams. The p value was labeled as follows: P<0.001 = ***, P<0.01 = ** *, and P<0.05 = *. The multivariate logistic model, including the risk score model and clinicopathological characteristics, was used to develop a nomogram to predict the survival probability of COAD patients. The 1-, 3-, and 5-year OS rates were used as the endpoints in the nomogram.

Correlations of the risk score model with hypoxia-related and immune-related factors

DEfrGenes in COAD tissues were identified using the R package limma. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses of DEfrGenes were performed using Metascape (http://metascape.org) with thresholds of min overlap 3, P<0.05, and min enrichment 3 [74]. Hypoxia-inducible factors included HIF1A, ARNT, EPAS1, ARNT2, HIF3A, and ARNTL. Hypoxia-related genes were retrieved from previous reports [75, 76]. Differential gene expression between the two risk score groups was conducted using the Wilcoxon signed-rank test. We employed several currently acknowledged algorithms to evaluate the correlations between the risk score and TIICs, including XCELL, TIMER, QUANTISEQ, MCPCOUNTER, EPIC, CIBERSORT−ABS, and CIBERSORT. Associations between the risk score groups and the expression of ICGs were assessed, and the results were presented as violin plots. The statistical significance was set at P<0.05.

Correlations between the risk score model and somatic variants

To evaluate the TMB, we counted the total number of non-synonymous mutations in COAD. TMB was defined as the number of somatic, coding, indel mutations, and base substitution per megabase (Mb) of the genome examined. To calculate the TMB score of each sample, the total number of mutations counted was divided by the exome size (38 Mb was used as the estimate of the exome size) [77]. The somatic variant status of reported prognostic-related genes in COAD was assessed in low- and high-risk groups. The included genes were APC, TP53, KRAS, NRAS, BRAF, FAT4, SMAD4, COL6A3, CDH10, DOCK2, TMEM132D, and VCANT, which were retrieved from a previous study [72].

Gene set enrichment analysis

Gene set enrichment analysis (high-risk score group vs. low-risk score group) (version 4.1.0; http://software.broadinstitute.org/gsea/index.jsp) was performed to investigate the potential molecular mechanisms by which the risk score model might act on tumor progression in COAD, as previously described [78, 79]. FDR <0.05 was used to identify significantly enriched pathways.

Statistical analysis

Perl software (version 5.32) was used to extract and structure the HTSeq FPKM and simple nucleotide variation data. The differentially expressed lncRNAs were identified using the Benjamini-Hochberg method based on the log fold change and FDR. Survival analyses of COAD patients based on the risk score model were assessed using the Kaplan-Meier method. Multivariate analysis was performed using the Cox regression model. R software version 4.0.3, with Bioconductor packages, was used to conduct the analyses.

Data statement

The data used for bioinformatics analyses in this study are freely available on the National Cancer Institute Genomic Data Commons (GDC) Data Portal (https://portal.gdc.cancer.gov/). The interpretation and reporting of these data are the sole responsibility of the authors.

Abbreviations

LncRNAs: long non-coding RNAs; FrlncRNA: ferroptosis-related lncRNA; COAD: colon adenocarcinoma; TCGA: The Cancer Genome Atlas; LASSO: least absolute shrinkage and selection operator; ROC: receiver operating characteristic; AUCs: area under ROC curves; OS: overall survival; HNSCC: head and neck squamous cell carcinoma; FrGenes: ferroptosis-related genes; TIICs: tumor immune infiltrating cells; ICGs: immune checkpoint genes; TMB: tumor mutation burden; Mb: megabase; HTSeq: high-throughput sequencing; FPKM: fragments per kilobase of transcript per million mapped reads; FDR: false discovery rate; GO: Gene Ontology; KEGG: Kyoto Encyclopedia of Genes and Genomes; GSEA: Gene set enrichment analysis.

Author Contributions

KTL, HL, LLL conceived and designed the study. KTL, LLL, TYH, MJ, ZZ performed the analyses. TYH, MJ, HZ, MY prepared all tables and figures. KTL, MJ and HL wrote the main manuscript. KTL, TYH, LLL revised the manuscript. All authors approved the final version of the manuscript.

Acknowledgments

We acknowledge the TCGA project group and all the R programming package developers.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Funding

This work was supported by the Natural Science Foundation of Zhejiang Province (No. LY21H160013); the Ningbo Clinical Medicine Research Center Project (No. 2019A21003); the Ningbo Natural Science Foundation (No.202003N4203); the Medical Science and Technology Project of Zhejiang Provincial Health Commission (No. 2020KY861).

References

  • 1. Siegel RL, Miller KD, Fuchs HE, Jemal A. Cancer Statistics, 2021. CA Cancer J Clin. 2021; 71:7–33. https://doi.org/10.3322/caac.21654 [PubMed]
  • 2. Chen W, Zheng R, Baade PD, Zhang S, Zeng H, Bray F, Jemal A, Yu XQ, He J. Cancer statistics in China, 2015. CA Cancer J Clin. 2016; 66:115–32. https://doi.org/10.3322/caac.21338 [PubMed]
  • 3. Barresi V, Reggiani Bonetti L, Ieni A, Caruso RA, Tuccari G. Histological grading in colorectal cancer: new insights and perspectives. Histol Histopathol. 2015; 30:1059–67. https://doi.org/10.14670/HH-11-633 [PubMed]
  • 4. Eadens MJ, Grothey A. Curable metastatic colorectal cancer. Curr Oncol Rep. 2011; 13:168–76. https://doi.org/10.1007/s11912-011-0157-0 [PubMed]
  • 5. Pita-Fernández S, González-Sáez L, López-Calviño B, Seoane-Pillado T, Rodríguez-Camacho E, Pazos-Sierra A, González-Santamaría P, Pértega-Díaz S. Effect of diagnostic delay on survival in patients with colorectal cancer: a retrospective cohort study. BMC Cancer. 2016; 16:664. https://doi.org/10.1186/s12885-016-2717-z [PubMed]
  • 6. Mercer TR, Dinger ME, Mattick JS. Long non-coding RNAs: insights into functions. Nat Rev Genet. 2009; 10:155–59. https://doi.org/10.1038/nrg2521 [PubMed]
  • 7. Statello L, Guo CJ, Chen LL, Huarte M. Gene regulation by long non-coding RNAs and its biological functions. Nat Rev Mol Cell Biol. 2021; 22:96–118. https://doi.org/10.1038/s41580-020-00315-9 [PubMed]
  • 8. Ulitsky I, Bartel DP. lincRNAs: genomics, evolution, and mechanisms. Cell. 2013; 154:26–46. https://doi.org/10.1016/j.cell.2013.06.020 [PubMed]
  • 9. Gao N, Li Y, Li J, Gao Z, Yang Z, Li Y, Liu H, Fan T. Long Non-Coding RNAs: The Regulatory Mechanisms, Research Strategies, and Future Directions in Cancers. Front Oncol. 2020; 10:598817. https://doi.org/10.3389/fonc.2020.598817 [PubMed]
  • 10. Guo XB, Hua Z, Li C, Peng LP, Wang JS, Wang B, Zhi QM. Biological significance of long non-coding RNA FTX expression in human colorectal cancer. Int J Clin Exp Med. 2015; 8:15591–600. [PubMed]
  • 11. Gupta RA, Shah N, Wang KC, Kim J, Horlings HM, Wong DJ, Tsai MC, Hung T, Argani P, Rinn JL, Wang Y, Brzoska P, Kong B, et al. Long non-coding RNA HOTAIR reprograms chromatin state to promote cancer metastasis. Nature. 2010; 464:1071–76. https://doi.org/10.1038/nature08975 [PubMed]
  • 12. Atianand MK, Caffrey DR, Fitzgerald KA. Immunobiology of Long Noncoding RNAs. Annu Rev Immunol. 2017; 35:177–98. https://doi.org/10.1146/annurev-immunol-041015-055459 [PubMed]
  • 13. Chen YG, Satpathy AT, Chang HY. Gene regulation in the immune system by long noncoding RNAs. Nat Immunol. 2017; 18:962–72. https://doi.org/10.1038/ni.3771 [PubMed]
  • 14. Dixon SJ. Ferroptosis: bug or feature? Immunol Rev. 2017; 277:150–57. https://doi.org/10.1111/imr.12533 [PubMed]
  • 15. Hassannia B, Vandenabeele P, Vanden Berghe T. Targeting Ferroptosis to Iron Out Cancer. Cancer Cell. 2019; 35:830–49. https://doi.org/10.1016/j.ccell.2019.04.002 [PubMed]
  • 16. Huang L, McClatchy DB, Maher P, Liang Z, Diedrich JK, Soriano-Castell D, Goldberg J, Shokhirev M, Yates JR 3rd, Schubert D, Currais A. Intracellular amyloid toxicity induces oxytosis/ferroptosis regulated cell death. Cell Death Dis. 2020; 11:828. https://doi.org/10.1038/s41419-020-03020-9 [PubMed]
  • 17. Mou Y, Wang J, Wu J, He D, Zhang C, Duan C, Li B. Ferroptosis, a new form of cell death: opportunities and challenges in cancer. J Hematol Oncol. 2019; 12:34. https://doi.org/10.1186/s13045-019-0720-y [PubMed]
  • 18. Xu T, Ding W, Ji X, Ao X, Liu Y, Yu W, Wang J. Molecular mechanisms of ferroptosis and its role in cancer therapy. J Cell Mol Med. 2019; 23:4900–12. https://doi.org/10.1111/jcmm.14511 [PubMed]
  • 19. Tang Y, Li C, Zhang YJ, Wu ZH. Ferroptosis-Related Long Non-Coding RNA signature predicts the prognosis of Head and neck squamous cell carcinoma. Int J Biol Sci. 2021; 17:702–11. https://doi.org/10.7150/ijbs.55552 [PubMed]
  • 20. Shen Y, Peng X, Shen C. Identification and validation of immune-related lncRNA prognostic signature for breast cancer. Genomics. 2020; 112:2640–46. https://doi.org/10.1016/j.ygeno.2020.02.015 [PubMed]
  • 21. Ghafouri-Fard S, Taheri M. Long non-coding RNA signature in gastric cancer. Exp Mol Pathol. 2020; 113:104365. https://doi.org/10.1016/j.yexmp.2019.104365 [PubMed]
  • 22. Wei J, Ge X, Tang Y, Qian Y, Lu W, Jiang K, Fang Y, Hwang M, Fu D, Xiao Q, Ding K. An Autophagy-Related Long Noncoding RNA Signature Contributes to Poor Prognosis in Colorectal Cancer. J Oncol. 2020; 2020:4728947. https://doi.org/10.1155/2020/4728947 [PubMed]
  • 23. Zeng JH, Liang L, He RQ, Tang RX, Cai XY, Chen JQ, Luo DZ, Chen G. Comprehensive investigation of a novel differentially expressed lncRNA expression profile signature to assess the survival of patients with colorectal adenocarcinoma. Oncotarget. 2017; 8:16811–28. https://doi.org/10.18632/oncotarget.15161 [PubMed]
  • 24. Hong W, Liang L, Gu Y, Qi Z, Qiu H, Yang X, Zeng W, Ma L, Xie J. Immune-Related lncRNA to Construct Novel Signature and Predict the Immune Landscape of Human Hepatocellular Carcinoma. Mol Ther Nucleic Acids. 2020; 22:937–47. https://doi.org/10.1016/j.omtn.2020.10.002 [PubMed]
  • 25. Danilova L, Ho WJ, Zhu Q, Vithayathil T, De Jesus-Acosta A, Azad NS, Laheru DA, Fertig EJ, Anders R, Jaffee EM, Yarchoan M. Programmed Cell Death Ligand-1 (PD-L1) and CD8 Expression Profiling Identify an Immunologic Subtype of Pancreatic Ductal Adenocarcinomas with Favorable Survival. Cancer Immunol Res. 2019; 7:886–95. https://doi.org/10.1158/2326-6066.CIR-18-0822 [PubMed]
  • 26. Manieri NA, Chiang EY, Grogan JL. TIGIT: A Key Inhibitor of the Cancer Immunity Cycle. Trends Immunol. 2017; 38:20–28. https://doi.org/10.1016/j.it.2016.10.002 [PubMed]
  • 27. Fourcade J, Sun Z, Chauvin JM, Ka M, Davar D, Pagliano O, Wang H, Saada S, Menna C, Amin R, Sander C, Kirkwood JM, Korman AJ, Zarour HM. CD226 opposes TIGIT to disrupt Tregs in melanoma. JCI Insight. 2018; 3:e121157. https://doi.org/10.1172/jci.insight.121157 [PubMed]
  • 28. Shen LD, Qi WH, Bai JJ, Zuo CY, Bai DL, Gao WD, Zong XL, Hao TT, Ma Y, Cao GC. Resibufogenin inhibited colorectal cancer cell growth and tumorigenesis through triggering ferroptosis and ROS production mediated by GPX4 inactivation. Anat Rec (Hoboken). 2021; 304:313–22. https://doi.org/10.1002/ar.24378 [PubMed]
  • 29. Sui X, Zhang R, Liu S, Duan T, Zhai L, Zhang M, Han X, Xiang Y, Huang X, Lin H, Xie T. RSL3 Drives Ferroptosis Through GPX4 Inactivation and ROS Production in Colorectal Cancer. Front Pharmacol. 2018; 9:1371. https://doi.org/10.3389/fphar.2018.01371 [PubMed]
  • 30. Park S, Oh J, Kim M, Jin EJ. Bromelain effectively suppresses Kras-mutant colorectal cancer by stimulating ferroptosis. Anim Cells Syst (Seoul). 2018; 22:334–40. https://doi.org/10.1080/19768354.2018.1512521 [PubMed]
  • 31. Qi X, Wang R, Lin Y, Yan D, Zuo J, Chen J, Shen B. A Ferroptosis-Related Gene Signature Identified as a Novel Prognostic Biomarker for Colon Cancer. Front Genet. 2021; 12:692426. https://doi.org/10.3389/fgene.2021.692426 [PubMed]
  • 32. Nie J, Shan D, Li S, Zhang S, Zi X, Xing F, Shi J, Liu C, Wang T, Sun X, Zhang Q, Zhou M, Luo S, et al. A Novel Ferroptosis Related Gene Signature for Prognosis Prediction in Patients With Colon Cancer. Front Oncol. 2021; 11:654076. https://doi.org/10.3389/fonc.2021.654076 [PubMed]
  • 33. Wang Y, Xia HB, Chen ZM, Meng L, Xu AM. Identification of a ferroptosis-related gene signature predictive model in colon cancer. World J Surg Oncol. 2021; 19:135. https://doi.org/10.1186/s12957-021-02244-z [PubMed]
  • 34. Jing X, Yang F, Shao C, Wei K, Xie M, Shen H, Shu Y. Role of hypoxia in cancer therapy by regulating the tumor microenvironment. Mol Cancer. 2019; 18:157. https://doi.org/10.1186/s12943-019-1089-9 [PubMed]
  • 35. Brahimi-Horn MC, Chiche J, Pouysségur J. Hypoxia and cancer. J Mol Med (Berl). 2007; 85:1301–07. https://doi.org/10.1007/s00109-007-0281-3 [PubMed]
  • 36. Semenza GL. Targeting HIF-1 for cancer therapy. Nat Rev Cancer. 2003; 3:721–32. https://doi.org/10.1038/nrc1187 [PubMed]
  • 37. Zhang Q, Qiao L, Liao J, Liu Q, Liu P, Liu L. A novel hypoxia gene signature indicates prognosis and immune microenvironments characters in patients with hepatocellular carcinoma. J Cell Mol Med. 2021; 25:3772–84. https://doi.org/10.1111/jcmm.16249 [PubMed]
  • 38. Zhang B, Tang B, Gao J, Li J, Kong L, Qin L. A hypoxia-related signature for clinically predicting diagnosis, prognosis and immune microenvironment of hepatocellular carcinoma patients. J Transl Med. 2020; 18:342. https://doi.org/10.1186/s12967-020-02492-9 [PubMed]
  • 39. Hu B, Yang XB, Sang XT. Development and Verification of the Hypoxia-Related and Immune-Associated Prognosis Signature for Hepatocellular Carcinoma. J Hepatocell Carcinoma. 2020; 7:315–30. https://doi.org/10.2147/JHC.S272109 [PubMed]
  • 40. Keysar SB, Trncic N, Larue SM, Fox MH. Hypoxia/reoxygenation-induced mutations in mammalian cells detected by the flow cytometry mutation assay and characterized by mutant spectrum. Radiat Res. 2010; 173:21–26. https://doi.org/10.1667/RR1838.1 [PubMed]
  • 41. Samanta D, Gilkes DM, Chaturvedi P, Xiang L, Semenza GL. Hypoxia-inducible factors are required for chemotherapy resistance of breast cancer stem cells. Proc Natl Acad Sci USA. 2014; 111:E5429–38. https://doi.org/10.1073/pnas.1421438111 [PubMed]
  • 42. Albadari N, Deng S, Li W. The transcriptional factors HIF-1 and HIF-2 and their novel inhibitors in cancer therapy. Expert Opin Drug Discov. 2019; 14:667–82. https://doi.org/10.1080/17460441.2019.1613370 [PubMed]
  • 43. Liang H, Xiao J, Zhou Z, Wu J, Ge F, Li Z, Zhang H, Sun J, Li F, Liu R, Chen C. Hypoxia induces miR-153 through the IRE1α-XBP1 pathway to fine tune the HIF1α/VEGFA axis in breast cancer angiogenesis. Oncogene. 2018; 37:1961–75. https://doi.org/10.1038/s41388-017-0089-8 [PubMed]
  • 44. Chen H, Feng J, Zhang Y, Shen A, Chen Y, Lin J, Lin W, Sferra TJ, Peng J. Pien Tze Huang Inhibits Hypoxia-Induced Angiogenesis via HIF-1 α /VEGF-A Pathway in Colorectal Cancer. Evid Based Complement Alternat Med. 2015; 2015:454279. https://doi.org/10.1155/2015/454279 [PubMed]
  • 45. Wei Y, Zhang Y, Meng Q, Cui L, Xu C. Hypoxia-induced circular RNA has_circRNA_403658 promotes bladder cancer cell growth through activation of LDHA. Am J Transl Res. 2019; 11:6838–49. [PubMed]
  • 46. Cui XG, Han ZT, He SH, Wu XD, Chen TR, Shao CH, Chen DL, Su N, Chen YM, Wang T, Wang J, Song DW, Yan WJ, et al. HIF1/2α mediates hypoxia-induced LDHA expression in human pancreatic cancer cells. Oncotarget. 2017; 8:24840–52. https://doi.org/10.18632/oncotarget.15266 [PubMed]
  • 47. Wang H, Peng R, Chen X, Jia R, Huang C, Huang Y, Xia L, Guo G. Effect of HK2, PKM2 and LDHA on Cetuximab efficacy in metastatic colorectal cancer. Oncol Lett. 2018; 15:5553–60. https://doi.org/10.3892/ol.2018.8005 [PubMed]
  • 48. Wang L, Bi R, Yin H, Liu H, Li L. ENO1 silencing impaires hypoxia-induced gemcitabine chemoresistance associated with redox modulation in pancreatic cancer cells. Am J Transl Res. 2019; 11:4470–80. [PubMed]
  • 49. Lei X, Lei Y, Li JK, Du WX, Li RG, Yang J, Li J, Li F, Tan HB. Immune cells within the tumor microenvironment: Biological functions and roles in cancer immunotherapy. Cancer Lett. 2020; 470:126–33. https://doi.org/10.1016/j.canlet.2019.11.009 [PubMed]
  • 50. Pauken KE, Wherry EJ. SnapShot: T Cell Exhaustion. Cell. 2015; 163:1038–38.e1. https://doi.org/10.1016/j.cell.2015.10.054 [PubMed]
  • 51. Helmy KY, Patel SA, Nahas GR, Rameshwar P. Cancer immunotherapy: accomplishments to date and future promise. Ther Deliv. 2013; 4:1307–20. https://doi.org/10.4155/tde.13.88 [PubMed]
  • 52. Hinshaw DC, Shevde LA. The Tumor Microenvironment Innately Modulates Cancer Progression. Cancer Res. 2019; 79:4557–66. https://doi.org/10.1158/0008-5472.CAN-18-3962 [PubMed]
  • 53. Iwai Y, Ishida M, Tanaka Y, Okazaki T, Honjo T, Minato N. Involvement of PD-L1 on tumor cells in the escape from host immune system and tumor immunotherapy by PD-L1 blockade. Proc Natl Acad Sci USA. 2002; 99:12293–97. https://doi.org/10.1073/pnas.192461099 [PubMed]
  • 54. Ward RJ, Crichton RR, Taylor DL, Della Corte L, Srai SK, Dexter DT. Iron and the immune system. J Neural Transm (Vienna). 2011; 118:315–28. https://doi.org/10.1007/s00702-010-0479-3 [PubMed]
  • 55. Porto G, De Sousa M. Iron overload and immunity. World J Gastroenterol. 2007; 13:4707–15. https://doi.org/10.3748/wjg.v13.i35.4707 [PubMed]
  • 56. Logtenberg ME, Scheeren FA, Schumacher TN. The CD47-SIRPα Immune Checkpoint. Immunity. 2020; 52:742–52. https://doi.org/10.1016/j.immuni.2020.04.011 [PubMed]
  • 57. Hayat SM, Bianconi V, Pirro M, Jaafari MR, Hatamipour M, Sahebkar A. CD47: role in the immune system and application to cancer therapy. Cell Oncol (Dordr). 2020; 43:19–30. https://doi.org/10.1007/s13402-019-00469-5 [PubMed]
  • 58. Bao R, Wang Y, Lai J, Zhu H, Zhao Y, Li S, Li N, Huang J, Yang Z, Wang F, Liu Z. Enhancing Anti-PD-1/PD-L1 Immune Checkpoint Inhibitory Cancer Therapy by CD276-Targeted Photodynamic Ablation of Tumor Cells and Tumor Vasculature. Mol Pharm. 2019; 16:339–48. https://doi.org/10.1021/acs.molpharmaceut.8b00997 [PubMed]
  • 59. Picarda E, Ohaegbulam KC, Zang X. Molecular Pathways: Targeting B7-H3 (CD276) for Human Cancer Immunotherapy. Clin Cancer Res. 2016; 22:3425–31. https://doi.org/10.1158/1078-0432.CCR-15-2428 [PubMed]
  • 60. Bekes I, Löb S, Holzheu I, Janni W, Baumann L, Wöckel A, Wulff C. Nectin-2 in ovarian cancer: How is it expressed and what might be its functional role? Cancer Sci. 2019; 110:1872–82. https://doi.org/10.1111/cas.13992 [PubMed]
  • 61. Chan TA, Yarchoan M, Jaffee E, Swanton C, Quezada SA, Stenzinger A, Peters S. Development of tumor mutation burden as an immunotherapy biomarker: utility for the oncology clinic. Ann Oncol. 2019; 30:44–56. https://doi.org/10.1093/annonc/mdy495 [PubMed]
  • 62. Snyder A, Makarov V, Merghoub T, Yuan J, Zaretsky JM, Desrichard A, Walsh LA, Postow MA, Wong P, Ho TS, Hollmann TJ, Bruggeman C, Kannan K, et al. Genetic basis for clinical response to CTLA-4 blockade in melanoma. N Engl J Med. 2014; 371:2189–99. https://doi.org/10.1056/NEJMoa1406498 [PubMed]
  • 63. Rizvi NA, Hellmann MD, Snyder A, Kvistborg P, Makarov V, Havel JJ, Lee W, Yuan J, Wong P, Ho TS, Miller ML, Rekhtman N, Moreira AL, et al. Cancer immunology. Mutational landscape determines sensitivity to PD-1 blockade in non-small cell lung cancer. Science. 2015; 348:124–28. https://doi.org/10.1126/science.aaa1348 [PubMed]
  • 64. Cao D, Xu H, Xu X, Guo T, Ge W. High tumor mutation burden predicts better efficacy of immunotherapy: a pooled analysis of 103078 cancer patients. Oncoimmunology. 2019; 8:e1629258. https://doi.org/10.1080/2162402X.2019.1629258 [PubMed]
  • 65. Lee DW, Han SW, Bae JM, Jang H, Han H, Kim H, Bang D, Jeong SY, Park KJ, Kang GH, Kim TY. Tumor Mutation Burden and Prognosis in Patients with Colorectal Cancer Treated with Adjuvant Fluoropyrimidine and Oxaliplatin. Clin Cancer Res. 2019; 25:6141–47. https://doi.org/10.1158/1078-0432.CCR-19-1105 [PubMed]
  • 66. Pai SG, Carneiro BA, Chae YK, Costa RL, Kalyan A, Shah HA, Helenowski I, Rademaker AW, Mahalingam D, Giles FJ. Correlation of tumor mutational burden and treatment outcomes in patients with colorectal cancer. J Gastrointest Oncol. 2017; 8:858–66. https://doi.org/10.21037/jgo.2017.06.20 [PubMed]
  • 67. Zhou Z, Xie X, Wang X, Zhang X, Li W, Sun T, Cai Y, Wu J, Dang C, Zhang H. Correlations Between Tumor Mutation Burden and Immunocyte Infiltration and Their Prognostic Value in Colon Cancer. Front Genet. 2021; 12:623424. https://doi.org/10.3389/fgene.2021.623424 [PubMed]
  • 68. McGrail DJ, Pilié PG, Rashid NU, Voorwerk L, Slagter M, Kok M, Jonasch E, Khasraw M, Heimberger AB, Lim B, Ueno NT, Litton JK, Ferrarotto R, et al. High tumor mutation burden fails to predict immune checkpoint blockade response across all cancer types. Ann Oncol. 2021; 32:661–72. https://doi.org/10.1016/j.annonc.2021.02.006 [PubMed]
  • 69. Nakayama M, Oshima M. Mutant p53 in colon cancer. J Mol Cell Biol. 2019; 11:267–76. https://doi.org/10.1093/jmcb/mjy075 [PubMed]
  • 70. Taieb J, Le Malicot K, Shi Q, Penault-Llorca F, Bouché O, Tabernero J, Mini E, Goldberg RM, Folprecht G, Luc Van Laethem J, Sargent DJ, Alberts SR, Emile JF, et al. Prognostic Value of BRAF and KRAS Mutations in MSI and MSS Stage III Colon Cancer. J Natl Cancer Inst. 2016; 109:djw272. https://doi.org/10.1093/jnci/djw272 [PubMed]
  • 71. Powell SM, Zilz N, Beazer-Barclay Y, Bryan TM, Hamilton SR, Thibodeau SN, Vogelstein B, Kinzler KW. APC mutations occur early during colorectal tumorigenesis. Nature. 1992; 359:235–37. https://doi.org/10.1038/359235a0 [PubMed]
  • 72. Yu J, Wu WK, Li X, He J, Li XX, Ng SS, Yu C, Gao Z, Yang J, Li M, Wang Q, Liang Q, Pan Y, et al. Novel recurrently mutated genes and a prognostic mutation signature in colorectal cancer. Gut. 2015; 64:636–45. https://doi.org/10.1136/gutjnl-2013-306620 [PubMed]
  • 73. Zhou N, Bao J. FerrDb: a manually curated resource for regulators and markers of ferroptosis and ferroptosis-disease associations. Database (Oxford). 2020; 2020:baaa021. https://doi.org/10.1093/database/baaa021 [PubMed]
  • 74. Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, Benner C, Chanda SK. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun. 2019; 10:1523. https://doi.org/10.1038/s41467-019-09234-6 [PubMed]
  • 75. Ye Y, Hu Q, Chen H, Liang K, Yuan Y, Xiang Y, Ruan H, Zhang Z, Song A, Zhang H, Liu L, Diao L, Lou Y, et al. Characterization of Hypoxia-associated Molecular Features to Aid Hypoxia-Targeted Therapy. Nat Metab. 2019; 1:431–44. https://doi.org/10.1038/s42255-019-0045-8 [PubMed]
  • 76. Buffa FM, Harris AL, West CM, Miller CJ. Large meta-analysis of multiple cancers reveals a common, compact and highly prognostic hypoxia metagene. Br J Cancer. 2010; 102:428–35. https://doi.org/10.1038/sj.bjc.6605450 [PubMed]
  • 77. Chalmers ZR, Connelly CF, Fabrizio D, Gay L, Ali SM, Ennis R, Schrock A, Campbell B, Shlien A, Chmielecki J, Huang F, He Y, Sun J, et al. Analysis of 100,000 human cancer genomes reveals the landscape of tumor mutational burden. Genome Med. 2017; 9:34. https://doi.org/10.1186/s13073-017-0424-2 [PubMed]
  • 78. Jin M, Li D, Liu W, Wang P, Xiang Z, Liu K. Pinin acts as a poor prognostic indicator for renal cell carcinoma by reducing apoptosis and promoting cell migration and invasion. J Cell Mol Med. 2021; 25:4340–48. https://doi.org/10.1111/jcmm.16495 [PubMed]
  • 79. Jin M, Zhang H, Yang J, Zheng Z, Liu K. Expression mode and prognostic value of FXYD family members in colon cancer. Aging (Albany NY). 2021; 13:18404–22. https://doi.org/10.18632/aging.203290 [PubMed]