Skip to main content
Download PDF

Identification and validation of a prognostic model for melanoma patients with 9 ferroptosis-related gene signature

BMC Genomics volume 23, Article number: 245 (2022) Cite this article

Abstract

Background

Melanoma is a highly heterogeneous and aggressive cutaneous malignancy. Ferroptosis, a new pathway of cell death depending on the intracellar iron, has been shown to be significantly associated with apoptosis of a number of tumors, including melanoma. Nevertheless, the relationship between ferroptosis-related genes (FRGs) and the melanoma patients’ prognosis needs to be explored.

Methods

Download expression profiles of FRGs and clinical data from The Cancer Genome Atlas (TCGA) database. 70% data were randomly selected from the TCGA database and utilized the univariate Cox analysis and the least absolute shrinkage and selection operator (LASSO) regression model to create a prognostic model, and the remaining 30% was used to validate the predictive power of the model. In addition, GSE65904 and GSE22153 date sets as the verification cohort to testify the predictive ability of the signature.

Results

We identified nine FRGs relating with melanoma patients’ overall survival (OS) and established a prognostic model based on their expression. During the research, patients were divided into group of high-risk and low-risk according to the results of LASSO regression analysis. Survival time was significantly longer in the low-risk group than that of in the high-risk group (P < 0.001). Enrichment analysis of different risk groups demonstrated that the reasons for the difference were related to immune-related pathways, and the degree of immune cell infiltration in the low-risk group was significantly higher than that in the high-risk group.

Conclusions

The FRG prognostic model we established can predict the prognosis of melanoma patients and may further guide subsequent treatment.

Peer Review reports

Introduction

Melanoma is a malignant tumor around the world which is associate with rapid growth, early metastasis, local recurrence and poor prognosis [1]. According to the data of Cancer statistics in 2020, the incidence of melanoma accounted for 7% and 4% of male and female patients respectively, ranking fifth and sixth respectively [2]. Worldwide, there are about 23,100 new reported cases of melanoma every year and about 55,500 deaths due to melanoma ranked as the sixth most common malignancy in the US. Therefore, early detection and recognition of melanoma are keys to improve survival rate. At present, the treatment of the skin melanoma mainly includes surgical resection, which is the standard therapy for the primary melanoma, and other treatments applied to treat advanced melanoma such as radiotherapy, chemotherapy and immunotherapy [3]. Although thorough surgical resection of the tumor can greatly improve the five-year survival rate of patients with melanoma [4], the prognosis of patients are still not satisfactory, with a 5 year survival rate of 40 − 50% especially for the lymphatic metastasis and gene mutations [4]. Specific biomarkers play an important role in the early screening, diagnosis, and prognosis of melanoma. Therefore, it is urgent to establish a more sensitive prognostic models to assess the patient’s current condition for monitoring recurrence and evaluating prognosis.

Ferroptosis is a non-apoptotic form of cell death dependent on intracellular iron, which is different from apoptosis, necrosis and autophagy [5], The occurrence of ferroptosis is closely related to the inhibition of glutathione peroxidase4 (GPX4) synthesis. The weakening of GPX4-dependent antioxidant defense system eventually leads to the accumulation of lipid ROS which is toxic to cells and the depletion of polyunsaturated fatty acids, finally resulting in cell death [6]. Besides, ferroptosis has been proved to be related to the prognosis of tumor [7]. Ubellacker et al. found that melanoma cells in lymph nodes may be swollen through incorporate oleic acid and other antioxidants to protect themselves from ferroptosis [8]. This may be one of the reasons for the earlier lymphatic metastasis of melanoma mentioned above. Therefore, it can be speculated that the relationship between the expression of FRGs and the prognosis of patients remains to be further explored and studied.

By downloading the mRNA expression profiles of some melanoma patients from the Cancer Genome Atlas (TCGA) and their corresponding clinical data, and combining with the FRGs obtained from the original pubmed literatures, we constructed a prognosis model related to ferroptosis gene expression, and used the mRNA expression profiles of the remaining melanoma patients in the TCGA database combined with their clinical data to verify the correctness of the model. Finally, the functional enrichment and tumor microenvironment were analyzed to explore the possible mechanism.

Materials and methods

The research flowchart is displayed in Fig. 1.

Fig. 1
figure 1

The procession flow diagram in the present study

Full size image

Data Collection and Preprocessing

First of all, we downloaded the mRNA expression profiles and corresponding clinical information of 471 melanoma patients from the TCGA database (https://portal.gdc.cancer.gov/), after applying the “scale” function in the “limma” R package (version 4.0.3) to normalize data among databases, we combined them together. After that, removing incomplete clinical data and 0-days follow-up duration from samples, 455 melanoma samples obtained and used for the primary cohort. By using the method of generating random numbers by Microsoft Office Excel (version 2019), 455 samples are randomly divided into parts accounting for 70% of the total number of samples (n = 318) to be the training cohort and 30% of the total number of samples (n = 137) to be the validation cohort. Table S1 lists the clinical characteristics of the patients above. After browsing the previous literatures [9,10,11,12], 60 FRGs were obtained. These genes will be shown in Supplementary Table S2.

Identification of FRGs Affect Prognosis

After associating the expression levels of FRGs and melanoma patients’ overall OS in the training cohort, and trough univariate Cox analysis, 12 genes with prognostic significance (p < 0.05) were regarded as serving as an independent predictor for OS were obtained. These genes were used to form a prognostic model.

Construction and Validation of the Prognostic Model

Through the least absolute shrinkage and selection operator (LASSO) Cox regression analysis with package “glmnet” in R software and associated 12 genes above and patients’ survival data to find the most significant 9 genes and their corresponding Cox regression coefficient to build model. This model offered a formula to calculate risk score of each patient: risk score=esum (each gene’s expression×corresponding coefficient). According to the median value of the risk score in the training cohort with R software, patients in the training cohort were divided into high-risk group and low-risk group. After getting the formula obtained by establishing the model through training cohort, using this formula to calculate each patient’s risk score in validation cohort from TCGA database and divide them into high and low risk groups for the validation of the model. The t-distributed stochastic neighbor embedding (t-SNE) and principal component analysis (PCA) was analyzed with “Rtsne” package and the prcomp function in the “stats” package to explore the distribution of high and low-risk groups. By comping the survival between the two groups above and evaluate the model’s predictive ability using the “survivalROC” package in R respectively with “timeROC” package, Kaplan–Meier survival curves and a time-dependent receiver operating characteristic curve (ROC) curve analysis were employed. Besides, we downloaded the GSE65904 and GSE22153 datasets to verify the predictive power of the model, the two datasets are embedded in GPL10558 (Illumina HumanHT-12 V4.0 Expression BeadChip) and GPL6102 (Illumina Human-6 V2.0 Expression Beadchip) platforms, respectively. After downloading the Series Matrix file of the data set, the clinical information of patients was extracted, and the probes were replaced with the gene ID by the annotation information of the corresponding platform. The samples with incomplete clinical information were eliminated after the matrix of the combination of patient gene expression profile and clinical information, and two datasets with 202 and 49 patient samples were obtained, respectively. After using the RemoveBatchEffect function in “limma” R package (Version 4.0.3) to remove the batch effect between TCGA data set and Gene Expression Omnibus (GEO) repository (https://www.ncbi.nlm.nih.gov/geo), the “scale” function in “limma” R Package (Version 4.0.3) was used to normalize the data. Finally, calculating patient risk scores by the formula obtained in the train cohort and applying survival analysis and analysis of AUC in ROC for risk signature at 1-, 2- and 3- year survival time in the same way mentioned above.

Functional analysis

With the “clusterProfiler” R package, the analysis of Kyoto Encyclopedia of Genes, Gene Ontology (GO) and kyoto encyclopedia of genes and genomes (KEGG; www.kegg.jp/kegg/kegg1.html ) between the groups of high- and low-risk was enriched. In order to investigate the pathways enriched in the subgroups of high- and low-risk and explore possible molecular mechanisms. Using the false discovery rate (FDR) method to adjust P-values and while the P-value < 0.05, the pathways were considered to be enriched significantly.

Immune score, calculation of immune-related pathways and immune cells infiltration between two groups

By using the “estimate” package of R software in order to estimate the expression of immune and stromal cells in malignant tumor tissue and obtain immune stromal component ratio in the tumor microenvironment (TME) [13]. The result of the “estimate” R package created three scores to evaluate the presence of stroma (Stromal Score), the level of immune cells infiltrations (Immune Score), and the sum of stromal score and immune score (Estimate Score) [14]. The scores go up with the proportion of corresponding condition in the TME. The single-sample geneset enrichment analysis (ssGSEA) [15] was finished by using “GSVA” package of R software, the result of the level of the 13 immune-related pathways expression and16 kinds of immune-related cells are received. Supplementary Table S3 illustrates the immune-related genes.

Statistical analysis

Multivariate and Univariate Cox regression analyses were used to estimate if the factor can be regarded as an independent predictor. Utilizing R package “timeROC” to predict overall survival. Using t-test and chi-square test of student to identify the difference of Stromal Score, Immune Score, and ESTIMATE Score between patients in different risk groups. R software (Version 4.05) were applied here for all statistical analyses. Statistically significant were regarded in case when p-value less than 0.05.

Quantitative real-time PCR

Human normal skin cell line (TE353.sk) and human melanoma cell line (SKMEL5) (purchased from YaJiBiological, China) for verification are prepared after culture in accordance with the instructions provided by the manufacturer. Before the 2 −ΔΔCt statistics were applied to calculate the gene expression level in the final step, the total RNA of two kinds of cells above was extracted by TRIzol reagent (Invitrogen, China), and PrimeScript RT kit (Takara, China) was used for reverse transcription and the SYBR PrimeScript RT-PCR kit (Takara, China) was used for quantitative reverse transcriptase polymerase chain reaction (qRT-PCR) analysis. Supplementary Table S4 includes the primer sequences involved in this study.

Results

Identification of prognostic FRGs

Univariate Cox regression analysis showed that 12 genes could be used as independent factors of OS in patients with melanoma. Among the 12 genes, 6 of them were protective genes (HR<1) and 6 were risk genes (HR > 1) (Fig. 2).

Fig. 2
figure 2

Forest plots to demonstrate the univariate Cox regression analysis results between gene expression and the OS of melanoma patients

Full size image

Prognostic model construction in the training cohort

Nine genes which are significant subject to the OS of melanoma patients were obtained after the utilizing of LASSO Cox regression analysis (i.e., ACSL4, ALOX5, ATP5MC3, CHAC1, CS, MT1G, ACACA, ZEB1, and ABCC1). A ferroptosis-related prognostic model was created on account of the best value of λ from LASSO Cox regression analysis. The score of risk can be figured out though the following formula: risk score=e (0.05785 * expression level of ACACA − 0.14149 * expression level of ACSL4 − 0.06033 * expression level of ALOX5 + 0.04117 * expression level of ATP5MC3 + 0.02936 * expression level of CHAC1 + 0.25160 * expression level of CS − 0.00617 * expression level of MT1G − 0.09239 * expression level of ZEB1 + 0.18878 * expression level of ABCC1). The risk values of patients in the train cohort were arranged from high to low, and the median value was found as the cut-off value to divide patients into two groups of low-risk (n = 159) and high -risk (n = 159) (Fig. 3(A)). PCA and t-SNE analysis showed the result that the patients in low-risk and high-risk groups were split into two directions (Fig. 3(B) and (C)). Figure 3(D) indicated that patients of high-risk group own a poor survival. The Kaplan–Meier survival analysis also confirmed the survival time of high-risk group was yielding reduced (Fig. 3(E), P < 0.001). An evaluation of the predictive performance was evaluated with the time-dependent ROC curves of the model and in Fig. 3(F), the area under curve (AUC) reached 0.598 at 1 year, 0.687 at 2 years, and 0.664 at 3 years.

Fig. 3
figure 3

The 9-gene model in the train cohort prognostic analysis. A The risk score distribution in the train cohort. B the train cohort PCA analysis. C t-SNE plot the train cohort. D The distribution of OS in the train cohort. E The Kaplan–Meier survival OS analysis between the high-risk and low-risk group of the train cohort. F The analysis of AUC in ROC for risk signature at 1-, 2- and 3-year survival time in the train cohort

Full size image

The prognostic model validation in the validation cohort

Calculating the score of risk for each patient in the validation cohort by using the formula obtained above to verify the accuracy of the prediction ability of the model. Using the method mentioned above, patients are split into high- (n = 68) and low-risk (n = 68) groups in the validation cohort according to risk score and median (Fig. 4(A)). In the validation cohort, a reliable clustering ability of risk score are also validated by PCA and t-SNE analysis (see Fig. 4(B) and 4(C)). It is obvious that the survival time of high-risk group patients is shorter than that of low-risk (Fig. 4(D), P < 0.001). Meanwhile, the Kaplan–Meier survival OS analysis between the high-risk group and low-risk group in the validation cohort shows there were also significant differences in survival between high and low-risk groups. (Figure 4(E)) and the values of AUC of year 1, 2 and 3 after ROC analysis are 0.759, 0.676 and 0.639 (Fig. 4(F)).

Fig. 4
figure 4

The prognostic model validation in the validation cohort. A The risk score distribution of in the validation cohort. B the validation cohort PCA analysis. C the validation cohort T-SNE plot. D The OS distribution in the validation cohort. E The Kaplan–Meier survival OS analysis between the high-risk group and low-risk group in the validation cohort. F The risk signature AUC in ROC analysis for 1-, 2- and 3-year survival time in the validation cohort

Full size image

Independent prognostic value of the risk score

In order to further explore whether risk score is an independent prognostic factor, we conducted univariate and multivariate Cox regression analysis. In the train cohort (HR= 4.383, 95% CI= 2.581 − 7.445, P < 0.001, Fig. 5(A)) and the validation cohort (HR= 6.672, 95% CI= 2.661 − 16.733, P < 0.001, Fig. 5(B)), univariate Cox regression analysis showed that the both risk score had a significant relationship with OS. Multivariate Cox regression analysis corrected for potential confounding factors also showed that risk score could be used as an independent predictor of OS (train cohort: HR= 4.429, 95% CI= 2.610 − 7.516, P < 0.001; validation cohort: HR= 7.032, 95% CI= 2.858 − 17.298, P < 0.001, Fig. 5(C) and 5(D)). In order to prove that our prognostic model is not a different way to quantify immune infiltration, we include the ESTIMATE score into the multivariate regression analysis model. The results show that the FRG based risk score still has predictive ability (P < 0.001 in train cohort and P < 0.05 in test cohort, Figure S1A and Figure S1B) although with a diminished hazard ratio.

Fig. 5
figure 5

Independent prognostic signature presented by risk score for melanoma. The univariate and multivariate Cox regression OS analyses results in the validation cohort (B/D) and in the train cohort (A/C)

Full size image

Functional analyses

Enrichment and KEGG pathway were carried out according to the differential expression genes (DEGs) between high-risk group and low-risk group in the train cohort and validation cohort to clarify the biological function and pathway related to risk score. According to the GO enrichment analysis result, the DEGs between risk groups from validation cohort and the train cohort were primarily enriched in immune response−activating cell surface receptor so as to signal pathway and immune response−activating signal transduction. (P. adjust <0.05, Fig. 6(A), (B)). It is can be seen that the cytokine−cytokine receptor interaction pathway was significantly enriched in both cohorts (P. adjust <0.05, Fig. 6(C), (D)) in KEGG pathway analysis.

Fig. 6
figure 6

The results of functional analyses. The barplot graph for the result of GO enrichment analysis for the train cohort (A) and the validation cohort (B). The barplot graph for the result of KEGG pathway analysis for the train cohort (C) and the validation cohort (D)

Full size image

The Immune-related pathways and immune cells infiltration analysis

In order to explore the difference of immune-related function between risk groups of high and low, each patient was scored by ssGSEA analysis, and then analyzed the enrichment of immune-related pathways and immune cells by the score of each patient. In the train cohort, B cells, D8 + T cells, NK cells, Neutrophils, pDCs, T helper cells, Th1 cells, Th2 cells, Tfh, TIL, T reg along with other immune cell subsets between the two groups were significantly different and down-regulated in the immune-related pathway, in the high-risk group (adjusted P < 0.05, see Fig. 7 (A)), Cytolytic activity, Inflammation-promoting was significantly down-regulated In the validation cohort, in the high-risk group (adjusted P < 0.05, Fig. 7 (B). the enrichment of Macrophage, Dendritic Cells (DCs) and Mast cells was different from that of train group, and the enrichment of remaining immune cells and the pathways of immune-related function was similar to that of train cohort (corrected P < 0.05, Fig. 7 (C) and (D)).

Fig. 7
figure 7

Differences of immune-related pathways and immune cells infiltration between groups of high- risk and low- risk. ssGSEA scores of 16 immune cells (A) and ssGSEA scores of 13 immune-related functions (B) in the train group. ssGSEA scores of 16 immune cells (C) and ssGSEA scores of 13 immune-related functions (D) in the validation group. ns: not significant; ∗: P < 0.05;∗∗: P < 0.01; ∗∗∗: P < 0.001

Full size image

Verification of risk score model in GEO repository

Data from GSE65904 and GSE22153 are used as validation cohorts 1 and validation cohorts 2. By using our FRG based riskScore, the samples of the two verification cohorts are divided into high and low risk groups. The Kaplan–Meier survival analysis shows the difference in survival time between low and high-risk groups still significant in two verification cohorts. (Figure 8(A), P < 0.001 and (B), P < 0.05). The result of time-dependent ROC curves was display in Fig. 8(C/D) to evaluate the predictive performance of the model. The area under curve (AUC) reached 0.701 at 1 year, 0.690 at 2 years, and 0.677 at 3 years for validation cohorts 1 and 0.977 at 1 year, 0.711 at 2 years, and 0.685 at 3 years for validation cohorts 2.

Fig. 8
figure 8

The Kaplan–Meier survival OS analysis between the high-risk and low-risk group of the validation cohort 1 (A) and validation cohort 2 (B). The analysis of AUC in ROC for risk signature at 1-, 2- and 3-year survival time in the validation cohort 1 (C) and validation cohort 2 (D)

Full size image

Genes expression levels in cell lines validation

In the 9 gene signatures, by qRT-PCR, we found that significant differences in the expression of ABCC1, ACSL4 and ALOX5 between normal skin cell lines and melanoma cell lines. Among them, the expression of ACSL4 and ALOX5 was up-regulated in normal skin cell lines (Fig. 9 (A) and (B), P < 0.05), along with the higher expression level of ABCC1 in melanocyte lines (Fig. 8 (C), P < 0.05).

Fig. 9
figure 9

Differences of Genes Expression Levels in Normal Skin and melanoma Cell Lines Validated by qRT-PCR. Relative expression levels of A ACSL4, B ALOX5, and C ABCC1 between normal skin cell line and melanoma cell line. ns: not significant; ∗P < 0.05; ∗∗P < 0.01; ∗∗∗P < 0.001; ∗∗∗∗P < 0.0001

Full size image

Discussion

In this study, we used the patients whose mRNA expression profiles and clinical data of 70% melanoma in the TCGA database as the train cohort and combined with 60 FRGs to construct a model containing 9 FRGs to predict the prognosis of patients of melanoma. The remaining 30% patients, data were used as the verification group to verify the model predictive function.

In the above two cohorts, the patients’ survival time in the high-risk group was significantly shorter than that in the group of low-risk. Analysis with Functional enrichment showed that the differences were mainly because of cytokine-cytokine receptor interaction and immune response-related pathways. Besides, immune-related functions and immune cell infiltration in the group of high-risk were remarkably lower than those in the group of low-risk. A signature model consisting of 9 ferroptosis genes, acyl-CoA synthetase long-chain family member 4(ACSL4), 5-Lipoxygenase (ALOX5), metallothionein (MT)-1G and Zinc finger E-box-binding homeobox 1 (ZEB1) as protective genes, ATP synthase membrane subunit c locus 3 (ATP5MC3), ChaC glutathione specific gamma-glutamylcyclotransferase 1 (CHAC1), citrate synthase (CS), acetyl-CoA carboxylase alpha (ACACA), ATP binding cassette subfamily C member 1 (ABCC1) as risk genes. Ferroptosis is a regulatory necrotic cell death controlled by glutathione peroxidase 4(GPX4), and the overexpression of ACSL4 will reduce the expression of GPX4. Moreover, ACSL4, a member of the long-chain acyl-CoA synthase family, can induce ferroptosis by oxidizing arachidonic acid [16, 17]. Therefore, ACSL4, a proferroptotic gene, plays a crucial role in cells ferroptosis process. Sebastian Doll et al. found that ACSL4 knockout cells showed significant resistance to ferroptosis, and reexpression of ACSL4 enabled cells to regain sensitivity to ferroptosis [18]. Jing Cheng et al. observed that solafenib increased cell viability by reducing siRNA-mediated ACSL4 silencing, suggesting that ACSL4 may protect glioma cells and inhibit their proliferation via activating a ferroptosis pathway [16]. It can be concluded that the high expression level of ACSL4 can lead more tumor cells to ferroptosis. As shown in our prognostic model, the HR value of ACSL4 is less than 1, which means that it is a protective gene, the higher its expression level brings the better the patients’ possible prognosis. It is worth mentioning that in José Pedro et al. pointed out that ASCL4 may activate and attract immune cells to clear ferroptosis tumor cells by sending signals such as “find me” and “eat me” to immune cells by participating in ferroptosis tumor cells, and put forward the related conjecture that the loss of these signals may lead to immune evasion of tumor cells [19]. ALOX5 is an iron-containing non-heme dioxygenase, which seems a key enzyme in the synthesis of leukotriene and also can be used to catalyze the peroxidation of polyunsaturated fatty acids [20]. Ferroptosis could be irritated by lipid peroxidation to mediate inflammation-related cell death [21], and ALOX5 plays a crucial role in both wire death and inflammation [22]. Previous study found that inhibition of ALOX5 expression could decrease ferroptosis in nerve cells derived from hemorrhagic stroke mice [22]. Faronato et al. [23] and Miess et al. [24] thought that the expression of ALOX5 in clear cell renal cell carcinomas deficient in von Hippel-Lindau (VHL) gene was greatly increased, which may be that this kind of cancer cell needs more eicosanoids synthesized through ALOX5 expression to promote local inflammatory response. In the prognostic model of melanoma patients with FRGs created by Xu and Chen [25], ALOX5 also formed their 5 gene signature as an indispensable member, and it was confirmed that its expression can be used as an independent prognostic factor to predict the OS of patients. Metallothioneins (MTS) is a protein that is highly expressed under the influence of different environmental stressors and is closely related to heavy metal detoxifications and antioxidants. The up-regulation of MT-1G expression can increase the drug resistance of tumor cells to sorafenib by inhibiting ferroptosis [26]. However, Sun, X et al. found that sorafenib can activate Nrf2 through the cystathionase pathway leading to the expression of MT-1G in liver cells. Nevertheless, inhibiting the expression of MT-1G will enhance the metastatic tumor activity of sorafenib against hepatoma [27]. In conclusion, MT-1G may become a new target for anti-tumor therapy in the future. ZEB1 is known as an epithelial marker to down-regulate the expression of e-cadherin to affect the epithelial-mesenchymal transition and thus participate in the invasion and metastasis of tumor cells [28]. Lee and his team tested the results of ferroptosis inducers induction in HNC cell lines or EMT inhibition and rat neoplasm graft models, found that overexpression of ZEB1 may enhance the neoplasm cells sensitivity to ferroptosis, and it absolutely was additionally confirmed in animal models that the neoplasm volume within the overexpression group was considerably reduced compared with the control group [29]. ATP5G3 encodes a fractional monetary unit of mitochondrial membrane ATP synthase, that catalyzes ATP synthesis throughout organic process. ATP5G3 was found to be up-regulated three days when the prevalence of secondary craniocerebral injury to accelerate ferroptosis of cells [30]. CHAC1, a γ-glutamyl cyclotransferase, which inhibits ferroptosis by enhancing the degradation of glutathione [31, 32]. In the research and exploration of Wang et al., finally found that Artesunate may increase the expression of CHAC1 through the ATF4-CHOP pathway, thereby rising the sensitivity of Burkitt’s cancer cells to ferroptosis [33]. CS will catalyse the synthesis of change state from oxalacetate, And gives material for carboxylic acid synthesis thus on provide needed macromolecule precursors for macromolecule peroxidation caused by ferroptosis [34]. Erastin is one of the small molecules that can induce ferroptosis [35]. Dixon et al. found that silencing CS can significantly reduce ferroptosis induced by Erastin [5]. ACACA mainly acts on the first stage of fatty acid synthesis, and is one of the rate-limiting enzymes that regulate fat and metabolism, and plays a vital role in the tumor cells survival [36]. In addition, knockout of ACACA could inhibit drug-induced cell ferroptosis [37]. Meanwhile, AMPK pathway could be activated to inhibit its downstream ACACA, subsequently slowing down lipid accumulation and ferroptosis [38]. The expression of ABCC1 can be positively regulated by the antioxidant transcription factor Nrf2 [39] to regulate the process of cell ferroptosis. Cao et al. down-regulated the expression of ATP binding cassette (ABC)-family transporter multidrug resistance protein 1 (MRP1). This prevents glutathione from flowing out of the cells and effectively inhibits ferroptosis [40]. In melanoma, the synergistic effect of ABCC1 and glutathione S-transferase M1 can also make tumor cells resistant to vincristine [41].

GO and KEGG enrichment analysis showed that the difference in expression between groups was mainly related to tumor microenvironment, and the enrichment of signaling pathway of immune response-activating cell surface receptor, transduction of immune response-activating signal and cytokine-cytokine receptor interaction pathway were the most significant. We can find that cAMP signal pathway is highly enriched in these differentially expressed genes. Arumugham et al. found that cAMP can regulate T cell activation and immune synaptic assembly to regulate the immune process [42]. These three immune-related pathways have also been found to be enriched in lung adenocarcinoma [43], testicular cancer [44], glioblastoma [45], and other tumors and can be used to predict the prognosis of patients. Although the specific mechanism of their role in melanoma needs to be further studied, we think that they have great potential as indicators for predicting the prognosis of patients with melanoma. We found that there was a significant difference in immune cell infiltration between the group of high- and low-risk. The immune cell infiltration in the group of low-risk was higher than that in the group of high-risk, including CD8+T cell and DCs. DCs were found to accumulate large amounts of lipids and polyunsaturated fatty acids in tumor patients, which led to a decline in their ability to present antigens and unable to fully stimulate activated CD8+T cell [46, 47]. In vitro experiments conducted by Matsushita et al. also confirmed that CD8+T cells can increase the specific lipid peroxidation of ferroptosis by releasing interferon-γ and increase the occurrence of ferroptosis, thus improving the effect of immunotherapy [48]. This phenomenon may be due to the releasing signal molecules such as interferon-γ by DCs and CD8+T cells to activate ferroptosis [49]. The level of immune cell infiltration in the high-risk group is significantly lower than that in the low-risk group, indicating that the tumor cells in the high-risk group are less sensitive to ferroptosis. In accordance, Friedmann Angeli et al. also found that tumor cells with ferroptosis release arachidonic acid as an immune activator to obtain drug resistance and immune evasion from other tumor cells [19].

In this study, the risk value calculated after the establishment of the model is a reliable independent prognostic index. Compared with the conventional prognostic indexes such as “age”, “gender”, “tumor stage” and “ESTIMATE score”, the risk score created by the above nine gene expression can better predict the survival of patients, which also confirms that the gene-based expression signal can accurately predict the prognosis of patients with melanoma. However, this study also has some limitations. First of all, the clinical information in TCGA database is incomplete, especially the lack of treatment-related information. Secondly, the prognostic prediction model constructed in this study is based on retrospective data, and no prospective clinical studies have been carried out to verify the model. Garg et al. proposed a prognostic signature consisting of 121 metastasis-related genes to predict the prognosis of patients with melanoma [50]. But unfortunately, none of them appeared in both their and our signatures. Ubellacker et al. found that a high level of ferroptosis inhibitor glutathione peroxidase 4 (GPX4) may be the reason for the earlier lymphatic metastasis of melanoma, but interestingly it appears also not in the metastasis-related gene signature constructed by Garg et al. Whether the genes or ferroptosis-related genes in our signature are related to the risk of metastasis is the focus of our next work and research. Since our prognostic model is based on ferroptosis-related genes, as mentioned in our Discussion section above, the ferroptosis behavior of tumors is closely related to immune infiltration and immune pathway, there may be some overlap between ESTIMATE score and riskScores to evaluate immunity. More specific forms of interaction between the two may need to be further verified and studied in other databases. However, at this stage, we still believe that the nine genes signature proposed in this research has the capacity to accurately predict the patients’ prognosis with melanoma and provide a new direction for new treatment strategies.

Conclusions

In conclusion, our research shows that the ferroptosis genes expression is related with the progression of patients of melanoma. We combined the ferroptosis gene expression of patients with clinical data to construct a signature containing 9 genes to accurately forecast the melanoma patients’ prognosis of. The research also revealed the close relationship between immune function and ferroptosis related genes through analysis of immune cell infiltration and immune-related functions.

Availability of data and materials

The datasets generated and analysed during the current study are available in the TCGA repository. https://tcga-data.nci.nih.gov/tcga. The independent validation cohort data is available in the GEO repository. https://www.ncbi.nlm.nih.gov/geo( GSE65904 and GSE22153).

Abbreviations

FRGs:

Ferroptosis-related genes

TCGA:

The cancer genome atlas

LASSO:

Least absolute shrinkage and selection operator

OS:

Overall survival

GPX4:

Glutathione peroxidase 4

PCA:

Principal component analysis

t-SNE:

T-distributed stochastic neighbor embedding

ROC:

Receiver operating characteristic curve

GO:

Gene ontology

KEGG:

Kyoto encyclopedia of genes and genomes

FDR:

False discovery rate

TME:

Tumor microenvironment

ssGSEA:

Single-sample geneset enrichment analysis

qRT-PCR:

Quantitative reverse transcriptase polymerase chain reaction

DEGs:

Differential expression genes

DCs:

Dendritic Cells; ACSL4:acyl-CoA synthetase long-chain family member 4

ALOX5:

5-Lipoxygenase

MT:

Metallothionein

ZEB1:

Zinc finger E-box-binding homeobox 1

ATP5MC3:

ATP synthase membrane subunit c locus 3

CHAC1:

ChaC glutathione specific gamma-glutamylcyclotransferase 1

CS:

Citrate synthase

ACACA:

Acetyl-CoA carboxylase alpha

ABCC1:

ATP binding cassette subfamily C member 1

MTs:

Metallothioneins

ABC:

ATP binding cassette

MRP1:

Multidrug resistance protein 1

GPX4:

Glutathione peroxidase4

References

  1. Gagliardi M, Saverio V, Monzani R, Ferrari E, Piacentini M, Corazzari M. Ferroptosis: a new unexpected chance to treat metastatic melanoma? Cell Cycle. 2020;19(19):2411–25.

    Article  CAS  Google Scholar 

  2. Siegel RL, Miller KD, Jemal A. Cancer statistics, 2020. CA Cancer J Clin 2020;70(1):7–30.

    Article  Google Scholar 

  3. Kozovska Z, Gabrisova V, Kucerova L. Malignant melanoma: diagnosis, treatment and cancer stem cells. Neoplasma. 2016;63(4):510–7.

    Article  CAS  Google Scholar 

  4. Schadendorf D, van Akkooi ACJ, Berking C, et al. Melanoma. Lancet. 2018;392(10151):971–84.

    Article  Google Scholar 

  5. Dixon SJ, Lemberg KM, Lamprecht MR, et al. Ferroptosis: an iron-dependent form of nonapoptotic cell death. Cell. 2012;149(5):1060–72.

    Article  CAS  Google Scholar 

  6. Cao JY, Dixon SJ. Mechanisms of ferroptosis. Cell Mol Life Sci. 2016;73(11–12):2195–209.

    Article  CAS  Google Scholar 

  7. Liu Z, Zhao Q, Zuo ZX, et al. Systematic Analysis of the Aberrances and Functional Implications of Ferroptosis in Cancer. iScience. 2020;23(7):101302.

  8. Ubellacker JM, Tasdogan A, Ramesh V, et al. Lymph protects metastasizing melanoma cells from ferroptosis. Nature.2020;585(7823):113-118.

  9. Stockwell BR, Friedmann Angeli JP, Bayir H, et al. Ferroptosis: A Regulated Cell Death Nexus Linking Metabolism, Redox Biology, and Disease. Cell. 2017;171(2):273–85.

    Article  CAS  Google Scholar 

  10. Hassannia B, Vandenabeele P, Vanden Berghe T. Targeting Ferroptosis to Iron Out Cancer. Cancer Cell. 2019;35(6):830–49.

    Article  CAS  Google Scholar 

  11. Bersuker K, Hendricks JM, Li Z, et al. The CoQ oxidoreductase FSP1 acts parallel to GPX4 to inhibit ferroptosis. Nature. 2019;575(7784):688–92.

    Article  CAS  Google Scholar 

  12. Doll S, Freitas FP, Shah R, et al. FSP1 is a glutathione-independent ferroptosis suppressor. Nature. 2019;575(7784):693–8.

    Article  CAS  Google Scholar 

  13. Jiang X, Yan Q, Xie L, et al. Construction and Validation of a Ferroptosis-Related Prognostic Model for Gastric Cancer. J Oncol. 2021;2021:6635526.

    PubMed  PubMed Central  Google Scholar 

  14. Yoshihara K, Shahmoradgoli M, Martinez E, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. 2013;4:2612.

    Article  Google Scholar 

  15. Rooney MS, Shukla SA, Wu CJ, Getz G, Hacohen N. Molecular and genetic properties of tumors associated with local immune cytolytic activity. Cell. 2015;160(1–2):48–61.

    Article  CAS  Google Scholar 

  16. Cheng J, Fan YQ, Liu BH, Zhou H, Wang JM, Chen QX. ACSL4 suppresses glioma cells proliferation via activating ferroptosis. Oncol Rep. 2020;43(1):147–58.

    CAS  PubMed  Google Scholar 

  17. Yuan H, Li X, Zhang X, Kang R, Tang D. Identification of ACSL4 as a biomarker and contributor of ferroptosis. Biochem Biophys Res Commun. 2016;478(3):1338–43.

    Article  CAS  Google Scholar 

  18. Doll S, Proneth B, Tyurina YY, et al. ACSL4 dictates ferroptosis sensitivity by shaping cellular lipid composition. Nat Chem Biol. 2017;13(1):91–8.

  19. Friedmann Angeli JP, Krysko DV, Conrad M. Ferroptosis at the crossroads of cancer-acquired drug resistance and immune evasion. Nat Rev Cancer. 2019;19(7):405–14.

    Article  CAS  Google Scholar 

  20. Sun QY, Zhou HH, Mao XY. Emerging Roles of 5-Lipoxygenase Phosphorylation in Inflammation and Cell Death. Oxid Med Cell Longev. 2019;2019:2749173.

    PubMed  PubMed Central  Google Scholar 

  21. Liu Y, Xu Z, Jin T, Xu K, Liu M, Xu H. Ferroptosis in Low-Grade Glioma: A New Marker for Diagnosis and Prognosis. Med Sci Monit. 2020;26.

  22. Karuppagounder SS, Alin L, Chen Y, et al. N-acetylcysteine targets 5 lipoxygenase-derived, toxic lipids and can synergize with prostaglandin E2 to inhibit ferroptosis and improve outcomes following hemorrhagic stroke in mice. Ann Neurol. 2018;84(6):854–72.

    Article  CAS  Google Scholar 

  23. Faronato M, Muzzonigro G, Milanese G, et al. Increased expression of 5-lipoxygenase is common in clear cell renal cell carcinoma. Histol Histopathol. 2007;22(10):1109-1118.

  24. Miess H, Dankworth B, Gouw AM, et al. The glutathione redox system is essential to prevent ferroptosis caused by impaired lipid metabolism in clear cell renal cell carcinoma. Oncogene. 2018;37(40):5435-5450.

  25. Xu C, Chen H. A Ferroptosis-Related Gene Model Predicts Prognosis and Immune Microenvironment for Cutaneous Melanoma. Front Genet. 2021;12:697043.

  26. Xie Y, Hou W, Song X, et al. Ferroptosis: process and function. Cell Death Differ. 2016;23(3):369-379.

  27. Sun X, Niu X, Chen R, et al. Metallothionein-1G facilitates sorafenib resistance through inhibition of ferroptosis. Hepatology. 2016;64(2):488-500.

  28. Sun L, Kokura K, Izumi V, et al. MPP8 and SIRT1 crosstalk in E-cadherin gene silencing and epithelial-mesenchymal transition. EMBO Rep. 2015;16(6):689-699.

  29. Lee J, You JH, Kim MS, Roh JL. Epigenetic reprogramming of epithelial-mesenchymal transition promotes ferroptosis of head and neck cancer. Redox Biol. 2020;37:101697.

  30. Huang L, He S, Cai Q, et al. Polydatin alleviates traumatic brain injury: Role of inhibiting ferroptosis. Biochem Biophys Res Commun. 2021;556:149-155.

  31. Li D, Liu S, Xu J, et al. Ferroptosis-related gene CHAC1 is a valid indicator for the poor prognosis of kidney renal clear cell carcinoma. J Cell Mol Med. 2021;25(7):3610-3621.

  32. Crawford RR, Prescott ET, Sylvester CF, et al. Human CHAC1 Protein Degrades Glutathione, and mRNA Induction Is Regulated by the Transcription Factors ATF4 and ATF3 and a Bipartite ATF/CRE Regulatory Element. J Biol Chem. 2015;290(25):15878-15891.

  33. Wang N, Zeng GZ, Yin JL, Bian ZX. Artesunate activates the ATF4-CHOP-CHAC1 pathway and affects ferroptosis in Burkitt’s Lymphoma. Biochem Biophys Res Commun. 2019;519(3):533-539.

  34. Wang H, Liu C, Zhao Y, Gao G. Mitochondria regulation in ferroptosis. Eur J Cell Biol. 2020;99(1):151058.

  35. Ratan RR. The Chemical Biology of Ferroptosis in the Central Nervous System. Cell Chem Biology. 2020;27(5):479-498.

  36. Brownsey RW, Boone AN, Elliott JE, Kulpa JE, Lee WM. Regulation of acetyl-CoA carboxylase. Biochem Soc Trans 2006;34(Pt 2):223-227.

  37. Dixon SJ, Winter GE, Musavi LS, et al. Human Haploid Cell Genetics Reveals Roles for Lipid Metabolism Genes in Nonapoptotic Cell Death. ACS Chem Biol 2015;10(7):1604-1609.

  38. Li C, Dong X, Du W, et al. LKB1-AMPK axis negatively regulates ferroptosis by inhibiting fatty acid synthesis. Signal Transduct Target Ther. 2020;5(1):187.

  39. Ji L, Li H, Gao P, et al. Nrf2 pathway regulates multidrug-resistance-associated protein 1 in small cell lung cancer. PLoS One. 2013;8(5):e63404.

  40. Cao JY, Poddar A, Magtanong L, et al. A Genome-wide Haploid Genetic Screen Identifies Regulators of Glutathione Abundance and Ferroptosis Sensitivity. Cell Rep. 2019;26(6):1544-1556 e1548.

  41. Depeille P, Cuq P, Passagne I, Evrard A, Vian L. Combined effects of GSTP1 and MRP1 in melanoma drug resistance. Brit J Cancer. 2005;93(2):216–23.

    Article  CAS  Google Scholar 

  42. Arumugham VB, Baldari CT. cAMP: a multifaceted modulator of immune synapse assembly and T cell activation. J Leukocyte Biol. 2017;101(6):1301–16.

    Article  CAS  Google Scholar 

  43. Zhang M, Zhu K, Pu H, et al. An Immune-Related Signature Predicts Survival in Patients With Lung Adenocarcinoma. Front Oncol. 2019;9:1314.

    Article  Google Scholar 

  44. Ke ZB, Wu YP, Huang P, et al. Identification of novel genes in testicular cancer microenvironment based on ESTIMATE algorithm-derived immune scores. J Cell Physiol. 2021;236(1):706–13.

    Article  CAS  Google Scholar 

  45. Zhou L, Tang H, Wang F, et al. Bioinformatics analyses of significant genes, related pathways and candidate prognostic biomarkers in glioblastoma. Mol Med Reports. 2018;18(5):4185–96.

    CAS  Google Scholar 

  46. Ramakrishnan R, Tyurin VA, Veglia F, et al. Oxidized lipids block antigen cross-presentation by dendritic cells in cancer. J Immunol. 2014;192(6):2920–31.

    Article  Google Scholar 

  47. Veglia F, Tyurin VA, Mohammadyani D, et al. Lipid bodies containing oxidatively truncated lipids block antigen cross-presentation by dendritic cells in cancer. Nat Commun. 2017;8(1):2122.

    Article  Google Scholar 

  48. Matsushita M, Freigang S, Schneider C, Conrad M, Bornkamm GW, Kopf M. T cell lipid peroxidation induces ferroptosis and prevents immunity to infection. J Experiment Med. 2015;212(4):555–68.

    Article  CAS  Google Scholar 

  49. Wang W, Green M, Choi JE, et al. CD8(+) T cells regulate tumour ferroptosis during cancer immunotherapy. Nature. 2019;569(7755):270–4.

    Article  CAS  Google Scholar 

  50. Garg M, Couturier DL, Nsengimana J, et al. Tumour gene expression signature in primary melanoma predicts long-term outcomes. Nat Commun. 2021;12(1):1137.

    Article  CAS  Google Scholar 

Download references

Acknowledgements

Not applicable.

Funding

National Natural Science Foundation of China (grant 81601701 to Dr J. Wang and 81901977 to Dr R. An) provides support on this work.

Author information

Author notes

  1. Yuxuan Chen and Linlin Guo contributed equally to this work.

Authors and Affiliations

  1. Department of Plastic Surgery, Union Hospital, Tongji Medical College, Huazhong University of Science and Technology, Wuhan, PR China

    Yuxuan Chen, Linlin Guo, Ran An & Jiecong Wang

  2. Department of Neurosurgery, Union Hospital, Tongji Medical College, Huazhong University of Science and Technology, Wuhan, PR China

    Zijie Zhou

Authors
  1. Yuxuan Chen
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Linlin Guo
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Zijie Zhou
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Ran An
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Jiecong Wang
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

YXC provide the idea of this manuscript, YXC and LLG wrote the manuscript, ZJZ drew the figures and tables, RA and JCW made some revisions of the review. All the authors read and approved the final version of the review.

Corresponding authors

Correspondence to Ran An or Jiecong Wang.

Ethics declarations

Ethics approval and consent to participate

The data of our study came from open access database and do not involve data of human participants. All methods were performed in accordance with the relevant guide lines and regulations.

Consent for publication

Not applicable.

Competing interests

The authors have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Chen, Y., Guo, L., Zhou, Z. et al. Identification and validation of a prognostic model for melanoma patients with 9 ferroptosis-related gene signature. BMC Genomics 23, 245 (2022). https://doi.org/10.1186/s12864-022-08475-y

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12864-022-08475-y

Keywords

BMC Genomics

ISSN: 1471-2164

Contact us